CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
tridi.f
1
3 module tridi_mod
4 contains
5
10 subroutine tridi1(l,n,cl,cm,cu,r1,au,a1)
11 !
12 use machine , only : kind_phys
13 implicit none
14 integer, parameter :: one = 1.0_kind_phys
15 integer k,n,l,i
16 real(kind=kind_phys) fk
17 !
18 real(kind=kind_phys) cl(l,2:n),cm(l,n),cu(l,n-1),r1(l,n), &
19 & au(l,n-1),a1(l,n)
20 !
21 do i=1,l
22 fk = one / cm(i,1)
23 au(i,1) = fk*cu(i,1)
24 a1(i,1) = fk*r1(i,1)
25 enddo
26 do k=2,n-1
27 do i=1,l
28 fk = one / (cm(i,k)-cl(i,k)*au(i,k-1))
29 au(i,k) = fk*cu(i,k)
30 a1(i,k) = fk*(r1(i,k)-cl(i,k)*a1(i,k-1))
31 enddo
32 enddo
33 do i=1,l
34 fk = one / (cm(i,n)-cl(i,n)*au(i,n-1))
35 a1(i,n) = fk*(r1(i,n)-cl(i,n)*a1(i,n-1))
36 enddo
37 do k=n-1,1,-1
38 do i=1,l
39 a1(i,k) = a1(i,k)-au(i,k)*a1(i,k+1)
40 enddo
41 enddo
42 !
43 return
44 end subroutine tridi1
45
46!-----------------------------------------------------------------------
49 subroutine tridi2(l,n,cl,cm,cu,r1,r2,au,a1,a2)
50!
51 use machine , only : kind_phys
52 implicit none
53 integer, parameter :: one = 1.0_kind_phys
54 integer k,n,l,i
55 real(kind=kind_phys) fk
56!
57 real(kind=kind_phys) cl(l,2:n),cm(l,n),cu(l,n-1),r1(l,n),r2(l,n), &
58 & au(l,n-1),a1(l,n),a2(l,n)
59!----------------------------------------------------------------------
60 do i=1,l
61 fk = one / cm(i,1)
62 au(i,1) = fk*cu(i,1)
63 a1(i,1) = fk*r1(i,1)
64 a2(i,1) = fk*r2(i,1)
65 enddo
66 do k=2,n-1
67 do i=1,l
68 fk = one / (cm(i,k)-cl(i,k)*au(i,k-1))
69 au(i,k) = fk*cu(i,k)
70 a1(i,k) = fk*(r1(i,k)-cl(i,k)*a1(i,k-1))
71 a2(i,k) = fk*(r2(i,k)-cl(i,k)*a2(i,k-1))
72 enddo
73 enddo
74 do i=1,l
75 fk = one / (cm(i,n)-cl(i,n)*au(i,n-1))
76 a1(i,n) = fk*(r1(i,n)-cl(i,n)*a1(i,n-1))
77 a2(i,n) = fk*(r2(i,n)-cl(i,n)*a2(i,n-1))
78 enddo
79 do k=n-1,1,-1
80 do i=1,l
81 a1(i,k) = a1(i,k)-au(i,k)*a1(i,k+1)
82 a2(i,k) = a2(i,k)-au(i,k)*a2(i,k+1)
83 enddo
84 enddo
85!-----------------------------------------------------------------------
86 return
87 end subroutine tridi2
88
89!-----------------------------------------------------------------------
94 subroutine tridin(l,n,nt,cl,cm,cu,r1,r2,au,a1,a2)
95!
96 use machine , only : kind_phys
97 implicit none
98 integer, parameter :: one = 1.0_kind_phys
99 integer is,k,kk,n,nt,l,i
100 real(kind=kind_phys) fk(l)
101!
102 real(kind=kind_phys) cl(l,2:n), cm(l,n), cu(l,n-1), &
103 & r1(l,n), r2(l,n*nt), &
104 & au(l,n-1), a1(l,n), a2(l,n*nt), &
105 & fkk(l,2:n-1)
106!-----------------------------------------------------------------------
107 do i=1,l
108 fk(i) = one / cm(i,1)
109 au(i,1) = fk(i)*cu(i,1)
110 a1(i,1) = fk(i)*r1(i,1)
111 enddo
112 do k = 1, nt
113 is = (k-1) * n
114 do i = 1, l
115 a2(i,1+is) = fk(i) * r2(i,1+is)
116 enddo
117 enddo
118 do k=2,n-1
119 do i=1,l
120 fkk(i,k) = one / (cm(i,k)-cl(i,k)*au(i,k-1))
121 au(i,k) = fkk(i,k)*cu(i,k)
122 a1(i,k) = fkk(i,k)*(r1(i,k)-cl(i,k)*a1(i,k-1))
123 enddo
124 enddo
125 do kk = 1, nt
126 is = (kk-1) * n
127 do k=2,n-1
128 do i=1,l
129 a2(i,k+is) = fkk(i,k)*(r2(i,k+is)-cl(i,k)*a2(i,k+is-1))
130 enddo
131 enddo
132 enddo
133 do i=1,l
134 fk(i) = one / (cm(i,n)-cl(i,n)*au(i,n-1))
135 a1(i,n) = fk(i)*(r1(i,n)-cl(i,n)*a1(i,n-1))
136 enddo
137 do k = 1, nt
138 is = (k-1) * n
139 do i = 1, l
140 a2(i,n+is) = fk(i)*(r2(i,n+is)-cl(i,n)*a2(i,n+is-1))
141 enddo
142 enddo
143 do k=n-1,1,-1
144 do i=1,l
145 a1(i,k) = a1(i,k) - au(i,k)*a1(i,k+1)
146 enddo
147 enddo
148 do kk = 1, nt
149 is = (kk-1) * n
150 do k=n-1,1,-1
151 do i=1,l
152 a2(i,k+is) = a2(i,k+is) - au(i,k)*a2(i,k+is+1)
153 enddo
154 enddo
155 enddo
156!-----------------------------------------------------------------------
157 return
158 end subroutine tridin
159
160!-----------------------------------------------------------------------
163 subroutine tridit(l,n,nt,cl,cm,cu,rt,au,at)
164!-----------------------------------------------------------------------
165!!
166 use machine , only : kind_phys
167 implicit none
168 integer, parameter :: one = 1.0_kind_phys
169 integer is,k,kk,n,nt,l,i
170 real(kind=kind_phys) fk(l)
171!!
172 real(kind=kind_phys) cl(l,2:n), cm(l,n), cu(l,n-1), &
173 & rt(l,n*nt), &
174 & au(l,n-1), at(l,n*nt), &
175 & fkk(l,2:n-1)
176!-----------------------------------------------------------------------
177 do i=1,l
178 fk(i) = one / cm(i,1)
179 au(i,1) = fk(i)*cu(i,1)
180 enddo
181 do k = 1, nt
182 is = (k-1) * n
183 do i = 1, l
184 at(i,1+is) = fk(i) * rt(i,1+is)
185 enddo
186 enddo
187 do k=2,n-1
188 do i=1,l
189 fkk(i,k) = one / (cm(i,k)-cl(i,k)*au(i,k-1))
190 au(i,k) = fkk(i,k)*cu(i,k)
191 enddo
192 enddo
193 do kk = 1, nt
194 is = (kk-1) * n
195 do k=2,n-1
196 do i=1,l
197 at(i,k+is) = fkk(i,k)*(rt(i,k+is)-cl(i,k)*at(i,k+is-1))
198 enddo
199 enddo
200 enddo
201 do i=1,l
202 fk(i) = one / (cm(i,n)-cl(i,n)*au(i,n-1))
203 enddo
204 do k = 1, nt
205 is = (k-1) * n
206 do i = 1, l
207 at(i,n+is) = fk(i)*(rt(i,n+is)-cl(i,n)*at(i,n+is-1))
208 enddo
209 enddo
210 do kk = 1, nt
211 is = (kk-1) * n
212 do k=n-1,1,-1
213 do i=1,l
214 at(i,k+is) = at(i,k+is) - au(i,k)*at(i,k+is+1)
215 enddo
216 enddo
217 enddo
218!-----------------------------------------------------------------------
219 return
220 end subroutine tridit
221 end module tridi_mod
subroutine tridi1(l, n, cl, cm, cu, r1, au, a1)
Routine to solve the tridiagonal system to calculate temperature and moisture at .
Definition tridi.f:11
subroutine tridit(l, n, nt, cl, cm, cu, rt, au, at)
This subroutine solves tridiagonal problem for TKE.
Definition tridi.f:164
subroutine tridin(l, n, nt, cl, cm, cu, r1, r2, au, a1, a2)
Routine to solve the tridiagonal system to calculate u- and v-momentum at .
Definition tridi.f:95
subroutine tridi2(l, n, cl, cm, cu, r1, r2, au, a1, a2)
This subroutine ..
Definition tridi.f:50