8 implicit none ;
private
20 real,
dimension(:,:),
intent(inout) :: a
21 real,
dimension(:),
intent(inout) :: b
22 real,
dimension(:),
intent(inout) :: x
23 integer,
intent(in) :: system_size
26 real,
parameter :: eps = 0.0
29 real :: swap_a, swap_b
30 logical :: found_pivot
33 do i = 1,system_size-1
43 do while ( ( .NOT. found_pivot ) .AND. ( k <= system_size ) )
45 if ( abs( a(k,i) ) > eps )
then
56 if ( .NOT. found_pivot )
then
58 call mom_error( fatal,
'The linear system is singular !' )
78 a(i,j) = a(i,j) / pivot
86 do k = (i+1),system_size
88 do j = (i+1),system_size
89 a(k,j) = a(k,j) - factor * a(i,j)
91 b(k) = b(k) - factor * b(i)
98 x(system_size) = b(system_size) / a(system_size,system_size)
99 do i = system_size-1,1,-1
101 do j = (i+1),system_size
102 x(i) = x(i) - a(i,j) * x(j)
114 real,
dimension(:),
intent(inout) :: ad
115 real,
dimension(:),
intent(inout) :: al
116 real,
dimension(:),
intent(inout) :: au
117 real,
dimension(:),
intent(inout) :: b
118 real,
dimension(:),
intent(inout) :: x
119 integer,
intent(in) :: system_size
128 al(k+1) = al(k+1) / ad(k)
129 ad(k+1) = ad(k+1) - al(k+1) * au(k)
134 b(k) = b(k) - al(k) * b(k-1)
140 x(k) = ( b(k) - au(k)*x(k+1) ) / ad(k)