BUGSrad
 All Classes Files Functions Variables
bandsolve.f
Go to the documentation of this file.
1 
2 
3 ! CVS: $Id: bandsolve.F,v 1.1 2001/04/30 08:43:42 norm Exp $
4 ! CVS: $Name: $
5 
6  module bandsolve
7 
8  implicit none
9 
10  public::bandec, banbks
11 
12  contains
13 
14  subroutine bandec(a, n, m1, m2, np, mp, al, mpl, indx, d)
15 
16  use kinds
17 
18  integer(kind=int_kind), intent(in) :: n, m1, m2, np, mp, mpl
19 
20  real(kind=dbl_kind), intent(inout), dimension(1:,1:) :: a
21 
22  integer(kind=int_kind), intent(out), dimension(1:) :: indx
23  real(kind=dbl_kind), intent(out), dimension(1:,1:) :: al
24  real(kind=dbl_kind), intent(out) :: d
25 
26  !Local variables
27  real(kind=dbl_kind), parameter:: tiny=1.0e-20_dbl_kind
28  integer(kind=int_kind) :: i, j, k, ell, mm
29  real(kind=dbl_kind) :: dum
30 
31  mm = m1 + m2 + 1
32 
33  ! if(mm > mp .or. m1 > mpl .or. n > np) then
34  ! print *, "Bad arguments in bandec"
35  ! !pause
36  ! endif
37 
38 
39  ell = m1
40 
41  do i = 1, m1
42  do j = m1+2-i, mm
43  a(i,j-ell) = a(i,j)
44  enddo
45  ell = ell - 1
46  do j = mm-ell, mm
47  a(i,j) = 0.0
48  enddo
49  enddo
50 
51  d = 1.0
52  ell = m1
53  do k = 1,n
54  dum = a(k,1)
55  i = k
56  if (ell < n) ell = ell+1
57  do j = k+1, ell
58  if (abs(a(j,1)) > abs(dum)) then
59  dum = a(j,1)
60  i = j
61  endif
62  enddo
63  indx(k) = i
64  if (dum == 0.0) a(k,1) = tiny !This was in num. rec.
65  if (i /= k) then
66  d = -d
67  do j = 1, mm
68  dum = a(k,j)
69  a(k,j) = a(i,j)
70  a(i,j) = dum
71  enddo
72  endif
73  do i = k+1, ell
74  dum = a(i,1)/a(k,1)
75  al(k,i-k) = dum
76  do j = 2, mm
77  a(i, j-1) = a(i,j) - dum*a(k,j)
78  enddo
79  a(i,mm) = 0.0
80  enddo
81  enddo
82  return
83  end subroutine bandec
84 
85 
86  subroutine banbks(a, n, m1, m2, np, mp, al, mpl, indx, b)
87 
88  use kinds
89 
90  integer(kind=int_kind), intent(in):: n, m1, m2, np, mp, mpl
91  integer(kind=int_kind), intent(in), dimension(:):: indx
92  real(kind=dbl_kind), intent(in), dimension(:,:):: a, al
93 
94  real(kind=dbl_kind), intent(inout), dimension(:):: b
95 
96  !Local variables
97 
98  integer(kind=int_kind):: i, k, ell, mm
99  real(kind=dbl_kind):: dum
100 
101  mm = m1 + m2 + 1
102 
103  if (mm > mp .or. m1 > mpl .or. n > np) then
104  print *, "Bad arguments in banbks"
105  endif
106 
107  ell = m1
108 
109  do k = 1, n
110  i = indx(k)
111  if (i /= k) then
112  dum = b(k)
113  b(k) = b(i)
114  b(i) = dum
115  endif
116  if (ell < n) ell = ell+1
117  do i = k+1, ell
118  b(i) = b(i) - al(k,i-k)*b(k)
119  enddo
120  enddo
121 
122  ell = 1
123  do i = n, 1, -1
124  dum = b(i)
125  do k = 2, ell
126  dum = dum - a(i,k)*b(k+i-1)
127  enddo
128  b(i) = dum/a(i,1)
129  if (ell < mm) ell = ell+1
130  enddo
131  return
132  end subroutine banbks
133 
134  end module bandsolve
Generated on Tue Aug 13 2013 10:38:25 for BUGSrad by   doxygen 1.8.3.1