Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer | :: | N | ||||
complex(kind=8), | dimension(N-1) | :: | sub_ | |||
complex(kind=8), | dimension(N) | :: | diag_ | |||
complex(kind=8), | dimension(N-1) | :: | over_ | |||
complex(kind=8), | dimension(N) | :: | Inv |
subroutine c_invert_tridiag_matrix(N,sub_,diag_,over_,Inv) integer :: i,j,N complex(8),dimension(N) :: diag_ complex(8),dimension(N-1) :: sub_ complex(8),dimension(N-1) :: over_ complex(8),dimension(N) :: Inv complex(8) :: Foo,Cleft,Cright complex(8),dimension(N) :: dleft complex(8),dimension(N) :: dright ! !DOWNWARD: dleft(1) = diag_(1) if(dleft(1)==0d0)stop "matrix is ill-conditioned: no inverse exists" do i=2,N foo = 1d0/dleft(i-1) cleft = sub_(i-1)*foo dleft(i) = diag_(i) - cleft*over_(i-1) !over_(i-1)/dleft(i-1)*sub_(i) if(dleft(i)==0d0)stop "matrix is ill-conditioned: no inverse exists" enddo ! !UPWARD: dright(N) = diag_(N) if(dright(N)==0d0)stop "matrix is ill-conditioned: no inverse exists" do i=N-1,1,-1 foo = 1d0/dright(i+1) cright = over_(i)*foo dright(i) = diag_(i) - cright*sub_(i) !sub_(i+1)/dright(i+1)*over_(i) if(dright(i)==0d0)stop "matrix is ill-conditioned: no inverse exists" enddo ! do i=1,N Foo = dleft(i) + dright(i) - diag_(i) Inv(i) = 1d0/Foo end do end subroutine c_invert_tridiag_matrix