module_small_step_em_ad.F
References to this file elsewhere.
1 ! DISCLAIMER
2 !
3 ! This file was generated by TAF version 1.7.18
4 !
5 ! FASTOPT DISCLAIMS ALL WARRANTIES, EXPRESS OR IMPLIED,
6 ! INCLUDING (WITHOUT LIMITATION) ALL IMPLIED WARRANTIES OF
7 ! MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, WITH
8 ! RESPECT TO THE SOFTWARE AND USER PROGRAMS. IN NO EVENT
9 ! SHALL FASTOPT BE LIABLE FOR ANY LOST OR ANTICIPATED PROF-
10 ! ITS, OR ANY INDIRECT, INCIDENTAL, EXEMPLARY, SPECIAL, OR
11 ! CONSEQUENTIAL DAMAGES, WHETHER OR NOT FASTOPT WAS ADVISED
12 ! OF THE POSSIBILITY OF SUCH DAMAGES.
13 !
14 ! Haftungsbeschraenkung
15 ! FastOpt gibt ausdruecklich keine Gewaehr, explizit oder indirekt,
16 ! bezueglich der Brauchbarkeit der Software fuer einen bestimmten
17 ! Zweck. Unter keinen Umstaenden ist FastOpt haftbar fuer
18 ! irgendeinen Verlust oder nicht eintretenden erwarteten Gewinn und
19 ! allen indirekten, zufaelligen, exemplarischen oder speziellen
20 ! Schaeden oder Folgeschaeden unabhaengig von einer eventuellen
21 ! Mitteilung darueber an FastOpt.
22 !
23 module a_module_small_step_em
24 !******************************************************************
25 !******************************************************************
26 !** This routine was generated by Automatic differentiation. **
27 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
28 !******************************************************************
29 !******************************************************************
30 !==============================================
31 ! referencing used modules
32 !==============================================
33 use module_configure
34 use module_model_constants
35 use module_small_step_em
36
37 !==============================================
38 ! all entries are defined explicitly
39 !==============================================
40 implicit none
41
42 contains
43 subroutine a_advance_mu_t( ww, a_ww, ww_1, a_ww_1, u, a_u, u_1, a_u_1, v, a_v, v_1, a_v_1, a_mu, a_mut, a_muave, a_muts, muu, &
44 &a_muu, muv, a_muv, a_mudf, a_t, t_1, a_t_1, a_t_ave, a_ft, mu_tend, a_mu_tend, rdx, rdy, dts, epssm, dnw, fnm, fnp, rdnw, msfu, &
45 &msfv, msft, config_flags, ids, ide, jds, jde, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte )
46 !******************************************************************
47 !******************************************************************
48 !** This routine was generated by Automatic differentiation. **
49 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
50 !******************************************************************
51 !******************************************************************
52 !==============================================
53 ! all entries are defined explicitly
54 !==============================================
55 implicit none
56
57 !==============================================
58 ! declare arguments
59 !==============================================
60 integer, intent(in) :: ime
61 integer, intent(in) :: ims
62 integer, intent(in) :: jme
63 integer, intent(in) :: jms
64 integer, intent(in) :: kme
65 integer, intent(in) :: kms
66 real, intent(inout) :: a_ft(ims:ime,kms:kme,jms:jme)
67 real, intent(inout) :: a_mu(ims:ime,jms:jme)
68 real, intent(inout) :: a_mu_tend(ims:ime,jms:jme)
69 real, intent(inout) :: a_muave(ims:ime,jms:jme)
70 real, intent(inout) :: a_mudf(ims:ime,jms:jme)
71 real, intent(inout) :: a_mut(ims:ime,jms:jme)
72 real, intent(inout) :: a_muts(ims:ime,jms:jme)
73 real, intent(inout) :: a_muu(ims:ime,jms:jme)
74 real, intent(inout) :: a_muv(ims:ime,jms:jme)
75 real, intent(inout) :: a_t(ims:ime,kms:kme,jms:jme)
76 real, intent(inout) :: a_t_1(ims:ime,kms:kme,jms:jme)
77 real, intent(inout) :: a_t_ave(ims:ime,kms:kme,jms:jme)
78 real, intent(inout) :: a_u(ims:ime,kms:kme,jms:jme)
79 real, intent(inout) :: a_u_1(ims:ime,kms:kme,jms:jme)
80 real, intent(inout) :: a_v(ims:ime,kms:kme,jms:jme)
81 real, intent(inout) :: a_v_1(ims:ime,kms:kme,jms:jme)
82 real, intent(inout) :: a_ww(ims:ime,kms:kme,jms:jme)
83 real, intent(inout) :: a_ww_1(ims:ime,kms:kme,jms:jme)
84 type (grid_config_rec_type), intent(in) :: config_flags
85 real, intent(in) :: dnw(kms:kme)
86 real, intent(in) :: dts
87 real, intent(in) :: epssm
88 real, intent(in) :: fnm(kms:kme)
89 real, intent(in) :: fnp(kms:kme)
90 integer, intent(in) :: ide
91 integer, intent(in) :: ids
92 integer, intent(in) :: ite
93 integer, intent(in) :: its
94 integer, intent(in) :: jde
95 integer, intent(in) :: jds
96 integer, intent(in) :: jte
97 integer, intent(in) :: jts
98 integer, intent(in) :: kde
99 integer, intent(in) :: kte
100 integer, intent(in) :: kts
101 real, intent(in) :: msft(ims:ime,jms:jme)
102 real, intent(in) :: msfu(ims:ime,jms:jme)
103 real, intent(in) :: msfv(ims:ime,jms:jme)
104 real, intent(in) :: mu_tend(ims:ime,jms:jme)
105 real, intent(in) :: muu(ims:ime,jms:jme)
106 real, intent(in) :: muv(ims:ime,jms:jme)
107 real, intent(in) :: rdnw(kms:kme)
108 real, intent(in) :: rdx
109 real, intent(in) :: rdy
110 real, intent(in) :: t_1(ims:ime,kms:kme,jms:jme)
111 real, intent(in) :: u(ims:ime,kms:kme,jms:jme)
112 real, intent(in) :: u_1(ims:ime,kms:kme,jms:jme)
113 real, intent(in) :: v(ims:ime,kms:kme,jms:jme)
114 real, intent(in) :: v_1(ims:ime,kms:kme,jms:jme)
115 real, intent(inout) :: ww(ims:ime,kms:kme,jms:jme)
116 real, intent(inout) :: ww_1(ims:ime,kms:kme,jms:jme)
117
118 !==============================================
119 ! declare local variables
120 !==============================================
121 real a_dmdt(its:ite)
122 real a_dvdxi(its:ite,kts:kte)
123 real a_wdtn(its:ite,kts:kte)
124 real dmdt(its:ite)
125 real dvdxi(its:ite,kts:kte)
126 integer i
127 integer i_end
128 integer i_start
129 integer j
130 integer j_end
131 integer j_start
132 integer k
133 integer k_end
134 integer k_start
135
136 !----------------------------------------------
137 ! RESET LOCAL ADJOINT VARIABLES
138 !----------------------------------------------
139 a_dmdt(:) = 0.
140 a_dvdxi(:,:) = 0.
141 a_wdtn(:,:) = 0.
142
143 !----------------------------------------------
144 ! ROUTINE BODY
145 !----------------------------------------------
146 i_start = its
147 ! recompute : i_start
148 i_end = ite
149 ! recompute : i_end
150 j_start = jts
151 ! recompute : j_start
152 j_end = jte
153 ! recompute : j_end
154 k_start = kts
155 ! recompute : k_start
156 k_end = kte-1
157 ! recompute : k_end
158 if (j_end .eq. jde) then
159 j_end = j_end-1
160 endif
161 ! recompute : j_end
162 if (i_end .eq. ide) then
163 i_end = i_end-1
164 endif
165 ! recompute : i_end
166 if ((config_flags%specified .or. config_flags%nested) .and. its .eq. ids) then
167 i_start = i_start+1
168 endif
169 ! recompute : i_start
170 if ((config_flags%specified .or. config_flags%nested) .and. ite .eq. ide) then
171 i_end = i_end-1
172 endif
173 ! recompute : i_end
174 if ((config_flags%specified .or. config_flags%nested) .and. jts .eq. jds) then
175 j_start = j_start+1
176 endif
177 ! recompute : j_start
178 if ((config_flags%specified .or. config_flags%nested) .and. jte .eq. jde) then
179 j_end = j_end-1
180 endif
181 ! recompute : j_end
182 do j = j_start, j_end
183 do i = i_start, i_end
184 dmdt(i) = 0.
185 end do
186 do k = k_start, k_end
187 do i = i_start, i_end
188 dvdxi(i,k) = msft(i,j)*msft(i,j)*(rdy*(v(i,k,j+1)+muv(i,j+1)*v_1(i,k,j+1)/msfv(i,j+1)-(v(i,k,j)+muv(i,j)*v_1(i,k,j)/msfv(i,j)&
189 &))+rdx*(u(i+1,k,j)+muu(i+1,j)*u_1(i+1,k,j)/msfu(i+1,j)-(u(i,k,j)+muu(i,j)*u_1(i,k,j)/msfu(i,j))))
190 dmdt(i) = dmdt(i)+dnw(k)*dvdxi(i,k)
191 end do
192 end do
193 do k = 2, k_end
194 do i = i_start, i_end
195 ww(i,k,j) = ww(i,k-1,j)-dnw(k-1)*(dmdt(i)+dvdxi(i,k-1)+mu_tend(i,j))/msft(i,j)
196 end do
197 end do
198 do k = 1, k_end
199 do i = i_start, i_end
200 ww(i,k,j) = ww(i,k,j)-ww_1(i,k,j)
201 end do
202 end do
203 end do
204 ! recompute : ww
205 do j = j_start, j_end
206 do k = 1, k_end
207 do i = i_start, i_end
208 a_t_1(i,k,j-1) = a_t_1(i,k,j-1)+0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdy*v(i,k,j)
209 a_t_1(i,k,j+1) = a_t_1(i,k,j+1)-0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdy*v(i,k,j+1)
210 a_t_1(i-1,k,j) = a_t_1(i-1,k,j)+0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdx*u(i,k,j)
211 a_t_1(i+1,k,j) = a_t_1(i+1,k,j)-0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdx*u(i+1,k,j)
212 a_t_1(i,k,j) = a_t_1(i,k,j)-a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*(0.5*rdy*(v(i,k,j+1)-v(i,k,j))+0.5*rdx*(u(i+1,k,j)-u(i,k,j)))
213 a_u(i+1,k,j) = a_u(i+1,k,j)-0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdx*(t_1(i+1,k,j)+t_1(i,k,j))
214 a_u(i,k,j) = a_u(i,k,j)+0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdx*(t_1(i,k,j)+t_1(i-1,k,j))
215 a_v(i,k,j+1) = a_v(i,k,j+1)-0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdy*(t_1(i,k,j+1)+t_1(i,k,j))
216 a_v(i,k,j) = a_v(i,k,j)+0.5*a_t(i,k,j)*dts*msft(i,j)*msft(i,j)*rdy*(t_1(i,k,j)+t_1(i,k,j-1))
217 a_wdtn(i,k+1) = a_wdtn(i,k+1)-a_t(i,k,j)*dts*msft(i,j)*rdnw(k)
218 a_wdtn(i,k) = a_wdtn(i,k)+a_t(i,k,j)*dts*msft(i,j)*rdnw(k)
219 end do
220 end do
221 do k = 2, k_end
222 do i = i_start, i_end
223 a_t_1(i,k-1,j) = a_t_1(i,k-1,j)+a_wdtn(i,k)*ww(i,k,j)*fnp(k)
224 a_t_1(i,k,j) = a_t_1(i,k,j)+a_wdtn(i,k)*ww(i,k,j)*fnm(k)
225 a_ww(i,k,j) = a_ww(i,k,j)+a_wdtn(i,k)*(fnm(k)*t_1(i,k,j)+fnp(k)*t_1(i,k-1,j))
226 a_wdtn(i,k) = 0.
227 end do
228 end do
229 do i = i_start, i_end
230 a_wdtn(i,kde) = 0.
231 a_wdtn(i,1) = 0.
232 end do
233 end do
234 do j = j_start, j_end
235 do k = 1, k_end
236 do i = i_start, i_end
237 a_ft(i,k,j) = a_ft(i,k,j)+a_t(i,k,j)*msft(i,j)*dts
238 a_t(i,k,j) = a_t(i,k,j)+a_t_ave(i,k,j)
239 a_t_ave(i,k,j) = 0.
240 end do
241 end do
242 end do
243 do j = j_start, j_end
244 do k = 1, k_end
245 do i = i_start, i_end
246 a_ww_1(i,k,j) = a_ww_1(i,k,j)-a_ww(i,k,j)
247 end do
248 end do
249 do k = k_end, 2, -1
250 do i = i_start, i_end
251 a_dmdt(i) = a_dmdt(i)-a_ww(i,k,j)*(dnw(k-1)/msft(i,j))
252 a_dvdxi(i,k-1) = a_dvdxi(i,k-1)-a_ww(i,k,j)*(dnw(k-1)/msft(i,j))
253 a_mu_tend(i,j) = a_mu_tend(i,j)-a_ww(i,k,j)*(dnw(k-1)/msft(i,j))
254 a_ww(i,k-1,j) = a_ww(i,k-1,j)+a_ww(i,k,j)
255 a_ww(i,k,j) = 0.
256 end do
257 end do
258 do i = i_start, i_end
259 a_mu(i,j) = a_mu(i,j)+0.5*a_muave(i,j)*(1+epssm)
260 a_muave(i,j) = 0.5*a_muave(i,j)*(1.-epssm)
261 a_mu(i,j) = a_mu(i,j)+a_muts(i,j)
262 a_mut(i,j) = a_mut(i,j)+a_muts(i,j)
263 a_muts(i,j) = 0.
264 a_dmdt(i) = a_dmdt(i)+a_mudf(i,j)
265 a_mu_tend(i,j) = a_mu_tend(i,j)+a_mudf(i,j)
266 a_mudf(i,j) = 0.
267 a_dmdt(i) = a_dmdt(i)+a_mu(i,j)*dts
268 a_mu_tend(i,j) = a_mu_tend(i,j)+a_mu(i,j)*dts
269 a_mu(i,j) = a_mu(i,j)+a_muave(i,j)
270 a_muave(i,j) = 0.
271 end do
272 do k = k_end, k_start, -1
273 do i = i_start, i_end
274 a_dvdxi(i,k) = a_dvdxi(i,k)+a_dmdt(i)*dnw(k)
275 a_muu(i+1,j) = a_muu(i+1,j)+a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdx*(u_1(i+1,k,j)/msfu(i+1,j))
276 a_muu(i,j) = a_muu(i,j)-a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdx*(u_1(i,k,j)/msfu(i,j))
277 a_muv(i,j+1) = a_muv(i,j+1)+a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdy*(v_1(i,k,j+1)/msfv(i,j+1))
278 a_muv(i,j) = a_muv(i,j)-a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdy*(v_1(i,k,j)/msfv(i,j))
279 a_u(i+1,k,j) = a_u(i+1,k,j)+a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdx
280 a_u(i,k,j) = a_u(i,k,j)-a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdx
281 a_u_1(i+1,k,j) = a_u_1(i+1,k,j)+a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdx*(muu(i+1,j)/msfu(i+1,j))
282 a_u_1(i,k,j) = a_u_1(i,k,j)-a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdx*(muu(i,j)/msfu(i,j))
283 a_v(i,k,j+1) = a_v(i,k,j+1)+a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdy
284 a_v(i,k,j) = a_v(i,k,j)-a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdy
285 a_v_1(i,k,j+1) = a_v_1(i,k,j+1)+a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdy*(muv(i,j+1)/msfv(i,j+1))
286 a_v_1(i,k,j) = a_v_1(i,k,j)-a_dvdxi(i,k)*msft(i,j)*msft(i,j)*rdy*(muv(i,j)/msfv(i,j))
287 a_dvdxi(i,k) = 0.
288 end do
289 end do
290 do i = i_start, i_end
291 a_dmdt(i) = 0.
292 end do
293 end do
294
295 end subroutine a_advance_mu_t
296
297
298 subroutine a_advance_uv( a_u, a_ru_tend, a_v, a_rv_tend, p, a_p, pb, ph, a_ph, php, a_php, alt, a_alt, al, a_al, mu, a_mu, muu, &
299 &a_muu, cqu, a_cqu, muv, a_muv, cqv, a_cqv, a_mudf, rdx, rdy, dts, cf1, cf2, cf3, fnm, fnp, emdiv, rdnw, config_flags, spec_zone, &
300 &non_hydrostatic, ids, ide, jds, jde, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte )
301 !******************************************************************
302 !******************************************************************
303 !** This routine was generated by Automatic differentiation. **
304 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
305 !******************************************************************
306 !******************************************************************
307 !==============================================
308 ! all entries are defined explicitly
309 !==============================================
310 implicit none
311
312 !==============================================
313 ! declare arguments
314 !==============================================
315 integer, intent(in) :: ime
316 integer, intent(in) :: ims
317 integer, intent(in) :: jme
318 integer, intent(in) :: jms
319 integer, intent(in) :: kme
320 integer, intent(in) :: kms
321 real, intent(inout) :: a_al(ims:ime,kms:kme,jms:jme)
322 real, intent(inout) :: a_alt(ims:ime,kms:kme,jms:jme)
323 real, intent(inout) :: a_cqu(ims:ime,kms:kme,jms:jme)
324 real, intent(inout) :: a_cqv(ims:ime,kms:kme,jms:jme)
325 real, intent(inout) :: a_mu(ims:ime,jms:jme)
326 real, intent(inout) :: a_mudf(ims:ime,jms:jme)
327 real, intent(inout) :: a_muu(ims:ime,jms:jme)
328 real, intent(inout) :: a_muv(ims:ime,jms:jme)
329 real, intent(inout) :: a_p(ims:ime,kms:kme,jms:jme)
330 real, intent(inout) :: a_ph(ims:ime,kms:kme,jms:jme)
331 real, intent(inout) :: a_php(ims:ime,kms:kme,jms:jme)
332 real, intent(inout) :: a_ru_tend(ims:ime,kms:kme,jms:jme)
333 real, intent(inout) :: a_rv_tend(ims:ime,kms:kme,jms:jme)
334 real, intent(inout) :: a_u(ims:ime,kms:kme,jms:jme)
335 real, intent(inout) :: a_v(ims:ime,kms:kme,jms:jme)
336 real, intent(in) :: al(ims:ime,kms:kme,jms:jme)
337 real, intent(in) :: alt(ims:ime,kms:kme,jms:jme)
338 real, intent(in) :: cf1
339 real, intent(in) :: cf2
340 real, intent(in) :: cf3
341 type (grid_config_rec_type), intent(in) :: config_flags
342 real, intent(in) :: cqu(ims:ime,kms:kme,jms:jme)
343 real, intent(in) :: cqv(ims:ime,kms:kme,jms:jme)
344 real, intent(in) :: dts
345 real, intent(in) :: emdiv
346 real, intent(in) :: fnm(kms:kme)
347 real, intent(in) :: fnp(kms:kme)
348 integer, intent(in) :: ide
349 integer, intent(in) :: ids
350 integer, intent(in) :: ite
351 integer, intent(in) :: its
352 integer, intent(in) :: jde
353 integer, intent(in) :: jds
354 integer, intent(in) :: jte
355 integer, intent(in) :: jts
356 integer, intent(in) :: kde
357 integer, intent(in) :: kte
358 integer, intent(in) :: kts
359 real, intent(in) :: mu(ims:ime,jms:jme)
360 real, intent(in) :: muu(ims:ime,jms:jme)
361 real, intent(in) :: muv(ims:ime,jms:jme)
362 logical, intent(in) :: non_hydrostatic
363 real, intent(in) :: p(ims:ime,kms:kme,jms:jme)
364 real, intent(in) :: pb(ims:ime,kms:kme,jms:jme)
365 real, intent(in) :: ph(ims:ime,kms:kme,jms:jme)
366 real, intent(in) :: php(ims:ime,kms:kme,jms:jme)
367 real, intent(in) :: rdnw(kms:kme)
368 real, intent(in) :: rdx
369 real, intent(in) :: rdy
370 integer, intent(in) :: spec_zone
371
372 !==============================================
373 ! declare local variables
374 !==============================================
375 real a_dpn(its:ite,kts:kte)
376 real a_dpxy(its:ite,kts:kte)
377 real a_mudf_xy(its:ite)
378 real dpn(its:ite,kts:kte)
379 real dpxy(its:ite,kts:kte)
380 real dx
381 real dy
382 integer i
383 integer i_end
384 integer i_end_up
385 integer i_endu
386 integer i_start
387 integer i_start_up
388 integer j
389 integer j_end
390 integer j_end_vp
391 integer j_endv
392 integer j_start
393 integer j_start_vp
394 integer k
395 integer k_end
396 integer k_start
397
398 !----------------------------------------------
399 ! RESET LOCAL ADJOINT VARIABLES
400 !----------------------------------------------
401 a_dpn(:,:) = 0.
402 a_dpxy(:,:) = 0.
403 a_mudf_xy(:) = 0.
404
405 !----------------------------------------------
406 ! ROUTINE BODY
407 !----------------------------------------------
408 if (config_flags%nested .or. config_flags%specified) then
409 i_start = max(its,ids+spec_zone)
410 i_end = min(ite,ide-spec_zone-1)
411 j_start = max(jts,jds+spec_zone)
412 j_end = min(jte,jde-spec_zone-1)
413 k_start = kts
414 k_end = min(kte,kde-1)
415 i_endu = min(ite,ide-spec_zone)
416 j_endv = min(jte,jde-spec_zone)
417 else
418 i_start = its
419 i_end = ite
420 j_start = jts
421 j_end = jte
422 k_start = kts
423 k_end = kte-1
424 i_endu = i_end
425 j_endv = j_end
426 if (j_end .eq. jde) then
427 j_end = j_end-1
428 endif
429 if (i_end .eq. ide) then
430 i_end = i_end-1
431 endif
432 endif
433 ! recompute : i_end,i_endu,i_start,j_end,j_endv,j_start,k_end,k_start
434 i_start_up = i_start
435 ! recompute : i_start_up
436 i_end_up = i_endu
437 ! recompute : i_end_up
438 j_start_vp = j_start
439 ! recompute : j_start_vp
440 j_end_vp = j_endv
441 ! recompute : j_end_vp
442 if ((config_flags%open_xs .or. config_flags%symmetric_xs) .and. its .eq. ids) then
443 i_start_up = i_start_up+1
444 endif
445 ! recompute : i_start_up
446 if ((config_flags%open_xe .or. config_flags%symmetric_xe) .and. ite .eq. ide) then
447 i_end_up = i_end_up-1
448 endif
449 ! recompute : i_end_up
450 if ((config_flags%open_ys .or. config_flags%symmetric_ys) .and. jts .eq. jds) then
451 j_start_vp = j_start_vp+1
452 endif
453 ! recompute : j_start_vp
454 if ((config_flags%open_ye .or. config_flags%symmetric_ye) .and. jte .eq. jde) then
455 j_end_vp = j_end_vp-1
456 endif
457 ! recompute : j_end_vp
458 dx = 1./rdx
459 ! recompute : dx
460 dy = 1./rdy
461 ! recompute : dy
462 a_v_outer_j_loop: do j = j_endv, j_start, -1
463 if (j .ge. j_start_vp .and. j .le. j_end_vp) then
464 do k = k_start, k_end
465 do i = i_start, i_end
466 dpxy(i,k) = 0.5*rdy*muv(i,j)*(ph(i,k+1,j)-ph(i,k+1,j-1)+ph(i,k,j)-ph(i,k,j-1)+(alt(i,k,j)+alt(i,k,j-1))*(p(i,k,j)-p(i,k,j-&
467 &1))+(al(i,k,j)+al(i,k,j-1))*(pb(i,k,j)-pb(i,k,j-1)))
468 end do
469 end do
470 ! recompute : dpxy
471 if (non_hydrostatic) then
472 do i = i_start, i_end
473 dpn(i,1) = 0.5*(cf1*(p(i,1,j)+p(i,1,j-1))+cf2*(p(i,2,j)+p(i,2,j-1))+cf3*(p(i,3,j)+p(i,3,j-1)))
474 dpn(i,kde) = 0.
475 end do
476 do k = k_start+1, k_end
477 do i = i_start, i_end
478 dpn(i,k) = 0.5*(fnm(k)*(p(i,k,j)+p(i,k,j-1))+fnp(k)*(p(i,k-1,j)+p(i,k-1,j-1)))
479 end do
480 end do
481 do k = k_start, k_end
482 do i = i_start, i_end
483 dpxy(i,k) = dpxy(i,k)+rdy*(php(i,k,j)-php(i,k,j-1))*(rdnw(k)*(dpn(i,k+1)-dpn(i,k))-0.5*(mu(i,j-1)+mu(i,j)))
484 end do
485 end do
486 endif
487 ! recompute : dpxy
488 do k = k_start, k_end
489 do i = i_start, i_end
490 a_cqv(i,k,j) = a_cqv(i,k,j)-a_v(i,k,j)*dts*dpxy(i,k)
491 a_dpxy(i,k) = a_dpxy(i,k)-a_v(i,k,j)*dts*cqv(i,k,j)
492 a_mudf_xy(i) = a_mudf_xy(i)+a_v(i,k,j)
493 end do
494 end do
495 if (non_hydrostatic) then
496 ! recompute : dpn
497 do k = k_start, k_end
498 do i = i_start, i_end
499 a_dpn(i,k+1) = a_dpn(i,k+1)+a_dpxy(i,k)*rdy*(php(i,k,j)-php(i,k,j-1))*rdnw(k)
500 a_dpn(i,k) = a_dpn(i,k)-a_dpxy(i,k)*rdy*(php(i,k,j)-php(i,k,j-1))*rdnw(k)
501 a_mu(i,j-1) = a_mu(i,j-1)-0.5*a_dpxy(i,k)*rdy*(php(i,k,j)-php(i,k,j-1))
502 a_mu(i,j) = a_mu(i,j)-0.5*a_dpxy(i,k)*rdy*(php(i,k,j)-php(i,k,j-1))
503 a_php(i,k,j-1) = a_php(i,k,j-1)-a_dpxy(i,k)*rdy*(rdnw(k)*(dpn(i,k+1)-dpn(i,k))-0.5*(mu(i,j-1)+mu(i,j)))
504 a_php(i,k,j) = a_php(i,k,j)+a_dpxy(i,k)*rdy*(rdnw(k)*(dpn(i,k+1)-dpn(i,k))-0.5*(mu(i,j-1)+mu(i,j)))
505 end do
506 end do
507 do k = k_start+1, k_end
508 do i = i_start, i_end
509 a_p(i,k-1,j-1) = a_p(i,k-1,j-1)+0.5*a_dpn(i,k)*fnp(k)
510 a_p(i,k-1,j) = a_p(i,k-1,j)+0.5*a_dpn(i,k)*fnp(k)
511 a_p(i,k,j-1) = a_p(i,k,j-1)+0.5*a_dpn(i,k)*fnm(k)
512 a_p(i,k,j) = a_p(i,k,j)+0.5*a_dpn(i,k)*fnm(k)
513 a_dpn(i,k) = 0.
514 end do
515 end do
516 do i = i_start, i_end
517 a_dpn(i,kde) = 0.
518 a_p(i,3,j-1) = a_p(i,3,j-1)+0.5*a_dpn(i,1)*cf3
519 a_p(i,3,j) = a_p(i,3,j)+0.5*a_dpn(i,1)*cf3
520 a_p(i,2,j-1) = a_p(i,2,j-1)+0.5*a_dpn(i,1)*cf2
521 a_p(i,2,j) = a_p(i,2,j)+0.5*a_dpn(i,1)*cf2
522 a_p(i,1,j-1) = a_p(i,1,j-1)+0.5*a_dpn(i,1)*cf1
523 a_p(i,1,j) = a_p(i,1,j)+0.5*a_dpn(i,1)*cf1
524 a_dpn(i,1) = 0.
525 end do
526 endif
527 do k = k_start, k_end
528 do i = i_start, i_end
529 a_al(i,k,j-1) = a_al(i,k,j-1)+0.5*a_dpxy(i,k)*rdy*muv(i,j)*(pb(i,k,j)-pb(i,k,j-1))
530 a_al(i,k,j) = a_al(i,k,j)+0.5*a_dpxy(i,k)*rdy*muv(i,j)*(pb(i,k,j)-pb(i,k,j-1))
531 a_alt(i,k,j-1) = a_alt(i,k,j-1)+0.5*a_dpxy(i,k)*rdy*muv(i,j)*(p(i,k,j)-p(i,k,j-1))
532 a_alt(i,k,j) = a_alt(i,k,j)+0.5*a_dpxy(i,k)*rdy*muv(i,j)*(p(i,k,j)-p(i,k,j-1))
533 a_muv(i,j) = a_muv(i,j)+0.5*a_dpxy(i,k)*rdy*(ph(i,k+1,j)-ph(i,k+1,j-1)+ph(i,k,j)-ph(i,k,j-1)+(alt(i,k,j)+alt(i,k,j-1))*&
534 &(p(i,k,j)-p(i,k,j-1))+(al(i,k,j)+al(i,k,j-1))*(pb(i,k,j)-pb(i,k,j-1)))
535 a_p(i,k,j-1) = a_p(i,k,j-1)-0.5*a_dpxy(i,k)*rdy*muv(i,j)*(alt(i,k,j)+alt(i,k,j-1))
536 a_p(i,k,j) = a_p(i,k,j)+0.5*a_dpxy(i,k)*rdy*muv(i,j)*(alt(i,k,j)+alt(i,k,j-1))
537 a_ph(i,k+1,j-1) = a_ph(i,k+1,j-1)-0.5*a_dpxy(i,k)*rdy*muv(i,j)
538 a_ph(i,k+1,j) = a_ph(i,k+1,j)+0.5*a_dpxy(i,k)*rdy*muv(i,j)
539 a_ph(i,k,j-1) = a_ph(i,k,j-1)-0.5*a_dpxy(i,k)*rdy*muv(i,j)
540 a_ph(i,k,j) = a_ph(i,k,j)+0.5*a_dpxy(i,k)*rdy*muv(i,j)
541 a_dpxy(i,k) = 0.
542 end do
543 end do
544 endif
545 do i = i_start, i_end
546 a_mudf(i,j-1) = a_mudf(i,j-1)+a_mudf_xy(i)*emdiv*dy
547 a_mudf(i,j) = a_mudf(i,j)-a_mudf_xy(i)*emdiv*dy
548 a_mudf_xy(i) = 0.
549 end do
550 do k = k_start, k_end
551 do i = i_start, i_end
552 a_rv_tend(i,k,j) = a_rv_tend(i,k,j)+a_v(i,k,j)*dts
553 end do
554 end do
555 end do a_v_outer_j_loop
556 do j = j_start, j_end
557 do k = k_start, k_end
558 do i = i_start_up, i_end_up
559 dpxy(i,k) = 0.5*rdx*muu(i,j)*(ph(i,k+1,j)-ph(i-1,k+1,j)+ph(i,k,j)-ph(i-1,k,j)+(alt(i,k,j)+alt(i-1,k,j))*(p(i,k,j)-p(i-1,k,j))&
560 &+(al(i,k,j)+al(i-1,k,j))*(pb(i,k,j)-pb(i-1,k,j)))
561 end do
562 end do
563 ! recompute : dpxy
564 if (non_hydrostatic) then
565 do i = i_start_up, i_end_up
566 dpn(i,1) = 0.5*(cf1*(p(i,1,j)+p(i-1,1,j))+cf2*(p(i,2,j)+p(i-1,2,j))+cf3*(p(i,3,j)+p(i-1,3,j)))
567 dpn(i,kde) = 0.
568 end do
569 do k = k_start+1, k_end
570 do i = i_start_up, i_end_up
571 dpn(i,k) = 0.5*(fnm(k)*(p(i,k,j)+p(i-1,k,j))+fnp(k)*(p(i,k-1,j)+p(i-1,k-1,j)))
572 end do
573 end do
574 do k = k_start, k_end
575 do i = i_start_up, i_end_up
576 dpxy(i,k) = dpxy(i,k)+rdx*(php(i,k,j)-php(i-1,k,j))*(rdnw(k)*(dpn(i,k+1)-dpn(i,k))-0.5*(mu(i-1,j)+mu(i,j)))
577 end do
578 end do
579 endif
580 ! recompute : dpxy
581 do k = k_start, k_end
582 do i = i_start_up, i_end_up
583 a_cqu(i,k,j) = a_cqu(i,k,j)-a_u(i,k,j)*dts*dpxy(i,k)
584 a_dpxy(i,k) = a_dpxy(i,k)-a_u(i,k,j)*dts*cqu(i,k,j)
585 a_mudf_xy(i) = a_mudf_xy(i)+a_u(i,k,j)
586 end do
587 end do
588 if (non_hydrostatic) then
589 ! recompute : dpn
590 do k = k_start, k_end
591 do i = i_start_up, i_end_up
592 a_dpn(i,k+1) = a_dpn(i,k+1)+a_dpxy(i,k)*rdx*(php(i,k,j)-php(i-1,k,j))*rdnw(k)
593 a_dpn(i,k) = a_dpn(i,k)-a_dpxy(i,k)*rdx*(php(i,k,j)-php(i-1,k,j))*rdnw(k)
594 a_mu(i-1,j) = a_mu(i-1,j)-0.5*a_dpxy(i,k)*rdx*(php(i,k,j)-php(i-1,k,j))
595 a_mu(i,j) = a_mu(i,j)-0.5*a_dpxy(i,k)*rdx*(php(i,k,j)-php(i-1,k,j))
596 a_php(i-1,k,j) = a_php(i-1,k,j)-a_dpxy(i,k)*rdx*(rdnw(k)*(dpn(i,k+1)-dpn(i,k))-0.5*(mu(i-1,j)+mu(i,j)))
597 a_php(i,k,j) = a_php(i,k,j)+a_dpxy(i,k)*rdx*(rdnw(k)*(dpn(i,k+1)-dpn(i,k))-0.5*(mu(i-1,j)+mu(i,j)))
598 end do
599 end do
600 do k = k_start+1, k_end
601 do i = i_start_up, i_end_up
602 a_p(i-1,k-1,j) = a_p(i-1,k-1,j)+0.5*a_dpn(i,k)*fnp(k)
603 a_p(i,k-1,j) = a_p(i,k-1,j)+0.5*a_dpn(i,k)*fnp(k)
604 a_p(i-1,k,j) = a_p(i-1,k,j)+0.5*a_dpn(i,k)*fnm(k)
605 a_p(i,k,j) = a_p(i,k,j)+0.5*a_dpn(i,k)*fnm(k)
606 a_dpn(i,k) = 0.
607 end do
608 end do
609 do i = i_start_up, i_end_up
610 a_dpn(i,kde) = 0.
611 a_p(i-1,3,j) = a_p(i-1,3,j)+0.5*a_dpn(i,1)*cf3
612 a_p(i,3,j) = a_p(i,3,j)+0.5*a_dpn(i,1)*cf3
613 a_p(i-1,2,j) = a_p(i-1,2,j)+0.5*a_dpn(i,1)*cf2
614 a_p(i,2,j) = a_p(i,2,j)+0.5*a_dpn(i,1)*cf2
615 a_p(i-1,1,j) = a_p(i-1,1,j)+0.5*a_dpn(i,1)*cf1
616 a_p(i,1,j) = a_p(i,1,j)+0.5*a_dpn(i,1)*cf1
617 a_dpn(i,1) = 0.
618 end do
619 endif
620 do k = k_start, k_end
621 do i = i_start_up, i_end_up
622 a_al(i-1,k,j) = a_al(i-1,k,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)*(pb(i,k,j)-pb(i-1,k,j))
623 a_al(i,k,j) = a_al(i,k,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)*(pb(i,k,j)-pb(i-1,k,j))
624 a_alt(i-1,k,j) = a_alt(i-1,k,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)*(p(i,k,j)-p(i-1,k,j))
625 a_alt(i,k,j) = a_alt(i,k,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)*(p(i,k,j)-p(i-1,k,j))
626 a_muu(i,j) = a_muu(i,j)+0.5*a_dpxy(i,k)*rdx*(ph(i,k+1,j)-ph(i-1,k+1,j)+ph(i,k,j)-ph(i-1,k,j)+(alt(i,k,j)+alt(i-1,k,j))*(p(i,&
627 &k,j)-p(i-1,k,j))+(al(i,k,j)+al(i-1,k,j))*(pb(i,k,j)-pb(i-1,k,j)))
628 a_p(i-1,k,j) = a_p(i-1,k,j)-0.5*a_dpxy(i,k)*rdx*muu(i,j)*(alt(i,k,j)+alt(i-1,k,j))
629 a_p(i,k,j) = a_p(i,k,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)*(alt(i,k,j)+alt(i-1,k,j))
630 a_ph(i-1,k+1,j) = a_ph(i-1,k+1,j)-0.5*a_dpxy(i,k)*rdx*muu(i,j)
631 a_ph(i,k+1,j) = a_ph(i,k+1,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)
632 a_ph(i-1,k,j) = a_ph(i-1,k,j)-0.5*a_dpxy(i,k)*rdx*muu(i,j)
633 a_ph(i,k,j) = a_ph(i,k,j)+0.5*a_dpxy(i,k)*rdx*muu(i,j)
634 a_dpxy(i,k) = 0.
635 end do
636 end do
637 do i = i_start_up, i_end_up
638 a_mudf(i-1,j) = a_mudf(i-1,j)+a_mudf_xy(i)*emdiv*dx
639 a_mudf(i,j) = a_mudf(i,j)-a_mudf_xy(i)*emdiv*dx
640 a_mudf_xy(i) = 0.
641 end do
642 do k = k_start, k_end
643 do i = i_start, i_endu
644 a_ru_tend(i,k,j) = a_ru_tend(i,k,j)+a_u(i,k,j)*dts
645 end do
646 end do
647 end do
648
649 end subroutine a_advance_uv
650
651
652 subroutine a_advance_w( w, a_w, rw_tend, a_rw_tend, ww, a_ww, u, a_u, v, a_v, mu1, a_mu1, mut, a_mut, muave, a_muave, muts, a_muts,&
653 & t_2ave, a_t_2ave, t_2, a_t_2, t_1, a_t_1, ph, a_ph, ph_1, a_ph_1, phb, ph_tend, a_ph_tend, ht, c2a, a_c2a, cqw, a_cqw, alt, &
654 &a_alt, alb, a, a_a, alpha, a_alpha, gamma, a_gamma, rdx, rdy, dts, t0, epssm, fnm, fnp, rdnw, rdn, cf1, cf2, cf3, msft, &
655 &config_flags, ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte )
656 !******************************************************************
657 !******************************************************************
658 !** This routine was generated by Automatic differentiation. **
659 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
660 !******************************************************************
661 !******************************************************************
662 !==============================================
663 ! all entries are defined explicitly
664 !==============================================
665 implicit none
666
667 !==============================================
668 ! declare arguments
669 !==============================================
670 integer, intent(in) :: ime
671 integer, intent(in) :: ims
672 integer, intent(in) :: jme
673 integer, intent(in) :: jms
674 integer, intent(in) :: kme
675 integer, intent(in) :: kms
676 real, intent(in) :: a(ims:ime,kms:kme,jms:jme)
677 real, intent(inout) :: a_a(ims:ime,kms:kme,jms:jme)
678 real, intent(inout) :: a_alpha(ims:ime,kms:kme,jms:jme)
679 real, intent(inout) :: a_alt(ims:ime,kms:kme,jms:jme)
680 real, intent(inout) :: a_c2a(ims:ime,kms:kme,jms:jme)
681 real, intent(inout) :: a_cqw(ims:ime,kms:kme,jms:jme)
682 real, intent(inout) :: a_gamma(ims:ime,kms:kme,jms:jme)
683 real, intent(inout) :: a_mu1(ims:ime,jms:jme)
684 real, intent(inout) :: a_muave(ims:ime,jms:jme)
685 real, intent(inout) :: a_mut(ims:ime,jms:jme)
686 real, intent(inout) :: a_muts(ims:ime,jms:jme)
687 real, intent(inout) :: a_ph(ims:ime,kms:kme,jms:jme)
688 real, intent(inout) :: a_ph_1(ims:ime,kms:kme,jms:jme)
689 real, intent(inout) :: a_ph_tend(ims:ime,kms:kme,jms:jme)
690 real, intent(inout) :: a_rw_tend(ims:ime,kms:kme,jms:jme)
691 real, intent(inout) :: a_t_1(ims:ime,kms:kme,jms:jme)
692 real, intent(inout) :: a_t_2(ims:ime,kms:kme,jms:jme)
693 real, intent(inout) :: a_t_2ave(ims:ime,kms:kme,jms:jme)
694 real, intent(inout) :: a_u(ims:ime,kms:kme,jms:jme)
695 real, intent(inout) :: a_v(ims:ime,kms:kme,jms:jme)
696 real, intent(inout) :: a_w(ims:ime,kms:kme,jms:jme)
697 real, intent(inout) :: a_ww(ims:ime,kms:kme,jms:jme)
698 real, intent(in) :: alb(ims:ime,kms:kme,jms:jme)
699 real, intent(in) :: alpha(ims:ime,kms:kme,jms:jme)
700 real, intent(in) :: alt(ims:ime,kms:kme,jms:jme)
701 real, intent(in) :: c2a(ims:ime,kms:kme,jms:jme)
702 real, intent(in) :: cf1
703 real, intent(in) :: cf2
704 real, intent(in) :: cf3
705 type (grid_config_rec_type), intent(in) :: config_flags
706 real, intent(in) :: cqw(ims:ime,kms:kme,jms:jme)
707 real, intent(in) :: dts
708 real, intent(in) :: epssm
709 real, intent(in) :: fnm(kms:kme)
710 real, intent(in) :: fnp(kms:kme)
711 real, intent(in) :: gamma(ims:ime,kms:kme,jms:jme)
712 real, intent(in) :: ht(ims:ime,jms:jme)
713 integer, intent(in) :: ide
714 integer, intent(in) :: ids
715 integer, intent(in) :: ite
716 integer, intent(in) :: its
717 integer, intent(in) :: jde
718 integer, intent(in) :: jds
719 integer, intent(in) :: jte
720 integer, intent(in) :: jts
721 integer, intent(in) :: kte
722 integer, intent(in) :: kts
723 real, intent(in) :: msft(ims:ime,jms:jme)
724 real, intent(in) :: mu1(ims:ime,jms:jme)
725 real, intent(in) :: muave(ims:ime,jms:jme)
726 real, intent(in) :: mut(ims:ime,jms:jme)
727 real, intent(in) :: muts(ims:ime,jms:jme)
728 real, intent(inout) :: ph(ims:ime,kms:kme,jms:jme)
729 real, intent(in) :: ph_1(ims:ime,kms:kme,jms:jme)
730 real, intent(in) :: ph_tend(ims:ime,kms:kme,jms:jme)
731 real, intent(in) :: phb(ims:ime,kms:kme,jms:jme)
732 real, intent(in) :: rdn(kms:kme)
733 real, intent(in) :: rdnw(kms:kme)
734 real, intent(in) :: rdx
735 real, intent(in) :: rdy
736 real, intent(in) :: rw_tend(ims:ime,kms:kme,jms:jme)
737 real, intent(in) :: t0
738 real, intent(in) :: t_1(ims:ime,kms:kme,jms:jme)
739 real, intent(in) :: t_2(ims:ime,kms:kme,jms:jme)
740 real, intent(inout) :: t_2ave(ims:ime,kms:kme,jms:jme)
741 real, intent(in) :: u(ims:ime,kms:kme,jms:jme)
742 real, intent(in) :: v(ims:ime,kms:kme,jms:jme)
743 real, intent(inout) :: w(ims:ime,kms:kme,jms:jme)
744 real, intent(in) :: ww(ims:ime,kms:kme,jms:jme)
745
746 !==============================================
747 ! declare local variables
748 !==============================================
749 real a_mut_inv(its:ite)
750 real a_rhs(its:ite,kts:kte)
751 real a_wdwn(its:ite,kts:kte)
752 integer i
753 integer i_end
754 integer i_start
755 integer j
756 integer j_end
757 integer j_start
758 integer k
759 integer k_end
760 real msft_inv(its:ite)
761 real mut_inv(its:ite)
762 real rhs(its:ite,kts:kte)
763 real t_2aveh(ims:ime,kms:kme,jms:jme)
764 real wdwn(its:ite,kts:kte)
765 real wh(ims:ime,kms:kme,jms:jme)
766 !mzz
767 real keep1_w(ims:ime,kms:kme,jms:jme)
768 real keep2_w(ims:ime,kms:kme,jms:jme)
769 real keep3_w(ims:ime,kms:kme,jms:jme)
770 !mzz
771
772 !----------------------------------------------
773 ! SAVE REQUIRED INPUT VARIABLES
774 !----------------------------------------------
775 wh(:,:,:) = w(:,:,:)
776 t_2aveh(:,:,:) = t_2ave(:,:,:)
777
778 !----------------------------------------------
779 ! RESET LOCAL ADJOINT VARIABLES
780 !----------------------------------------------
781 a_mut_inv(:) = 0.
782 a_rhs(:,:) = 0.
783 a_wdwn(:,:) = 0.
784
785 !----------------------------------------------
786 ! ROUTINE BODY
787 !----------------------------------------------
788 i_start = its
789 ! recompute : i_start
790 i_end = ite
791 ! recompute : i_end
792 j_start = jts
793 ! recompute : j_start
794 j_end = jte
795 ! recompute : j_end
796 k_end = kte-1
797 ! recompute : k_end
798 if (j_end .eq. jde) then
799 j_end = j_end-1
800 endif
801 ! recompute : j_end
802 if (i_end .eq. ide) then
803 i_end = i_end-1
804 endif
805 ! recompute : i_end
806 if ((config_flags%specified .or. config_flags%nested) .and. its .eq. ids) then
807 i_start = i_start+1
808 endif
809 ! recompute : i_start
810 if ((config_flags%specified .or. config_flags%nested) .and. ite .eq. ide) then
811 i_end = i_end-1
812 endif
813 ! recompute : i_end
814 if ((config_flags%specified .or. config_flags%nested) .and. jts .eq. jds) then
815 j_start = j_start+1
816 endif
817 ! recompute : j_start
818 if ((config_flags%specified .or. config_flags%nested) .and. jte .eq. jde) then
819 j_end = j_end-1
820 endif
821 ! recompute : j_end
822 do j = j_start, j_end
823 do i = i_start, i_end
824 mut_inv(i) = 1./mut(i,j)
825 msft_inv(i) = 1./msft(i,j)
826 end do
827 ! recompute : msft_inv,mut_inv
828 do k = 1, k_end
829 do i = i_start, i_end
830 t_2ave(i,k,j) = 0.5*((1.+epssm)*t_2(i,k,j)+(1.-epssm)*t_2ave(i,k,j))
831 t_2ave(i,k,j) = (t_2ave(i,k,j)-mu1(i,j)*t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j)))
832 end do
833 end do
834 ! recompute : t_2ave
835 do k = 2, k_end+1
836 do i = i_start, i_end
837 wdwn(i,k) = 0.5*(ww(i,k,j)+ww(i,k-1,j))*rdnw(k-1)*(ph_1(i,k,j)-ph_1(i,k-1,j)+phb(i,k,j)-phb(i,k-1,j))
838 rhs(i,k) = dts*(ph_tend(i,k,j)+0.5*g*(1.-epssm)*w(i,k,j))
839 end do
840 end do
841 ! recompute : rhs,wdwn
842 do k = 2, k_end
843 do i = i_start, i_end
844 rhs(i,k) = rhs(i,k)-dts*(fnm(k)*wdwn(i,k+1)+fnp(k)*wdwn(i,k))
845 end do
846 end do
847 ! recompute : rhs
848 do k = 2, k_end+1
849 do i = i_start, i_end
850 rhs(i,k) = ph(i,k,j)+msft(i,j)*rhs(i,k)*mut_inv(i)
851 end do
852 end do
853 ! recompute : rhs
854 do i = i_start, i_end
855 w(i,1,j) = 0.5*rdy*((ht(i,j+1)-ht(i,j))*(cf1*v(i,1,j+1)+cf2*v(i,2,j+1)+cf3*v(i,3,j+1))+(ht(i,j)-ht(i,j-1))*(cf1*v(i,1,j)+cf2*&
856 &v(i,2,j)+cf3*v(i,3,j)))+0.5*rdx*((ht(i+1,j)-ht(i,j))*(cf1*u(i+1,1,j)+cf2*u(i+1,2,j)+cf3*u(i+1,3,j))+(ht(i,j)-ht(i-1,j))*(cf1*&
857 &u(i,1,j)+cf2*u(i,2,j)+cf3*u(i,3,j)))
858 end do
859 ! recompute : w
860 do k = 2, k_end
861 do i = i_start, i_end
862 w(i,k,j) = w(i,k,j)+dts*rw_tend(i,k,j)+msft_inv(i)*cqw(i,k,j)*0.5*dts*g*mut_inv(i)*rdn(k)*(c2a(i,k,j)*rdnw(k)*((1.+epssm)*&
863 &(rhs(i,k+1)-rhs(i,k))+(1.-epssm)*(ph(i,k+1,j)-ph(i,k,j)))-c2a(i,k-1,j)*rdnw(k-1)*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-&
864 &epssm)*(ph(i,k,j)-ph(i,k-1,j))))+dts*g*msft_inv(i)*(rdn(k)*(c2a(i,k,j)*alt(i,k,j)*t_2ave(i,k,j)-c2a(i,k-1,j)*alt(i,k-1,j)*&
865 &t_2ave(i,k-1,j))+(rdn(k)*(c2a(i,k,j)*alb(i,k,j)-c2a(i,k-1,j)*alb(i,k-1,j))*mut_inv(i)-1.)*muave(i,j))
866 end do
867 end do
868 ! recompute : w
869 k = k_end+1
870 ! recompute : k
871 do i = i_start, i_end
872 w(i,k,j) = w(i,k,j)+dts*rw_tend(i,k,j)+msft_inv(i)*((-(0.5*dts*g*mut_inv(i)*rdnw(k-1)**2*2.*c2a(i,k-1,j)*((1.+epssm)*(rhs(i,k)-&
873 &rhs(i,k-1))+(1.-epssm)*(ph(i,k,j)-ph(i,k-1,j)))))-dts*g*(2.*rdnw(k-1)*c2a(i,k-1,j)*alt(i,k-1,j)*t_2ave(i,k-1,j)+(1.+2.*rdnw(k-&
874 &1)*c2a(i,k-1,j)*alb(i,k-1,j)*mut_inv(i))*muave(i,j)))
875 end do
876 ! recompute : w
877 do k = 2, k_end+1
878 do i = i_start, i_end
879 !mzz
880 keep1_w(i,k-1,j) =w(i,k-1,j)
881 keep2_w(i,k,j) =w(i,k,j)
882 !mzz
883 w(i,k,j) = (w(i,k,j)-a(i,k,j)*w(i,k-1,j))*alpha(i,k,j)
884 end do
885 end do
886 ! recompute : w
887 do k = k_end, 2, -1
888 do i = i_start, i_end
889 !mzz
890 keep3_w(i,k+1,j) =w(i,k+1,j)
891 !mzz
892 w(i,k,j) = w(i,k,j)-gamma(i,k,j)*w(i,k+1,j)
893 end do
894 end do
895 ! recompute : w
896 do k = 2, k_end+1
897 do i = i_start, i_end
898 a_muts(i,j) = a_muts(i,j)-a_ph(i,k,j)*(0.5*msft(i,j)*dts*g*(1.+epssm)*w(i,k,j)/(muts(i,j)*muts(i,j)))
899 a_rhs(i,k) = a_rhs(i,k)+a_ph(i,k,j)
900 a_w(i,k,j) = a_w(i,k,j)+a_ph(i,k,j)*(0.5*msft(i,j)*dts*g*(1.+epssm)/muts(i,j))
901 a_ph(i,k,j) = 0.
902 end do
903 end do
904 ! recdepend vars : cf1,cf2,cf3,ht,i_end,i_start,j,rdx,rdy,u,v
905 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1361
906 ! recompute vars : w
907
908 do i = i_start, i_end
909 w(i,1,j) = 0.5*rdy*((ht(i,j+1)-ht(i,j))*(cf1*v(i,1,j+1)+cf2*v(i,2,j+1)+cf3*v(i,3,j+1))+(ht(i,j)-ht(i,j-1))*(cf1*v(i,1,j)+cf2*&
910 &v(i,2,j)+cf3*v(i,3,j)))+0.5*rdx*((ht(i+1,j)-ht(i,j))*(cf1*u(i+1,1,j)+cf2*u(i+1,2,j)+cf3*u(i+1,3,j))+(ht(i,j)-ht(i-1,j))*(cf1*&
911 &u(i,1,j)+cf2*u(i,2,j)+cf3*u(i,3,j)))
912 end do
913 ! recompute vars : w
914 ! recdepend vars : alb,alt,c2a,cqw,dts,epssm,g,i_end,i_start,j,k_end,ms
915 ! ft_inv,muave,mut_inv,ph,rdn,rdnw,rhs,rw_tend,t_2ave,w
916 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1382
917 ! recompute vars : w
918 do k = 2, k_end
919 do i = i_start, i_end
920 w(i,k,j) = w(i,k,j)+dts*rw_tend(i,k,j)+msft_inv(i)*cqw(i,k,j)*0.5*dts*g*mut_inv(i)*rdn(k)*(c2a(i,k,j)*rdnw(k)*((1.+epssm)*&
921 &(rhs(i,k+1)-rhs(i,k))+(1.-epssm)*(ph(i,k+1,j)-ph(i,k,j)))-c2a(i,k-1,j)*rdnw(k-1)*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-&
922 &epssm)*(ph(i,k,j)-ph(i,k-1,j))))+dts*g*msft_inv(i)*(rdn(k)*(c2a(i,k,j)*alt(i,k,j)*t_2ave(i,k,j)-c2a(i,k-1,j)*alt(i,k-1,j)*&
923 &t_2ave(i,k-1,j))+(rdn(k)*(c2a(i,k,j)*alb(i,k,j)-c2a(i,k-1,j)*alb(i,k-1,j))*mut_inv(i)-1.)*muave(i,j))
924 end do
925 end do
926 ! recompute vars : w
927 ! recdepend vars : k_end,w
928 ! recompute pos : ASSIGN_STMT module_small_step_em.f90:1404
929 ! recompute vars : k
930 k = k_end+1
931 ! recompute vars : k
932 ! recdepend vars : alb,alt,c2a,dts,epssm,g,i_end,i_start,j,k,msft_inv,m
933 ! uave,mut_inv,ph,rdnw,rhs,rw_tend,t_2ave,w
934 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1406
935 ! recompute vars : w
936 do i = i_start, i_end
937 w(i,k,j) = w(i,k,j)+dts*rw_tend(i,k,j)+msft_inv(i)*((-(0.5*dts*g*mut_inv(i)*rdnw(k-1)**2*2.*c2a(i,k-1,j)*((1.+epssm)*(rhs(i,k)-&
938 &rhs(i,k-1))+(1.-epssm)*(ph(i,k,j)-ph(i,k-1,j)))))-dts*g*(2.*rdnw(k-1)*c2a(i,k-1,j)*alt(i,k-1,j)*t_2ave(i,k-1,j)+(1.+2.*rdnw(k-&
939 &1)*c2a(i,k-1,j)*alb(i,k-1,j)*mut_inv(i))*muave(i,j)))
940 end do
941 ! recompute vars : w
942 ! recdepend vars : a,alpha,i_end,i_start,j,k_end,w
943 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1418
944 ! recompute vars : w
945 do k = 2, k_end+1
946 do i = i_start, i_end
947 w(i,k,j) = (w(i,k,j)-a(i,k,j)*w(i,k-1,j))*alpha(i,k,j)
948 end do
949 end do
950 ! recompute vars : w
951 do k = k_end, 3, -1
952 do i = i_start, i_end
953 w(i,k,j) = w(i,k,j)-gamma(i,k,j)*w(i,k+1,j)
954 end do
955 end do
956 do k = 2, k_end
957 do i = i_start, i_end
958 !mzz
959 w(i,k+1,j) =keep3_w(i,k+1,j)
960 !mzz
961 a_gamma(i,k,j) = a_gamma(i,k,j)-a_w(i,k,j)*w(i,k+1,j)
962 a_w(i,k+1,j) = a_w(i,k+1,j)-a_w(i,k,j)*gamma(i,k,j)
963 end do
964 end do
965 ! recdepend vars : cf1,cf2,cf3,ht,i_end,i_start,j,rdx,rdy,u,v
966 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1361
967 ! recompute vars : w
968 do i = i_start, i_end
969 w(i,1,j) = 0.5*rdy*((ht(i,j+1)-ht(i,j))*(cf1*v(i,1,j+1)+cf2*v(i,2,j+1)+cf3*v(i,3,j+1))+(ht(i,j)-ht(i,j-1))*(cf1*v(i,1,j)+cf2*&
970 &v(i,2,j)+cf3*v(i,3,j)))+0.5*rdx*((ht(i+1,j)-ht(i,j))*(cf1*u(i+1,1,j)+cf2*u(i+1,2,j)+cf3*u(i+1,3,j))+(ht(i,j)-ht(i-1,j))*(cf1*&
971 &u(i,1,j)+cf2*u(i,2,j)+cf3*u(i,3,j)))
972 end do
973 ! recompute vars : w
974 ! recdepend vars : alb,alt,c2a,cqw,dts,epssm,g,i_end,i_start,j,k_end,ms
975 ! ft_inv,muave,mut_inv,ph,rdn,rdnw,rhs,rw_tend,t_2ave,w
976 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1382
977 ! recompute vars : w
978 do k = 2, k_end
979 do i = i_start, i_end
980 w(i,k,j) = w(i,k,j)+dts*rw_tend(i,k,j)+msft_inv(i)*cqw(i,k,j)*0.5*dts*g*mut_inv(i)*rdn(k)*(c2a(i,k,j)*rdnw(k)*((1.+epssm)*&
981 &(rhs(i,k+1)-rhs(i,k))+(1.-epssm)*(ph(i,k+1,j)-ph(i,k,j)))-c2a(i,k-1,j)*rdnw(k-1)*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-&
982 &epssm)*(ph(i,k,j)-ph(i,k-1,j))))+dts*g*msft_inv(i)*(rdn(k)*(c2a(i,k,j)*alt(i,k,j)*t_2ave(i,k,j)-c2a(i,k-1,j)*alt(i,k-1,j)*&
983 &t_2ave(i,k-1,j))+(rdn(k)*(c2a(i,k,j)*alb(i,k,j)-c2a(i,k-1,j)*alb(i,k-1,j))*mut_inv(i)-1.)*muave(i,j))
984 end do
985 end do
986 ! recompute vars : w
987 ! recdepend vars : k_end,w
988 ! recompute pos : ASSIGN_STMT module_small_step_em.f90:1404
989 ! recompute vars : k
990 k = k_end+1
991 ! recompute vars : k
992 ! recdepend vars : alb,alt,c2a,dts,epssm,g,i_end,i_start,j,k,msft_inv,m
993 ! uave,mut_inv,ph,rdnw,rhs,rw_tend,t_2ave,w
994 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1406
995 ! recompute vars : w
996 do i = i_start, i_end
997 w(i,k,j) = w(i,k,j)+dts*rw_tend(i,k,j)+msft_inv(i)*((-(0.5*dts*g*mut_inv(i)*rdnw(k-1)**2*2.*c2a(i,k-1,j)*((1.+epssm)*(rhs(i,k)-&
998 &rhs(i,k-1))+(1.-epssm)*(ph(i,k,j)-ph(i,k-1,j)))))-dts*g*(2.*rdnw(k-1)*c2a(i,k-1,j)*alt(i,k-1,j)*t_2ave(i,k-1,j)+(1.+2.*rdnw(k-&
999 &1)*c2a(i,k-1,j)*alb(i,k-1,j)*mut_inv(i))*muave(i,j)))
1000 end do
1001 ! recompute vars : w
1002 do k = 2, k_end
1003 do i = i_start, i_end
1004 w(i,k,j) = (w(i,k,j)-a(i,k,j)*w(i,k-1,j))*alpha(i,k,j)
1005 end do
1006 end do
1007 do k = k_end+1, 2, -1
1008 do i = i_start, i_end
1009 !mzz
1010 w(i,k-1,j)=keep1_w(i,k-1,j)
1011 w(i,k,j)=keep2_w(i,k,j)
1012 !mzz
1013 a_a(i,k,j) = a_a(i,k,j)-a_w(i,k,j)*w(i,k-1,j)*alpha(i,k,j)
1014 a_alpha(i,k,j) = a_alpha(i,k,j)+a_w(i,k,j)*(w(i,k,j)-a(i,k,j)*w(i,k-1,j))
1015 a_w(i,k-1,j) = a_w(i,k-1,j)-a_w(i,k,j)*a(i,k,j)*alpha(i,k,j)
1016 a_w(i,k,j) = a_w(i,k,j)*alpha(i,k,j)
1017 end do
1018 end do
1019
1020 ! recdepend vars : k_end
1021 ! recompute pos : ASSIGN_STMT module_small_step_em.f90:1404
1022 ! recompute vars : k
1023 k = k_end+1
1024 ! recompute vars : k
1025 do i = i_start, i_end
1026 a_alt(i,k-1,j) = a_alt(i,k-1,j)-2*a_w(i,k,j)*msft_inv(i)*dts*g*rdnw(k-1)*c2a(i,k-1,j)*t_2ave(i,k-1,j)
1027 a_c2a(i,k-1,j) = a_c2a(i,k-1,j)-a_w(i,k,j)*msft_inv(i)*dts*g*(mut_inv(i)*rdnw(k-1)**2*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-&
1028 &epssm)*(ph(i,k,j)-ph(i,k-1,j)))+2*rdnw(k-1)*alt(i,k-1,j)*t_2ave(i,k-1,j)+2*rdnw(k-1)*alb(i,k-1,j)*mut_inv(i)*muave(i,j))
1029 a_muave(i,j) = a_muave(i,j)-a_w(i,k,j)*msft_inv(i)*dts*g*(1+2.*rdnw(k-1)*c2a(i,k-1,j)*alb(i,k-1,j)*mut_inv(i))
1030 a_mut_inv(i) = a_mut_inv(i)-a_w(i,k,j)*msft_inv(i)*(dts*g*rdnw(k-1)**2*c2a(i,k-1,j)*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-&
1031 &epssm)*(ph(i,k,j)-ph(i,k-1,j)))+2*dts*g*rdnw(k-1)*c2a(i,k-1,j)*alb(i,k-1,j)*muave(i,j))
1032 a_ph(i,k-1,j) = a_ph(i,k-1,j)+a_w(i,k,j)*msft_inv(i)*dts*g*mut_inv(i)*rdnw(k-1)**2*c2a(i,k-1,j)*(1.-epssm)
1033 a_ph(i,k,j) = a_ph(i,k,j)-a_w(i,k,j)*msft_inv(i)*dts*g*mut_inv(i)*rdnw(k-1)**2*c2a(i,k-1,j)*(1.-epssm)
1034 a_rhs(i,k-1) = a_rhs(i,k-1)+a_w(i,k,j)*msft_inv(i)*dts*g*mut_inv(i)*rdnw(k-1)**2*c2a(i,k-1,j)*(1.+epssm)
1035 a_rhs(i,k) = a_rhs(i,k)-a_w(i,k,j)*msft_inv(i)*dts*g*mut_inv(i)*rdnw(k-1)**2*c2a(i,k-1,j)*(1+epssm)
1036 a_rw_tend(i,k,j) = a_rw_tend(i,k,j)+a_w(i,k,j)*dts
1037 a_t_2ave(i,k-1,j) = a_t_2ave(i,k-1,j)-2*a_w(i,k,j)*msft_inv(i)*dts*g*rdnw(k-1)*c2a(i,k-1,j)*alt(i,k-1,j)
1038 end do
1039 do k = 2, k_end
1040 do i = i_start, i_end
1041 a_alt(i,k-1,j) = a_alt(i,k-1,j)-a_w(i,k,j)*dts*g*msft_inv(i)*rdn(k)*c2a(i,k-1,j)*t_2ave(i,k-1,j)
1042 a_alt(i,k,j) = a_alt(i,k,j)+a_w(i,k,j)*dts*g*msft_inv(i)*rdn(k)*c2a(i,k,j)*t_2ave(i,k,j)
1043 a_c2a(i,k-1,j) = a_c2a(i,k-1,j)+a_w(i,k,j)*((-(0.5*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*rdnw(k-1)*((1.+epssm)*&
1044 &(rhs(i,k)-rhs(i,k-1))+(1.-epssm)*(ph(i,k,j)-ph(i,k-1,j)))))+dts*g*msft_inv(i)*((-(rdn(k)*alt(i,k-1,j)*t_2ave(i,k-1,j)))-&
1045 &rdn(k)*alb(i,k-1,j)*mut_inv(i)*muave(i,j)))
1046 a_c2a(i,k,j) = a_c2a(i,k,j)+a_w(i,k,j)*(0.5*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*rdnw(k)*((1.+epssm)*(rhs(i,k+1)-&
1047 &rhs(i,k))+(1.-epssm)*(ph(i,k+1,j)-ph(i,k,j)))+dts*g*msft_inv(i)*rdn(k)*(alt(i,k,j)*t_2ave(i,k,j)+alb(i,k,j)*mut_inv(i)*muave(i,j)))
1048 a_cqw(i,k,j) = a_cqw(i,k,j)+0.5*a_w(i,k,j)*msft_inv(i)*dts*g*mut_inv(i)*rdn(k)*(c2a(i,k,j)*rdnw(k)*((1.+epssm)*(rhs(i,k+1)-&
1049 &rhs(i,k))+(1.-epssm)*(ph(i,k+1,j)-ph(i,k,j)))-c2a(i,k-1,j)*rdnw(k-1)*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-epssm)*(ph(i,k,j)&
1050 &-ph(i,k-1,j))))
1051 a_muave(i,j) = a_muave(i,j)+a_w(i,k,j)*dts*g*msft_inv(i)*((-1)+rdn(k)*(c2a(i,k,j)*alb(i,k,j)-c2a(i,k-1,j)*alb(i,k-1,j))*&
1052 &mut_inv(i))
1053 a_mut_inv(i) = a_mut_inv(i)+a_w(i,k,j)*(0.5*msft_inv(i)*cqw(i,k,j)*dts*g*rdn(k)*(c2a(i,k,j)*rdnw(k)*((1.+epssm)*(rhs(i,k+1)-&
1054 &rhs(i,k))+(1.-epssm)*(ph(i,k+1,j)-ph(i,k,j)))-c2a(i,k-1,j)*rdnw(k-1)*((1.+epssm)*(rhs(i,k)-rhs(i,k-1))+(1.-epssm)*(ph(i,k,j)&
1055 &-ph(i,k-1,j))))+dts*g*msft_inv(i)*rdn(k)*(c2a(i,k,j)*alb(i,k,j)-c2a(i,k-1,j)*alb(i,k-1,j))*muave(i,j))
1056 a_ph(i,k-1,j) = a_ph(i,k-1,j)+0.5*a_w(i,k,j)*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*c2a(i,k-1,j)*rdnw(k-1)*(1.-epssm)
1057 a_ph(i,k+1,j) = a_ph(i,k+1,j)+0.5*a_w(i,k,j)*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*c2a(i,k,j)*rdnw(k)*(1.-epssm)
1058 a_ph(i,k,j) = a_ph(i,k,j)-0.5*a_w(i,k,j)*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*(c2a(i,k,j)*rdnw(k)*(1.-epssm)+c2a(i,&
1059 &k-1,j)*rdnw(k-1)*(1.-epssm))
1060 a_rhs(i,k-1) = a_rhs(i,k-1)+0.5*a_w(i,k,j)*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*c2a(i,k-1,j)*rdnw(k-1)*(1.+epssm)
1061 a_rhs(i,k+1) = a_rhs(i,k+1)+0.5*a_w(i,k,j)*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*c2a(i,k,j)*rdnw(k)*(1+epssm)
1062 a_rhs(i,k) = a_rhs(i,k)-0.5*a_w(i,k,j)*msft_inv(i)*cqw(i,k,j)*dts*g*mut_inv(i)*rdn(k)*(c2a(i,k,j)*rdnw(k)*(1.+epssm)+c2a(i,k-&
1063 &1,j)*rdnw(k-1)*(1+epssm))
1064 a_rw_tend(i,k,j) = a_rw_tend(i,k,j)+a_w(i,k,j)*dts
1065 a_t_2ave(i,k-1,j) = a_t_2ave(i,k-1,j)-a_w(i,k,j)*dts*g*msft_inv(i)*rdn(k)*c2a(i,k-1,j)*alt(i,k-1,j)
1066 a_t_2ave(i,k,j) = a_t_2ave(i,k,j)+a_w(i,k,j)*dts*g*msft_inv(i)*rdn(k)*c2a(i,k,j)*alt(i,k,j)
1067 end do
1068 end do
1069 do i = i_start, i_end
1070 a_u(i+1,3,j) = a_u(i+1,3,j)+0.5*a_w(i,1,j)*rdx*(ht(i+1,j)-ht(i,j))*cf3
1071 a_u(i,3,j) = a_u(i,3,j)+0.5*a_w(i,1,j)*rdx*(ht(i,j)-ht(i-1,j))*cf3
1072 a_u(i+1,2,j) = a_u(i+1,2,j)+0.5*a_w(i,1,j)*rdx*(ht(i+1,j)-ht(i,j))*cf2
1073 a_u(i,2,j) = a_u(i,2,j)+0.5*a_w(i,1,j)*rdx*(ht(i,j)-ht(i-1,j))*cf2
1074 a_u(i+1,1,j) = a_u(i+1,1,j)+0.5*a_w(i,1,j)*rdx*(ht(i+1,j)-ht(i,j))*cf1
1075 a_u(i,1,j) = a_u(i,1,j)+0.5*a_w(i,1,j)*rdx*(ht(i,j)-ht(i-1,j))*cf1
1076 a_v(i,3,j+1) = a_v(i,3,j+1)+0.5*a_w(i,1,j)*rdy*(ht(i,j+1)-ht(i,j))*cf3
1077 a_v(i,3,j) = a_v(i,3,j)+0.5*a_w(i,1,j)*rdy*(ht(i,j)-ht(i,j-1))*cf3
1078 a_v(i,2,j+1) = a_v(i,2,j+1)+0.5*a_w(i,1,j)*rdy*(ht(i,j+1)-ht(i,j))*cf2
1079 a_v(i,2,j) = a_v(i,2,j)+0.5*a_w(i,1,j)*rdy*(ht(i,j)-ht(i,j-1))*cf2
1080 a_v(i,1,j+1) = a_v(i,1,j+1)+0.5*a_w(i,1,j)*rdy*(ht(i,j+1)-ht(i,j))*cf1
1081 a_v(i,1,j) = a_v(i,1,j)+0.5*a_w(i,1,j)*rdy*(ht(i,j)-ht(i,j-1))*cf1
1082 a_w(i,1,j) = 0.
1083 end do
1084 w(:,:,:) = wh(:,:,:)
1085 ! recdepend vars : dts,epssm,g,i_end,i_start,j,k_end,ph_tend,w
1086 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1335
1087 ! recompute vars : rhs
1088 do k = 2, k_end+1
1089 do i = i_start, i_end
1090 rhs(i,k) = dts*(ph_tend(i,k,j)+0.5*g*(1.-epssm)*w(i,k,j))
1091 end do
1092 end do
1093 ! recompute vars : rhs
1094 ! recdepend vars : dts,fnm,fnp,i_end,i_start,k_end,rhs,wdwn
1095 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:1343
1096 ! recompute vars : rhs
1097 do k = 2, k_end
1098 do i = i_start, i_end
1099 rhs(i,k) = rhs(i,k)-dts*(fnm(k)*wdwn(i,k+1)+fnp(k)*wdwn(i,k))
1100 end do
1101 end do
1102 ! recompute vars : rhs
1103 do k = 2, k_end+1
1104 do i = i_start, i_end
1105 a_mut_inv(i) = a_mut_inv(i)+a_rhs(i,k)*msft(i,j)*rhs(i,k)
1106 a_ph(i,k,j) = a_ph(i,k,j)+a_rhs(i,k)
1107 a_rhs(i,k) = a_rhs(i,k)*msft(i,j)*mut_inv(i)
1108 end do
1109 end do
1110 do k = 2, k_end
1111 do i = i_start, i_end
1112 a_wdwn(i,k+1) = a_wdwn(i,k+1)-a_rhs(i,k)*dts*fnm(k)
1113 a_wdwn(i,k) = a_wdwn(i,k)-a_rhs(i,k)*dts*fnp(k)
1114 end do
1115 end do
1116 do k = 2, k_end+1
1117 do i = i_start, i_end
1118 a_ph_tend(i,k,j) = a_ph_tend(i,k,j)+a_rhs(i,k)*dts
1119 a_w(i,k,j) = a_w(i,k,j)+0.5*a_rhs(i,k)*dts*g*(1.-epssm)
1120 a_rhs(i,k) = 0.
1121 a_ph_1(i,k-1,j) = a_ph_1(i,k-1,j)-0.5*a_wdwn(i,k)*(ww(i,k,j)+ww(i,k-1,j))*rdnw(k-1)
1122 a_ph_1(i,k,j) = a_ph_1(i,k,j)+0.5*a_wdwn(i,k)*(ww(i,k,j)+ww(i,k-1,j))*rdnw(k-1)
1123 a_ww(i,k-1,j) = a_ww(i,k-1,j)+0.5*a_wdwn(i,k)*rdnw(k-1)*(ph_1(i,k,j)-ph_1(i,k-1,j)+phb(i,k,j)-phb(i,k-1,j))
1124 a_ww(i,k,j) = a_ww(i,k,j)+0.5*a_wdwn(i,k)*rdnw(k-1)*(ph_1(i,k,j)-ph_1(i,k-1,j)+phb(i,k,j)-phb(i,k-1,j))
1125 a_wdwn(i,k) = 0.
1126 end do
1127 end do
1128 t_2ave(:,:,:) = t_2aveh(:,:,:)
1129 do k = 1, k_end
1130 do i = i_start, i_end
1131 t_2ave(i,k,j) = 0.5*((1.+epssm)*t_2(i,k,j)+(1.-epssm)*t_2ave(i,k,j))
1132 ! recompute : t_2ave
1133 a_mu1(i,j) = a_mu1(i,j)-a_t_2ave(i,k,j)*(t_1(i,k,j)/(muts(i,j)*(t0+t_1(i,k,j))))
1134 a_muts(i,j) = a_muts(i,j)-a_t_2ave(i,k,j)*((t_2ave(i,k,j)-mu1(i,j)*t_1(i,k,j))*(t0+t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j))*&
1135 &muts(i,j)*(t0+t_1(i,k,j))))
1136 a_t_1(i,k,j) = a_t_1(i,k,j)-a_t_2ave(i,k,j)*(mu1(i,j)/(muts(i,j)*(t0+t_1(i,k,j)))+(t_2ave(i,k,j)-mu1(i,j)*t_1(i,k,j))*muts(i,&
1137 &j)/(muts(i,j)*(t0+t_1(i,k,j))*muts(i,j)*(t0+t_1(i,k,j))))
1138 a_t_2ave(i,k,j) = a_t_2ave(i,k,j)/(muts(i,j)*(t0+t_1(i,k,j)))
1139 a_t_2(i,k,j) = a_t_2(i,k,j)+0.5*a_t_2ave(i,k,j)*(1+epssm)
1140 a_t_2ave(i,k,j) = 0.5*a_t_2ave(i,k,j)*(1.-epssm)
1141 end do
1142 end do
1143 do i = i_start, i_end
1144 a_mut(i,j) = a_mut(i,j)-a_mut_inv(i)/(mut(i,j)*mut(i,j))
1145 a_mut_inv(i) = 0.
1146 end do
1147 end do
1148
1149 !----------------------------------------------
1150 ! FREE DYNAMIC MEMORY
1151 !----------------------------------------------
1152
1153 end subroutine a_advance_w
1154
1155
1156 subroutine a_calc_coef_w( a, a_a, alpha, a_alpha, gamma, a_gamma, mut, a_mut, cqw, a_cqw, rdn, rdnw, c2a, a_c2a, dts, g, epssm, &
1157 &ide, jde, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte )
1158 !******************************************************************
1159 !******************************************************************
1160 !** This routine was generated by Automatic differentiation. **
1161 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
1162 !******************************************************************
1163 !******************************************************************
1164 !==============================================
1165 ! all entries are defined explicitly
1166 !==============================================
1167 implicit none
1168
1169 !==============================================
1170 ! declare arguments
1171 !==============================================
1172 integer, intent(in) :: ime
1173 integer, intent(in) :: ims
1174 integer, intent(in) :: jme
1175 integer, intent(in) :: jms
1176 integer, intent(in) :: kme
1177 integer, intent(in) :: kms
1178 real, intent(inout) :: a(ims:ime,kms:kme,jms:jme)
1179 real, intent(inout) :: a_a(ims:ime,kms:kme,jms:jme)
1180 real, intent(inout) :: a_alpha(ims:ime,kms:kme,jms:jme)
1181 real, intent(inout) :: a_c2a(ims:ime,kms:kme,jms:jme)
1182 real, intent(inout) :: a_cqw(ims:ime,kms:kme,jms:jme)
1183 real, intent(inout) :: a_gamma(ims:ime,kms:kme,jms:jme)
1184 real, intent(inout) :: a_mut(ims:ime,jms:jme)
1185 real, intent(inout) :: alpha(ims:ime,kms:kme,jms:jme)
1186 real, intent(in) :: c2a(ims:ime,kms:kme,jms:jme)
1187 real, intent(in) :: cqw(ims:ime,kms:kme,jms:jme)
1188 real, intent(in) :: dts
1189 real, intent(in) :: epssm
1190 real, intent(in) :: g
1191 real, intent(inout) :: gamma(ims:ime,kms:kme,jms:jme)
1192 integer, intent(in) :: ide
1193 integer, intent(in) :: ite
1194 integer, intent(in) :: its
1195 integer, intent(in) :: jde
1196 integer, intent(in) :: jte
1197 integer, intent(in) :: jts
1198 integer, intent(in) :: kde
1199 real, intent(in) :: mut(ims:ime,jms:jme)
1200 real, intent(in) :: rdn(kms:kme)
1201 real, intent(in) :: rdnw(kms:kme)
1202
1203 !==============================================
1204 ! declare local variables
1205 !==============================================
1206 real a_b
1207 real a_c
1208 real a_cof(ims:ime)
1209 real b
1210 real c
1211 real cof(ims:ime)
1212 integer i
1213 integer i_end
1214 integer i_start
1215 integer j
1216 integer j_end
1217 integer j_start
1218 integer k
1219 integer k1
1220
1221 !----------------------------------------------
1222 ! RESET LOCAL ADJOINT VARIABLES
1223 !----------------------------------------------
1224 a_b = 0.
1225 a_c = 0.
1226 a_cof(:) = 0.
1227
1228 !----------------------------------------------
1229 ! ROUTINE BODY
1230 !----------------------------------------------
1231 i_start = its
1232 ! recompute : i_start
1233 i_end = ite
1234 ! recompute : i_end
1235 j_start = jts
1236 ! recompute : j_start
1237 j_end = jte
1238 ! recompute : j_end
1239 if (j_end .eq. jde) then
1240 j_end = j_end-1
1241 endif
1242 ! recompute : j_end
1243 if (i_end .eq. ide) then
1244 i_end = i_end-1
1245 endif
1246 ! recompute : i_end
1247 do j = j_start, j_end
1248 a_b = 0.
1249 a_c = 0.
1250 do i = i_start, i_end
1251 cof(i) = (0.5*dts*g*(1.+epssm)/mut(i,j))**2
1252 a(i,2,j) = 0.
1253 a(i,kde,j) = -(2.*cof(i)*rdnw(kde-1)**2*c2a(i,kde-1,j))
1254 gamma(i,1,j) = 0.
1255 end do
1256 ! recompute : cof,gamma
1257 do k = 3, kde-1
1258 do i = i_start, i_end
1259 a(i,k,j) = -(cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k-1)*c2a(i,k-1,j))
1260 end do
1261 end do
1262 ! recompute : a
1263 do k = 2, kde-1
1264 do i = i_start, i_end
1265 b = 1.+cqw(i,k,j)*cof(i)*rdn(k)*(rdnw(k)*c2a(i,k,j)+rdnw(k-1)*c2a(i,k-1,j))
1266 c = -(cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k)*c2a(i,k,j))
1267 alpha(i,k,j) = 1./(b-a(i,k,j)*gamma(i,k-1,j))
1268 gamma(i,k,j) = c*alpha(i,k,j)
1269 end do
1270 end do
1271 ! recompute : gamma
1272 do i = i_start, i_end
1273 a_b = 0.
1274 a_c = 0.
1275 b = 1.+2.*cof(i)*rdnw(kde-1)**2*c2a(i,kde-1,j)
1276 ! recompute : b
1277 c = 0.
1278 ! recompute : c
1279 alpha(i,kde,j) = 1./(b-a(i,kde,j)*gamma(i,kde-1,j))
1280 ! recompute : alpha
1281 a_alpha(i,kde,j) = a_alpha(i,kde,j)+a_gamma(i,kde,j)*c
1282 a_c = a_c+a_gamma(i,kde,j)*alpha(i,kde,j)
1283 a_gamma(i,kde,j) = 0.
1284 a_a(i,kde,j) = a_a(i,kde,j)+a_alpha(i,kde,j)*(1.*gamma(i,kde-1,j)/((b-a(i,kde,j)*gamma(i,kde-1,j))*(b-a(i,kde,j)*gamma(i,kde-1,&
1285 &j))))
1286 a_b = a_b-a_alpha(i,kde,j)/((b-a(i,kde,j)*gamma(i,kde-1,j))*(b-a(i,kde,j)*gamma(i,kde-1,j)))
1287 a_gamma(i,kde-1,j) = a_gamma(i,kde-1,j)+a_alpha(i,kde,j)*(1.*a(i,kde,j)/((b-a(i,kde,j)*gamma(i,kde-1,j))*(b-a(i,kde,j)*gamma(i,&
1288 &kde-1,j))))
1289 a_alpha(i,kde,j) = 0.
1290 a_c = 0.
1291 a_c2a(i,kde-1,j) = a_c2a(i,kde-1,j)+2*a_b*cof(i)*rdnw(kde-1)**2
1292 a_cof(i) = a_cof(i)+2*a_b*rdnw(kde-1)**2*c2a(i,kde-1,j)
1293 a_b = 0.
1294 end do
1295 do k = kde-1, 2, -1
1296 ! recdepend vars : i_end,i_start,j
1297 ! recompute pos : DOLOOP_STMT module_small_step_em.f90:658
1298 ! recompute vars : gamma
1299 do i = i_start, i_end
1300 gamma(i,1,j) = 0.
1301 end do
1302 ! recompute vars : gamma
1303 do k1 = 2, k-1
1304 do i = i_start, i_end
1305 b = 1.+cqw(i,k1,j)*cof(i)*rdn(k1)*(rdnw(k1)*c2a(i,k1,j)+rdnw(k1-1)*c2a(i,k1-1,j))
1306 c = -(cqw(i,k1,j)*cof(i)*rdn(k1)*rdnw(k1)*c2a(i,k1,j))
1307 alpha(i,k1,j) = 1./(b-a(i,k1,j)*gamma(i,k1-1,j))
1308 gamma(i,k1,j) = c*alpha(i,k1,j)
1309 end do
1310 end do
1311 do i = i_start, i_end
1312 a_b = 0.
1313 a_c = 0.
1314 b = 1.+cqw(i,k,j)*cof(i)*rdn(k)*(rdnw(k)*c2a(i,k,j)+rdnw(k-1)*c2a(i,k-1,j))
1315 ! recompute : b
1316 c = -(cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k)*c2a(i,k,j))
1317 ! recompute : c
1318 alpha(i,k,j) = 1./(b-a(i,k,j)*gamma(i,k-1,j))
1319 ! recompute : alpha
1320 a_alpha(i,k,j) = a_alpha(i,k,j)+a_gamma(i,k,j)*c
1321 a_c = a_c+a_gamma(i,k,j)*alpha(i,k,j)
1322 a_gamma(i,k,j) = 0.
1323 a_a(i,k,j) = a_a(i,k,j)+a_alpha(i,k,j)*(1.*gamma(i,k-1,j)/((b-a(i,k,j)*gamma(i,k-1,j))*(b-a(i,k,j)*gamma(i,k-1,j))))
1324 a_b = a_b-a_alpha(i,k,j)/((b-a(i,k,j)*gamma(i,k-1,j))*(b-a(i,k,j)*gamma(i,k-1,j)))
1325 a_gamma(i,k-1,j) = a_gamma(i,k-1,j)+a_alpha(i,k,j)*(1.*a(i,k,j)/((b-a(i,k,j)*gamma(i,k-1,j))*(b-a(i,k,j)*gamma(i,k-1,j))))
1326 a_alpha(i,k,j) = 0.
1327 a_c2a(i,k,j) = a_c2a(i,k,j)-a_c*cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k)
1328 a_cof(i) = a_cof(i)-a_c*cqw(i,k,j)*rdn(k)*rdnw(k)*c2a(i,k,j)
1329 a_cqw(i,k,j) = a_cqw(i,k,j)-a_c*cof(i)*rdn(k)*rdnw(k)*c2a(i,k,j)
1330 a_c = 0.
1331 a_c2a(i,k-1,j) = a_c2a(i,k-1,j)+a_b*cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k-1)
1332 a_c2a(i,k,j) = a_c2a(i,k,j)+a_b*cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k)
1333 a_cof(i) = a_cof(i)+a_b*cqw(i,k,j)*rdn(k)*(rdnw(k)*c2a(i,k,j)+rdnw(k-1)*c2a(i,k-1,j))
1334 a_cqw(i,k,j) = a_cqw(i,k,j)+a_b*cof(i)*rdn(k)*(rdnw(k)*c2a(i,k,j)+rdnw(k-1)*c2a(i,k-1,j))
1335 a_b = 0.
1336 end do
1337 end do
1338 do k = 3, kde-1
1339 do i = i_start, i_end
1340 a_c2a(i,k-1,j) = a_c2a(i,k-1,j)-a_a(i,k,j)*cqw(i,k,j)*cof(i)*rdn(k)*rdnw(k-1)
1341 a_cof(i) = a_cof(i)-a_a(i,k,j)*cqw(i,k,j)*rdn(k)*rdnw(k-1)*c2a(i,k-1,j)
1342 a_cqw(i,k,j) = a_cqw(i,k,j)-a_a(i,k,j)*cof(i)*rdn(k)*rdnw(k-1)*c2a(i,k-1,j)
1343 a_a(i,k,j) = 0.
1344 end do
1345 end do
1346 do i = i_start, i_end
1347 cof(i) = (0.5*dts*g*(1.+epssm)/mut(i,j))**2
1348 ! recompute : cof
1349 a_gamma(i,1,j) = 0.
1350 a_c2a(i,kde-1,j) = a_c2a(i,kde-1,j)-2*a_a(i,kde,j)*cof(i)*rdnw(kde-1)**2
1351 a_cof(i) = a_cof(i)-2*a_a(i,kde,j)*rdnw(kde-1)**2*c2a(i,kde-1,j)
1352 a_a(i,kde,j) = 0.
1353 a_a(i,2,j) = 0.
1354 a_mut(i,j) = a_mut(i,j)-2*a_cof(i)*0.5*dts*g*(1.+epssm)/(mut(i,j)*mut(i,j))*(0.5*dts*g*(1.+epssm)/mut(i,j))
1355 a_cof(i) = 0.
1356 end do
1357 end do
1358
1359 end subroutine a_calc_coef_w
1360
1361
1362 subroutine a_calc_p_rho( al, a_al, p, a_p, ph, a_ph, alt, a_alt, t_2, a_t_2, t_1, a_t_1, c2a, a_c2a, a_pm1, mu, a_mu, muts, a_muts,&
1363 & znu, t0, rdnw, dnw, smdiv, non_hydrostatic, step, ide, jde, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte )
1364 !******************************************************************
1365 !******************************************************************
1366 !** This routine was generated by Automatic differentiation. **
1367 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
1368 !******************************************************************
1369 !******************************************************************
1370 !==============================================
1371 ! all entries are defined explicitly
1372 !==============================================
1373 implicit none
1374
1375 !==============================================
1376 ! declare arguments
1377 !==============================================
1378 integer, intent(in) :: ime
1379 integer, intent(in) :: ims
1380 integer, intent(in) :: jme
1381 integer, intent(in) :: jms
1382 integer, intent(in) :: kme
1383 integer, intent(in) :: kms
1384 real, intent(inout) :: a_al(ims:ime,kms:kme,jms:jme)
1385 real, intent(inout) :: a_alt(ims:ime,kms:kme,jms:jme)
1386 real, intent(inout) :: a_c2a(ims:ime,kms:kme,jms:jme)
1387 real, intent(inout) :: a_mu(ims:ime,jms:jme)
1388 real, intent(inout) :: a_muts(ims:ime,jms:jme)
1389 real, intent(inout) :: a_p(ims:ime,kms:kme,jms:jme)
1390 real, intent(inout) :: a_ph(ims:ime,kms:kme,jms:jme)
1391 real, intent(inout) :: a_pm1(ims:ime,kms:kme,jms:jme)
1392 real, intent(inout) :: a_t_1(ims:ime,kms:kme,jms:jme)
1393 real, intent(inout) :: a_t_2(ims:ime,kms:kme,jms:jme)
1394 real, intent(out) :: al(ims:ime,kms:kme,jms:jme)
1395 real, intent(in) :: alt(ims:ime,kms:kme,jms:jme)
1396 real, intent(in) :: c2a(ims:ime,kms:kme,jms:jme)
1397 real, intent(in) :: dnw(kms:kme)
1398 integer, intent(in) :: ide
1399 integer, intent(in) :: ite
1400 integer, intent(in) :: its
1401 integer, intent(in) :: jde
1402 integer, intent(in) :: jte
1403 integer, intent(in) :: jts
1404 integer, intent(in) :: kde
1405 integer, intent(in) :: kte
1406 integer, intent(in) :: kts
1407 real, intent(in) :: mu(ims:ime,jms:jme)
1408 real, intent(in) :: muts(ims:ime,jms:jme)
1409 logical, intent(in) :: non_hydrostatic
1410 real, intent(out) :: p(ims:ime,kms:kme,jms:jme)
1411 real, intent(inout) :: ph(ims:ime,kms:kme,jms:jme)
1412 real, intent(in) :: rdnw(kms:kme)
1413 real, intent(in) :: smdiv
1414 integer, intent(in) :: step
1415 real, intent(in) :: t0
1416 real, intent(in) :: t_1(ims:ime,kms:kme,jms:jme)
1417 real, intent(in) :: t_2(ims:ime,kms:kme,jms:jme)
1418 real, intent(in) :: znu(kms:kme)
1419
1420 !==============================================
1421 ! declare local variables
1422 !==============================================
1423 real a_ptmp
1424 integer i
1425 integer i_end
1426 integer i_start
1427 integer j
1428 integer j_end
1429 integer j_start
1430 integer k
1431 integer k_end
1432 integer k_start
1433
1434 !----------------------------------------------
1435 ! RESET LOCAL ADJOINT VARIABLES
1436 !----------------------------------------------
1437 a_ptmp = 0.
1438
1439 !----------------------------------------------
1440 ! ROUTINE BODY
1441 !----------------------------------------------
1442 i_start = its
1443 ! recompute : i_start
1444 i_end = ite
1445 ! recompute : i_end
1446 j_start = jts
1447 ! recompute : j_start
1448 j_end = jte
1449 ! recompute : j_end
1450 k_start = kts
1451 ! recompute : k_start
1452 k_end = min(kte,kde-1)
1453 ! recompute : k_end
1454 if (i_end .eq. ide) then
1455 i_end = i_end-1
1456 endif
1457 ! recompute : i_end
1458 if (j_end .eq. jde) then
1459 j_end = j_end-1
1460 endif
1461 ! recompute : j_end
1462 if (step .eq. 0) then
1463 do j = j_start, j_end
1464 do k = k_start, k_end
1465 do i = i_start, i_end
1466 a_p(i,k,j) = a_p(i,k,j)+a_pm1(i,k,j)
1467 a_pm1(i,k,j) = 0.
1468 end do
1469 end do
1470 end do
1471 else
1472 do j = j_start, j_end
1473 a_ptmp = 0.
1474 do k = k_start, k_end
1475 a_ptmp = 0.
1476 do i = i_start, i_end
1477 a_ptmp = 0.
1478 a_ptmp = a_ptmp+a_pm1(i,k,j)
1479 a_pm1(i,k,j) = 0.
1480 a_pm1(i,k,j) = a_pm1(i,k,j)-a_p(i,k,j)*smdiv
1481 a_p(i,k,j) = a_p(i,k,j)*(1+smdiv)
1482 a_p(i,k,j) = a_p(i,k,j)+a_ptmp
1483 a_ptmp = 0.
1484 end do
1485 end do
1486 end do
1487 endif
1488 if (non_hydrostatic) then
1489 do j = j_start, j_end
1490 do k = k_start, k_end
1491 do i = i_start, i_end
1492 al(i,k,j) = -(1./muts(i,j)*(alt(i,k,j)*mu(i,j)+rdnw(k)*(ph(i,k+1,j)-ph(i,k,j))))
1493 ! recompute : al
1494 a_al(i,k,j) = a_al(i,k,j)-a_p(i,k,j)*c2a(i,k,j)
1495 a_alt(i,k,j) = a_alt(i,k,j)+a_p(i,k,j)*c2a(i,k,j)*((t_2(i,k,j)-mu(i,j)*t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j))))
1496 a_c2a(i,k,j) = a_c2a(i,k,j)+a_p(i,k,j)*(alt(i,k,j)*(t_2(i,k,j)-mu(i,j)*t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j)))-al(i,k,j))
1497 a_mu(i,j) = a_mu(i,j)-a_p(i,k,j)*c2a(i,k,j)*(alt(i,k,j)*t_1(i,k,j)/(muts(i,j)*(t0+t_1(i,k,j))))
1498 a_muts(i,j) = a_muts(i,j)-a_p(i,k,j)*c2a(i,k,j)*(alt(i,k,j)*(t_2(i,k,j)-mu(i,j)*t_1(i,k,j))*(t0+t_1(i,k,j))/(muts(i,j)*(t0+&
1499 &t_1(i,k,j))*muts(i,j)*(t0+t_1(i,k,j))))
1500 a_t_1(i,k,j) = a_t_1(i,k,j)-a_p(i,k,j)*c2a(i,k,j)*(alt(i,k,j)*mu(i,j)/(muts(i,j)*(t0+t_1(i,k,j)))+alt(i,k,j)*(t_2(i,k,j)-&
1501 &mu(i,j)*t_1(i,k,j))*muts(i,j)/(muts(i,j)*(t0+t_1(i,k,j))*muts(i,j)*(t0+t_1(i,k,j))))
1502 a_t_2(i,k,j) = a_t_2(i,k,j)+a_p(i,k,j)*c2a(i,k,j)*(alt(i,k,j)/(muts(i,j)*(t0+t_1(i,k,j))))
1503 a_p(i,k,j) = 0.
1504 a_alt(i,k,j) = a_alt(i,k,j)-a_al(i,k,j)*1./muts(i,j)*mu(i,j)
1505 a_mu(i,j) = a_mu(i,j)-a_al(i,k,j)*1./muts(i,j)*alt(i,k,j)
1506 a_muts(i,j) = a_muts(i,j)+a_al(i,k,j)/(muts(i,j)*muts(i,j))*(alt(i,k,j)*mu(i,j)+rdnw(k)*(ph(i,k+1,j)-ph(i,k,j)))
1507 a_ph(i,k+1,j) = a_ph(i,k+1,j)-a_al(i,k,j)*1./muts(i,j)*rdnw(k)
1508 a_ph(i,k,j) = a_ph(i,k,j)+a_al(i,k,j)*1./muts(i,j)*rdnw(k)
1509 a_al(i,k,j) = 0.
1510 end do
1511 end do
1512 end do
1513 else
1514 do j = j_start, j_end
1515 do k = k_end, k_start, -1
1516 do i = i_start, i_end
1517 p(i,k,j) = mu(i,j)*znu(k)
1518 ! recompute : p
1519 al(i,k,j) = alt(i,k,j)*(t_2(i,k,j)-mu(i,j)*t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j)))-p(i,k,j)/c2a(i,k,j)
1520 ! recompute : al
1521 a_al(i,k,j) = a_al(i,k,j)-a_ph(i,k+1,j)*dnw(k)*muts(i,j)
1522 a_alt(i,k,j) = a_alt(i,k,j)-a_ph(i,k+1,j)*dnw(k)*mu(i,j)
1523 a_mu(i,j) = a_mu(i,j)-a_ph(i,k+1,j)*dnw(k)*alt(i,k,j)
1524 a_muts(i,j) = a_muts(i,j)-a_ph(i,k+1,j)*dnw(k)*al(i,k,j)
1525 a_ph(i,k,j) = a_ph(i,k,j)+a_ph(i,k+1,j)
1526 a_ph(i,k+1,j) = 0.
1527 a_alt(i,k,j) = a_alt(i,k,j)+a_al(i,k,j)*((t_2(i,k,j)-mu(i,j)*t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j))))
1528 a_c2a(i,k,j) = a_c2a(i,k,j)+a_al(i,k,j)*(p(i,k,j)/(c2a(i,k,j)*c2a(i,k,j)))
1529 a_mu(i,j) = a_mu(i,j)-a_al(i,k,j)*(alt(i,k,j)*t_1(i,k,j)/(muts(i,j)*(t0+t_1(i,k,j))))
1530 a_muts(i,j) = a_muts(i,j)-a_al(i,k,j)*(alt(i,k,j)*(t_2(i,k,j)-mu(i,j)*t_1(i,k,j))*(t0+t_1(i,k,j))/(muts(i,j)*(t0+t_1(i,k,j)&
1531 &)*muts(i,j)*(t0+t_1(i,k,j))))
1532 a_p(i,k,j) = a_p(i,k,j)-a_al(i,k,j)/c2a(i,k,j)
1533 a_t_1(i,k,j) = a_t_1(i,k,j)-a_al(i,k,j)*(alt(i,k,j)*mu(i,j)/(muts(i,j)*(t0+t_1(i,k,j)))+alt(i,k,j)*(t_2(i,k,j)-mu(i,j)*&
1534 &t_1(i,k,j))*muts(i,j)/(muts(i,j)*(t0+t_1(i,k,j))*muts(i,j)*(t0+t_1(i,k,j))))
1535 a_t_2(i,k,j) = a_t_2(i,k,j)+a_al(i,k,j)*(alt(i,k,j)/(muts(i,j)*(t0+t_1(i,k,j))))
1536 a_al(i,k,j) = 0.
1537 a_mu(i,j) = a_mu(i,j)+a_p(i,k,j)*znu(k)
1538 a_p(i,k,j) = 0.
1539 end do
1540 end do
1541 end do
1542 endif
1543
1544 end subroutine a_calc_p_rho
1545
1546 subroutine a_small_step_finish( u_2, a_u_2, v_2, a_v_2, w_2, a_w_2, t_2, a_t_2, a_ph_2, a_mu_2, mut, a_mut, muts, a_muts, muu, &
1547 &a_muu, muus, a_muus, muv, a_muv, muvs, a_muvs, u_save, a_u_save, v_save, a_v_save, w_save, a_w_save, t_save, a_t_save, a_ph_save, &
1548 &a_mu_save, msfu, msfv, msft, ide, jde, kds, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte )
1549
1550 !******************************************************************
1551 !******************************************************************
1552 !** This routine was generated by Automatic differentiation. **
1553 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
1554 !******************************************************************
1555 !******************************************************************
1556 !==============================================
1557 ! all entries are defined explicitly
1558 !==============================================
1559 implicit none
1560
1561 !==============================================
1562 ! declare arguments
1563 !==============================================
1564 integer, intent(in) :: ime
1565 integer, intent(in) :: ims
1566 integer, intent(in) :: jme
1567 integer, intent(in) :: jms
1568 real, intent(inout) :: a_mu_2(ims:ime,jms:jme)
1569 real, intent(inout) :: a_mu_save(ims:ime,jms:jme)
1570 real, intent(inout) :: a_mut(ims:ime,jms:jme)
1571 real, intent(inout) :: a_muts(ims:ime,jms:jme)
1572 real, intent(inout) :: a_muu(ims:ime,jms:jme)
1573 real, intent(inout) :: a_muus(ims:ime,jms:jme)
1574 real, intent(inout) :: a_muv(ims:ime,jms:jme)
1575 real, intent(inout) :: a_muvs(ims:ime,jms:jme)
1576 integer, intent(in) :: kme
1577 integer, intent(in) :: kms
1578 real, intent(inout) :: a_ph_2(ims:ime,kms:kme,jms:jme)
1579 real, intent(inout) :: a_ph_save(ims:ime,kms:kme,jms:jme)
1580 real, intent(inout) :: a_t_2(ims:ime,kms:kme,jms:jme)
1581 real, intent(inout) :: a_t_save(ims:ime,kms:kme,jms:jme)
1582 real, intent(inout) :: a_u_2(ims:ime,kms:kme,jms:jme)
1583 real, intent(inout) :: a_u_save(ims:ime,kms:kme,jms:jme)
1584 real, intent(inout) :: a_v_2(ims:ime,kms:kme,jms:jme)
1585 real, intent(inout) :: a_v_save(ims:ime,kms:kme,jms:jme)
1586 real, intent(inout) :: a_w_2(ims:ime,kms:kme,jms:jme)
1587 real, intent(inout) :: a_w_save(ims:ime,kms:kme,jms:jme)
1588 integer, intent(in) :: ide
1589 integer, intent(in) :: ite
1590 integer, intent(in) :: its
1591 integer, intent(in) :: jde
1592 integer, intent(in) :: jte
1593 integer, intent(in) :: jts
1594 integer, intent(in) :: kde
1595 integer, intent(in) :: kds
1596 real, intent(in) :: msft(ims:ime,jms:jme)
1597 real, intent(in) :: msfu(ims:ime,jms:jme)
1598 real, intent(in) :: msfv(ims:ime,jms:jme)
1599 real, intent(inout) :: mut(ims:ime,jms:jme)
1600 real, intent(inout) :: muts(ims:ime,jms:jme)
1601 real, intent(inout) :: muu(ims:ime,jms:jme)
1602 real, intent(inout) :: muus(ims:ime,jms:jme)
1603 real, intent(inout) :: muv(ims:ime,jms:jme)
1604 real, intent(inout) :: muvs(ims:ime,jms:jme)
1605 real, intent(inout) :: t_2(ims:ime,kms:kme,jms:jme)
1606 real, intent(in) :: t_save(ims:ime,kms:kme,jms:jme)
1607 real, intent(inout) :: u_2(ims:ime,kms:kme,jms:jme)
1608 real, intent(in) :: u_save(ims:ime,kms:kme,jms:jme)
1609 real, intent(inout) :: v_2(ims:ime,kms:kme,jms:jme)
1610 real, intent(in) :: v_save(ims:ime,kms:kme,jms:jme)
1611 real, intent(inout) :: w_2(ims:ime,kms:kme,jms:jme)
1612 real, intent(in) :: w_save(ims:ime,kms:kme,jms:jme)
1613
1614
1615 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) ::ww1,a_ww1,a_ww
1616 REAL, DIMENSION(ims:ime, kms:kme, jms:jme) :: h_diabatic
1617 INTEGER :: number_of_small_timesteps
1618 REAL :: dts
1619
1620
1621
1622 !==============================================
1623 ! declare local variables
1624 !==============================================
1625 integer i
1626 integer i_end
1627 integer i_endu
1628 integer i_start
1629 integer j
1630 integer j_end
1631 integer j_endv
1632 integer j_start
1633 integer k
1634
1635 !----------------------------------------------
1636 ! ROUTINE BODY
1637 !----------------------------------------------
1638 i_start = its
1639 ! recompute : i_start
1640 i_end = ite
1641 ! recompute : i_end
1642 j_start = jts
1643 ! recompute : j_start
1644 j_end = jte
1645 ! recompute : j_end
1646 i_endu = i_end
1647 ! recompute : i_endu
1648 j_endv = j_end
1649 ! recompute : j_endv
1650 if (i_end .eq. ide) then
1651 i_end = i_end-1
1652 endif
1653 ! recompute : i_end
1654 if (j_end .eq. jde) then
1655 j_end = j_end-1
1656 endif
1657 ! recompute : j_end
1658 do j = j_start, j_end
1659 do i = i_start, i_end
1660 a_mu_save(i,j) = a_mu_save(i,j)+a_mu_2(i,j)
1661 end do
1662 end do
1663 do j = j_start, j_end
1664 do k = kds, kde-1
1665 do i = i_start, i_end
1666 !zzma a_mut(i,j) = a_mut(i,j)+a_t_2(i,k,j)*(t_save(i,k,j)/muts(i,j))
1667 !zzma a_muts(i,j) = a_muts(i,j)-a_t_2(i,k,j)*((t_2(i,k,j)+t_save(i,k,j)*mut(i,j))/(muts(i,j)*muts(i,j)))
1668 !zzma a_t_save(i,k,j) = a_t_save(i,k,j)+a_t_2(i,k,j)*(mut(i,j)/muts(i,j))
1669 !zzma a_t_2(i,k,j) = a_t_2(i,k,j)/muts(i,j)
1670
1671 a_mut(i,j) = a_mut(i,j)+a_t_2(i,k,j)*(- dts*number_of_small_timesteps*h_diabatic(i,k,j) +t_save(i,k,j)/muts(i,j))
1672 a_muts(i,j) = a_muts(i,j)-a_t_2(i,k,j)*((t_2(i,k,j)- dts*number_of_small_timesteps*mut(i,j)*h_diabatic(i,k,j) &
1673 +t_save(i,k,j)*mut(i,j))/(muts(i,j)*muts(i,j)))
1674 a_t_save(i,k,j) = a_t_save(i,k,j)+a_t_2(i,k,j)*(mut(i,j)/muts(i,j))
1675 a_t_2(i,k,j) = a_t_2(i,k,j)/muts(i,j)
1676
1677 end do
1678 end do
1679 end do
1680 do j = j_start, j_end
1681 do k = kds, kde
1682 do i = i_start, i_end
1683 !------------------------- added in May 9 -----------------
1684 a_ww1(i,k,j) =a_ww1(i,k,j) +a_ww(i,k,j)
1685 !------------------------- added in May 9 -----------------
1686
1687 a_ph_save(i,k,j) = a_ph_save(i,k,j)+a_ph_2(i,k,j)
1688 a_mut(i,j) = a_mut(i,j)+a_w_2(i,k,j)*(w_save(i,k,j)/muts(i,j))
1689 a_muts(i,j) = a_muts(i,j)-a_w_2(i,k,j)*((msft(i,j)*w_2(i,k,j)+w_save(i,k,j)*mut(i,j))/(muts(i,j)*muts(i,j)))
1690 a_w_save(i,k,j) = a_w_save(i,k,j)+a_w_2(i,k,j)*(mut(i,j)/muts(i,j))
1691 a_w_2(i,k,j) = a_w_2(i,k,j)*(msft(i,j)/muts(i,j))
1692 end do
1693 end do
1694 end do
1695 do j = j_start, j_end
1696 do k = kds, kde-1
1697 do i = i_start, i_endu
1698 a_muu(i,j) = a_muu(i,j)+a_u_2(i,k,j)*(u_save(i,k,j)/muus(i,j))
1699 a_muus(i,j) = a_muus(i,j)-a_u_2(i,k,j)*((msfu(i,j)*u_2(i,k,j)+u_save(i,k,j)*muu(i,j))/(muus(i,j)*muus(i,j)))
1700 a_u_save(i,k,j) = a_u_save(i,k,j)+a_u_2(i,k,j)*(muu(i,j)/muus(i,j))
1701 a_u_2(i,k,j) = a_u_2(i,k,j)*(msfu(i,j)/muus(i,j))
1702 end do
1703 end do
1704 end do
1705 do j = j_start, j_endv
1706 do k = kds, kde-1
1707 do i = i_start, i_end
1708 a_muv(i,j) = a_muv(i,j)+a_v_2(i,k,j)*(v_save(i,k,j)/muvs(i,j))
1709 a_muvs(i,j) = a_muvs(i,j)-a_v_2(i,k,j)*((msfv(i,j)*v_2(i,k,j)+v_save(i,k,j)*muv(i,j))/(muvs(i,j)*muvs(i,j)))
1710 a_v_save(i,k,j) = a_v_save(i,k,j)+a_v_2(i,k,j)*(muv(i,j)/muvs(i,j))
1711 a_v_2(i,k,j) = a_v_2(i,k,j)*(msfv(i,j)/muvs(i,j))
1712 end do
1713 end do
1714 end do
1715
1716 end subroutine a_small_step_finish
1717
1718
1719 subroutine a_small_step_prep( u_1, a_u_1, u_2, a_u_2, v_1, a_v_1, v_2, a_v_2, w_1, a_w_1, w_2, a_w_2, t_1, a_t_1, t_2, a_t_2, &
1720 &a_ph_1, a_ph_2, mub, mu_1, a_mu_1, mu_2, a_mu_2, muu, a_muu, muus, a_muus, muv, a_muv, muvs, a_muvs, mut, a_mut, muts, a_muts, &
1721 &a_mudf, a_u_save, a_v_save, a_w_save, a_t_save, a_ph_save, a_mu_save, a_ww, a_ww_save, a_c2a, pb, p, a_p, alt, a_alt, msfu, msfv, &
1722 &msft, rk_step, ide, jde, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte )
1723 !******************************************************************
1724 !******************************************************************
1725 !** This routine was generated by Automatic differentiation. **
1726 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
1727 !******************************************************************
1728 !******************************************************************
1729 !==============================================
1730 ! all entries are defined explicitly
1731 !==============================================
1732 implicit none
1733
1734 !==============================================
1735 ! declare arguments
1736 !==============================================
1737 integer, intent(in) :: ime
1738 integer, intent(in) :: ims
1739 integer, intent(in) :: jme
1740 integer, intent(in) :: jms
1741 integer, intent(in) :: kme
1742 integer, intent(in) :: kms
1743 real, intent(inout) :: a_alt(ims:ime,kms:kme,jms:jme)
1744 real, intent(inout) :: a_c2a(ims:ime,kms:kme,jms:jme)
1745 real, intent(inout) :: a_mu_1(ims:ime,jms:jme)
1746 real, intent(inout) :: a_mu_2(ims:ime,jms:jme)
1747 real, intent(inout) :: a_mu_save(ims:ime,jms:jme)
1748 real, intent(inout) :: a_mudf(ims:ime,jms:jme)
1749 real, intent(inout) :: a_mut(ims:ime,jms:jme)
1750 real, intent(inout) :: a_muts(ims:ime,jms:jme)
1751 real, intent(inout) :: a_muu(ims:ime,jms:jme)
1752 real, intent(inout) :: a_muus(ims:ime,jms:jme)
1753 real, intent(inout) :: a_muv(ims:ime,jms:jme)
1754 real, intent(inout) :: a_muvs(ims:ime,jms:jme)
1755 real, intent(inout) :: a_p(ims:ime,kms:kme,jms:jme)
1756 real, intent(inout) :: a_ph_1(ims:ime,kms:kme,jms:jme)
1757 real, intent(inout) :: a_ph_2(ims:ime,kms:kme,jms:jme)
1758 real, intent(inout) :: a_ph_save(ims:ime,kms:kme,jms:jme)
1759 real, intent(inout) :: a_t_1(ims:ime,kms:kme,jms:jme)
1760 real, intent(inout) :: a_t_2(ims:ime,kms:kme,jms:jme)
1761 real, intent(inout) :: a_t_save(ims:ime,kms:kme,jms:jme)
1762 real, intent(inout) :: a_u_1(ims:ime,kms:kme,jms:jme)
1763 real, intent(inout) :: a_u_2(ims:ime,kms:kme,jms:jme)
1764 real, intent(inout) :: a_u_save(ims:ime,kms:kme,jms:jme)
1765 real, intent(inout) :: a_v_1(ims:ime,kms:kme,jms:jme)
1766 real, intent(inout) :: a_v_2(ims:ime,kms:kme,jms:jme)
1767 real, intent(inout) :: a_v_save(ims:ime,kms:kme,jms:jme)
1768 real, intent(inout) :: a_w_1(ims:ime,kms:kme,jms:jme)
1769 real, intent(inout) :: a_w_2(ims:ime,kms:kme,jms:jme)
1770 real, intent(inout) :: a_w_save(ims:ime,kms:kme,jms:jme)
1771 real, intent(inout) :: a_ww(ims:ime,kms:kme,jms:jme)
1772 real, intent(inout) :: a_ww_save(ims:ime,kms:kme,jms:jme)
1773 real, intent(in) :: alt(ims:ime,kms:kme,jms:jme)
1774 integer, intent(in) :: ide
1775 integer, intent(in) :: ite
1776 integer, intent(in) :: its
1777 integer, intent(in) :: jde
1778 integer, intent(in) :: jte
1779 integer, intent(in) :: jts
1780 integer, intent(in) :: kde
1781 integer, intent(in) :: kte
1782 integer, intent(in) :: kts
1783 real, intent(inout) :: msft(ims:ime,jms:jme)
1784 real, intent(inout) :: msfu(ims:ime,jms:jme)
1785 real, intent(inout) :: msfv(ims:ime,jms:jme)
1786 real, intent(inout) :: mu_1(ims:ime,jms:jme)
1787 real, intent(inout) :: mu_2(ims:ime,jms:jme)
1788 real, intent(inout) :: mub(ims:ime,jms:jme)
1789 real, intent(inout) :: mut(ims:ime,jms:jme)
1790 real, intent(out) :: muts(ims:ime,jms:jme)
1791 real, intent(inout) :: muu(ims:ime,jms:jme)
1792 real, intent(out) :: muus(ims:ime,jms:jme)
1793 real, intent(inout) :: muv(ims:ime,jms:jme)
1794 real, intent(out) :: muvs(ims:ime,jms:jme)
1795 real, intent(in) :: p(ims:ime,kms:kme,jms:jme)
1796 real, intent(in) :: pb(ims:ime,kms:kme,jms:jme)
1797 integer, intent(in) :: rk_step
1798 real, intent(inout) :: t_1(ims:ime,kms:kme,jms:jme)
1799 real, intent(inout) :: t_2(ims:ime,kms:kme,jms:jme)
1800 real, intent(inout) :: u_1(ims:ime,kms:kme,jms:jme)
1801 real, intent(inout) :: u_2(ims:ime,kms:kme,jms:jme)
1802 real, intent(inout) :: v_1(ims:ime,kms:kme,jms:jme)
1803 real, intent(inout) :: v_2(ims:ime,kms:kme,jms:jme)
1804 real, intent(inout) :: w_1(ims:ime,kms:kme,jms:jme)
1805 real, intent(inout) :: w_2(ims:ime,kms:kme,jms:jme)
1806
1807 !==============================================
1808 ! declare local variables
1809 !==============================================
1810 integer i
1811 integer i_end
1812 integer i_endu
1813 integer i_start
1814 integer j
1815 integer j_end
1816 integer j_endv
1817 integer j_start
1818 integer k
1819 integer k_end
1820 integer k_start
1821
1822 !----------------------------------------------
1823 ! ROUTINE BODY
1824 !----------------------------------------------
1825 i_start = its
1826 ! recompute : i_start
1827 i_end = ite
1828 ! recompute : i_end
1829 j_start = jts
1830 ! recompute : j_start
1831 j_end = jte
1832 ! recompute : j_end
1833 k_start = kts
1834 ! recompute : k_start
1835 k_end = min(kte,kde-1)
1836 ! recompute : k_end
1837 i_endu = i_end
1838 ! recompute : i_endu
1839 j_endv = j_end
1840 ! recompute : j_endv
1841 if (i_end .eq. ide) then
1842 i_end = i_end-1
1843 endif
1844 ! recompute : i_end
1845 if (j_end .eq. jde) then
1846 j_end = j_end-1
1847 endif
1848
1849 ! recompute : j_end
1850 if (rk_step .eq. 1) then
1851
1852 !--------- added by zzma -----------
1853 DO j=j_start, j_end
1854 DO i=i_start, i_end
1855 mu_1(i,j)=mu_2(i,j)
1856 ENDDO
1857 ENDDO
1858 !--------- added by zzma -----------
1859
1860 do j = j_start, j_end
1861 do k = k_start, k_end
1862 do i = i_start, i_endu
1863 u_1(i,k,j) = u_2(i,k,j)
1864 end do
1865 end do
1866 end do
1867 do j = j_start, j_endv
1868 do k = k_start, k_end
1869 do i = i_start, i_end
1870 v_1(i,k,j) = v_2(i,k,j)
1871 end do
1872 end do
1873 end do
1874 do j = j_start, j_end
1875 do k = k_start, k_end
1876 do i = i_start, i_end
1877 t_1(i,k,j) = t_2(i,k,j)
1878 end do
1879 end do
1880 end do
1881 do j = j_start, j_end
1882 do k = k_start, min(kde,kte)
1883 do i = i_start, i_end
1884 w_1(i,k,j) = w_2(i,k,j)
1885
1886 !--------- added by zzma -----------
1887 ! ph_1(i,k,j) = ph_2(i,k,j)
1888 !--------- added by zzma -----------
1889
1890 end do
1891 end do
1892 end do
1893 do j = j_start, j_end
1894 do i = i_start, i_end
1895 muts(i,j) = mub(i,j)+mu_2(i,j)
1896 end do
1897 do i = i_start, i_endu
1898 muus(i,j) = muu(i,j)
1899 end do
1900 end do
1901 do j = j_start, j_endv
1902 do i = i_start, i_end
1903 muvs(i,j) = muv(i,j)
1904 end do
1905 end do
1906 else
1907 do j = j_start, j_end
1908 do i = i_start, i_end
1909 muts(i,j) = mub(i,j)+mu_1(i,j)
1910 end do
1911 do i = i_start, i_endu
1912 muus(i,j) = 0.5*(mub(i,j)+mu_1(i,j)+mub(i-1,j)+mu_1(i-1,j))
1913 end do
1914 end do
1915 do j = j_start, j_endv
1916 do i = i_start, i_end
1917 muvs(i,j) = 0.5*(mub(i,j)+mu_1(i,j)+mub(i,j-1)+mu_1(i,j-1))
1918 end do
1919 end do
1920
1921 !--------- added by zzma -----------
1922 DO j=j_start, j_end
1923 DO i=i_start, i_end
1924 mu_2(i,j)=mu_1(i,j)-mu_2(i,j)
1925 ENDDO
1926 ENDDO
1927 !--------- added by zzma -----------
1928
1929 endif
1930 ! recompute : muts,muus,muvs,t_1,u_1,v_1,w_1
1931 do j = j_start, j_end
1932 do k = k_start, kde
1933 do i = i_start, i_end
1934 a_ww(i,k,j) = a_ww(i,k,j)+a_ww_save(i,k,j)
1935 a_ww_save(i,k,j) = 0.
1936 end do
1937 end do
1938 end do
1939 do j = j_start, j_end
1940 do k = k_start, kde
1941 do i = i_start, i_end
1942 a_ph_1(i,k,j) = a_ph_1(i,k,j)+a_ph_2(i,k,j)
1943 a_ph_2(i,k,j) = -a_ph_2(i,k,j)
1944 a_ph_2(i,k,j) = a_ph_2(i,k,j)+a_ph_save(i,k,j)
1945 a_ph_save(i,k,j) = 0.
1946 a_mut(i,j) = a_mut(i,j)-a_w_2(i,k,j)*(w_2(i,k,j)/msft(i,j))
1947 a_muts(i,j) = a_muts(i,j)+a_w_2(i,k,j)*(w_1(i,k,j)/msft(i,j))
1948 a_w_1(i,k,j) = a_w_1(i,k,j)+a_w_2(i,k,j)*(muts(i,j)/msft(i,j))
1949 a_w_2(i,k,j) = -(a_w_2(i,k,j)*(mut(i,j)/msft(i,j)))
1950 a_w_2(i,k,j) = a_w_2(i,k,j)+a_w_save(i,k,j)
1951 a_w_save(i,k,j) = 0.
1952 end do
1953 end do
1954 end do
1955 do j = j_start, j_end
1956 do k = k_start, k_end
1957 do i = i_start, i_end
1958 a_mut(i,j) = a_mut(i,j)-a_t_2(i,k,j)*t_2(i,k,j)
1959 a_muts(i,j) = a_muts(i,j)+a_t_2(i,k,j)*t_1(i,k,j)
1960 a_t_1(i,k,j) = a_t_1(i,k,j)+a_t_2(i,k,j)*muts(i,j)
1961 a_t_2(i,k,j) = -(a_t_2(i,k,j)*mut(i,j))
1962 a_t_2(i,k,j) = a_t_2(i,k,j)+a_t_save(i,k,j)
1963 a_t_save(i,k,j) = 0.
1964 end do
1965 end do
1966 end do
1967 do j = j_start, j_endv
1968 do k = k_start, k_end
1969 do i = i_start, i_end
1970 a_muv(i,j) = a_muv(i,j)-a_v_2(i,k,j)*(v_2(i,k,j)/msfv(i,j))
1971 a_muvs(i,j) = a_muvs(i,j)+a_v_2(i,k,j)*(v_1(i,k,j)/msfv(i,j))
1972 a_v_1(i,k,j) = a_v_1(i,k,j)+a_v_2(i,k,j)*(muvs(i,j)/msfv(i,j))
1973 a_v_2(i,k,j) = -(a_v_2(i,k,j)*(muv(i,j)/msfv(i,j)))
1974 a_v_2(i,k,j) = a_v_2(i,k,j)+a_v_save(i,k,j)
1975 a_v_save(i,k,j) = 0.
1976 end do
1977 end do
1978 end do
1979 do j = j_start, j_end
1980 do k = k_start, k_end
1981 do i = i_start, i_endu
1982 a_muu(i,j) = a_muu(i,j)-a_u_2(i,k,j)*(u_2(i,k,j)/msfu(i,j))
1983 a_muus(i,j) = a_muus(i,j)+a_u_2(i,k,j)*(u_1(i,k,j)/msfu(i,j))
1984 a_u_1(i,k,j) = a_u_1(i,k,j)+a_u_2(i,k,j)*(muus(i,j)/msfu(i,j))
1985 a_u_2(i,k,j) = -(a_u_2(i,k,j)*(muu(i,j)/msfu(i,j)))
1986 a_u_2(i,k,j) = a_u_2(i,k,j)+a_u_save(i,k,j)
1987 a_u_save(i,k,j) = 0.
1988 end do
1989 end do
1990 end do
1991 do j = j_start, j_end
1992 do k = k_start, k_end
1993 do i = i_start, i_end
1994 a_alt(i,k,j) = a_alt(i,k,j)-a_c2a(i,k,j)*(cpovcv*(pb(i,k,j)+p(i,k,j))/(alt(i,k,j)*alt(i,k,j)))
1995 a_p(i,k,j) = a_p(i,k,j)+a_c2a(i,k,j)*(cpovcv/alt(i,k,j))
1996 a_c2a(i,k,j) = 0.
1997 end do
1998 end do
1999 end do
2000 do j = j_start, j_end
2001 do i = i_start, i_end
2002 a_ww_save(i,1,j) = 0.
2003 a_ww_save(i,kde,j) = 0.
2004 end do
2005 end do
2006 if (rk_step .eq. 1) then
2007 do j = j_start, j_end
2008 do i = i_start, i_end
2009 a_mu_2(i,j) = 0.
2010 a_mu_2(i,j) = a_mu_2(i,j)+a_mu_save(i,j)
2011 a_mu_save(i,j) = 0.
2012 end do
2013 end do
2014 do j = j_start, j_endv
2015 do i = i_start, i_end
2016 a_muv(i,j) = a_muv(i,j)+a_muvs(i,j)
2017 a_muvs(i,j) = 0.
2018 end do
2019 end do
2020 do j = j_start, j_end
2021 do i = i_start, i_endu
2022 a_muu(i,j) = a_muu(i,j)+a_muus(i,j)
2023 a_muus(i,j) = 0.
2024 end do
2025 do i = i_start, i_end
2026 a_mu_2(i,j) = a_mu_2(i,j)+a_muts(i,j)
2027 a_muts(i,j) = 0.
2028 end do
2029 end do
2030 do j = j_start, j_end
2031 do k = k_start, min(kde,kte)
2032 do i = i_start, i_end
2033 a_ph_2(i,k,j) = a_ph_2(i,k,j)+a_ph_1(i,k,j)
2034 a_ph_1(i,k,j) = 0.
2035 a_w_2(i,k,j) = a_w_2(i,k,j)+a_w_1(i,k,j)
2036 a_w_1(i,k,j) = 0.
2037 end do
2038 end do
2039 end do
2040 do j = j_start, j_end
2041 do k = k_start, k_end
2042 do i = i_start, i_end
2043 a_t_2(i,k,j) = a_t_2(i,k,j)+a_t_1(i,k,j)
2044 a_t_1(i,k,j) = 0.
2045 end do
2046 end do
2047 end do
2048 do j = j_start, j_endv
2049 do k = k_start, k_end
2050 do i = i_start, i_end
2051 a_v_2(i,k,j) = a_v_2(i,k,j)+a_v_1(i,k,j)
2052 a_v_1(i,k,j) = 0.
2053 end do
2054 end do
2055 end do
2056 do j = j_start, j_end
2057 do k = k_start, k_end
2058 do i = i_start, i_endu
2059 a_u_2(i,k,j) = a_u_2(i,k,j)+a_u_1(i,k,j)
2060 a_u_1(i,k,j) = 0.
2061 end do
2062 end do
2063 end do
2064 do j = j_start, j_end
2065 do i = i_start, i_end
2066 a_mudf(i,j) = 0.
2067 a_ww_save(i,1,j) = 0.
2068 a_ww_save(i,kde,j) = 0.
2069 a_mu_2(i,j) = a_mu_2(i,j)+a_mu_1(i,j)
2070 a_mu_1(i,j) = 0.
2071 end do
2072 end do
2073 else
2074 !---------- zzma -----------
2075 do j = j_start, j_end
2076 do i = i_start, i_end
2077 ! a_mu_1(i,j) = 0.0
2078 end do
2079 end do
2080 do j = j_start, j_endv
2081 do i = i_start, i_end
2082 ! a_mu_1(i,j) = 0.0
2083 end do
2084 end do
2085 do j = j_start, j_endv
2086 do i = i_start, i_endu
2087 ! a_mu_1(i,j) = 0.0
2088 end do
2089 end do
2090
2091 ! a_mu_1=0
2092 ! a_u_1=0
2093 ! a_v_1=0
2094 ! a_w_1=0
2095 ! a_t_1=0
2096 ! a_ph_1=0
2097 ! a_muts =0
2098 ! a_muvs=0
2099 ! a_muus=0
2100 ! a_mu_save=0
2101
2102 !---------- zzma -----------
2103 do j = j_start, j_end
2104 do i = i_start, i_end
2105 a_mu_1(i,j) = a_mu_1(i,j)+a_mu_2(i,j)
2106 a_mu_2(i,j) = -a_mu_2(i,j)
2107 a_mu_2(i,j) = a_mu_2(i,j)+a_mu_save(i,j)
2108 a_mu_save(i,j) = 0.
2109 end do
2110 end do
2111 do j = j_start, j_endv
2112 do i = i_start, i_end
2113 a_mu_1(i,j-1) = a_mu_1(i,j-1)+0.5*a_muvs(i,j)
2114 a_mu_1(i,j) = a_mu_1(i,j)+0.5*a_muvs(i,j)
2115 a_muvs(i,j) = 0.
2116 end do
2117 end do
2118 do j = j_start, j_end
2119 do i = i_start, i_endu
2120 a_mu_1(i-1,j) = a_mu_1(i-1,j)+0.5*a_muus(i,j)
2121 a_mu_1(i,j) = a_mu_1(i,j)+0.5*a_muus(i,j)
2122 a_muus(i,j) = 0.
2123 end do
2124 do i = i_start, i_end
2125 a_mu_1(i,j) = a_mu_1(i,j)+a_muts(i,j)
2126 a_muts(i,j) = 0.
2127 end do
2128 end do
2129 endif
2130
2131 end subroutine a_small_step_prep
2132
2133
2134 subroutine a_sumflux( a_ru, a_rv, a_ww, u_lin, a_u_lin, v_lin, a_v_lin, a_ww_lin, muu, a_muu, muv, a_muv, a_ru_m, a_rv_m, a_ww_m, &
2135 &msfu, msfv, iteration, number_of_small_timesteps, ide, jde, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte )
2136 !******************************************************************
2137 !******************************************************************
2138 !** This routine was generated by Automatic differentiation. **
2139 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
2140 !******************************************************************
2141 !******************************************************************
2142 !==============================================
2143 ! all entries are defined explicitly
2144 !==============================================
2145 implicit none
2146
2147 !==============================================
2148 ! declare arguments
2149 !==============================================
2150 integer, intent(in) :: ime
2151 integer, intent(in) :: ims
2152 integer, intent(in) :: jme
2153 integer, intent(in) :: jms
2154 real, intent(inout) :: a_muu(ims:ime,jms:jme)
2155 real, intent(inout) :: a_muv(ims:ime,jms:jme)
2156 integer, intent(in) :: kme
2157 integer, intent(in) :: kms
2158 real, intent(inout) :: a_ru(ims:ime,kms:kme,jms:jme)
2159 real, intent(inout) :: a_ru_m(ims:ime,kms:kme,jms:jme)
2160 real, intent(inout) :: a_rv(ims:ime,kms:kme,jms:jme)
2161 real, intent(inout) :: a_rv_m(ims:ime,kms:kme,jms:jme)
2162 real, intent(inout) :: a_u_lin(ims:ime,kms:kme,jms:jme)
2163 real, intent(inout) :: a_v_lin(ims:ime,kms:kme,jms:jme)
2164 real, intent(inout) :: a_ww(ims:ime,kms:kme,jms:jme)
2165 real, intent(inout) :: a_ww_lin(ims:ime,kms:kme,jms:jme)
2166 real, intent(inout) :: a_ww_m(ims:ime,kms:kme,jms:jme)
2167 integer, intent(in) :: ide
2168 integer, intent(in) :: ite
2169 integer, intent(in) :: iteration
2170 integer, intent(in) :: its
2171 integer, intent(in) :: jde
2172 integer, intent(in) :: jte
2173 integer, intent(in) :: jts
2174 integer, intent(in) :: kde
2175 integer, intent(in) :: kte
2176 integer, intent(in) :: kts
2177 real, intent(in) :: msfu(ims:ime,jms:jme)
2178 real, intent(in) :: msfv(ims:ime,jms:jme)
2179 real, intent(in) :: muu(ims:ime,jms:jme)
2180 real, intent(in) :: muv(ims:ime,jms:jme)
2181 integer, intent(in) :: number_of_small_timesteps
2182 real, intent(in) :: u_lin(ims:ime,kms:kme,jms:jme)
2183 real, intent(in) :: v_lin(ims:ime,kms:kme,jms:jme)
2184
2185 !==============================================
2186 ! declare local variables
2187 !==============================================
2188 integer i
2189 integer j
2190 integer k
2191
2192 !----------------------------------------------
2193 ! ROUTINE BODY
2194 !----------------------------------------------
2195 if (iteration .eq. number_of_small_timesteps) then
2196 do j = jts, min(jde-1,jte)
2197 do k = kts, kte
2198 do i = its, min(ide-1,ite)
2199 a_ww_lin(i,k,j) = a_ww_lin(i,k,j)+a_ww_m(i,k,j)
2200 a_ww_m(i,k,j) = a_ww_m(i,k,j)/float(number_of_small_timesteps)
2201 end do
2202 end do
2203 end do
2204 do j = jts, jte
2205 do k = kts, min(kde-1,kte)
2206 do i = its, min(ide-1,ite)
2207 a_muv(i,j) = a_muv(i,j)+a_rv_m(i,k,j)*(v_lin(i,k,j)/msfv(i,j))
2208 a_v_lin(i,k,j) = a_v_lin(i,k,j)+a_rv_m(i,k,j)*(muv(i,j)/msfv(i,j))
2209 a_rv_m(i,k,j) = a_rv_m(i,k,j)/float(number_of_small_timesteps)
2210 end do
2211 end do
2212 end do
2213 do j = jts, min(jde-1,jte)
2214 do k = kts, min(kde-1,kte)
2215 do i = its, ite
2216 a_muu(i,j) = a_muu(i,j)+a_ru_m(i,k,j)*(u_lin(i,k,j)/msfu(i,j))
2217 a_u_lin(i,k,j) = a_u_lin(i,k,j)+a_ru_m(i,k,j)*(muu(i,j)/msfu(i,j))
2218 a_ru_m(i,k,j) = a_ru_m(i,k,j)/float(number_of_small_timesteps)
2219 end do
2220 end do
2221 end do
2222 endif
2223 do j = jts, min(jde-1,jte)
2224 do k = kts, kte
2225 do i = its, min(ide-1,ite)
2226 a_ww(i,k,j) = a_ww(i,k,j)+a_ww_m(i,k,j)
2227 end do
2228 end do
2229 end do
2230 do j = jts, jte
2231 do k = kts, min(kde-1,kte)
2232 do i = its, min(ide-1,ite)
2233 a_rv(i,k,j) = a_rv(i,k,j)+a_rv_m(i,k,j)
2234 end do
2235 end do
2236 end do
2237 do j = jts, min(jde-1,jte)
2238 do k = kts, min(kde-1,kte)
2239 do i = its, ite
2240 a_ru(i,k,j) = a_ru(i,k,j)+a_ru_m(i,k,j)
2241 end do
2242 end do
2243 end do
2244 if (iteration .eq. 1) then
2245 do j = jts, jte
2246 do k = kts, kte
2247 do i = its, ite
2248 a_ww_m(i,k,j) = 0.
2249 a_rv_m(i,k,j) = 0.
2250 a_ru_m(i,k,j) = 0.
2251 end do
2252 end do
2253 end do
2254 endif
2255
2256 end subroutine a_sumflux
2257
2258
2259 end module a_module_small_step_em
2260
2261