Skip to content

Instantly share code, notes, and snippets.

@shane5ul
Created March 5, 2014 19:59
Show Gist options
  • Save shane5ul/9375343 to your computer and use it in GitHub Desktop.
Save shane5ul/9375343 to your computer and use it in GitHub Desktop.
Prettified f90ppr example
!
!--------------------------------------------------------------------------------------------------
subroutine tridiag (a, b, c, r, u, n, code)
! Solves for a vector u of length n the tridiagonal linear set ! from numerical recipes
! m u = r, where a, b and c are the three main diagonals of matrix !
! m(n,n), the other terms are 0. r is the right side vector. !
!--------------------------------------------------------------------------------------------------
integer, parameter :: nmax = 100
integer, intent (in) :: n
double precision, intent (in) :: a (n), b (n), c (n), r (n)
double precision, intent (out) :: u (n)
integer, intent (out) :: code
!
double precision :: bet, gam (nmax)
integer :: j
!
if (b(1) .eq. 0.d0) then
code = 1 ! error unless code=0
return
end if
!
bet = b (1)
u (1) = r (1) / bet
do j = 2, n ! decomposition and forward substitution
gam (j) = c (j-1) / bet
bet = b (j) - a (j) * gam (j)
if (bet .eq. 0.d0) then ! algorithm fails
code = 2
return
end if
u (j) = (r(j)-a(j)*u(j-1)) / bet
end do
!
do j = n - 1, 1, - 1 !back substitution
u (j) = u (j) - gam (j+1) * u (j+1)
end do
!
code = 0
return
!
end subroutine tridiag
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment