solve_em_ad.F
References to this file elsewhere.
1 ! DISCLAIMER
2 !
3 ! This file was generated by TAF version 1.7.18
4 !***********************************************************
5 !***********************************************************
6
7 ! revised by Qingnong Xiao and Zaizhong Ma in June 08,2005
8
9 !***********************************************************
10 !***********************************************************
11 !
12 ! FASTOPT DISCLAIMS ALL WARRANTIES, EXPRESS OR IMPLIED,
13 ! INCLUDING (WITHOUT LIMITATION) ALL IMPLIED WARRANTIES OF
14 ! MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, WITH
15 ! RESPECT TO THE SOFTWARE AND USER PROGRAMS. IN NO EVENT
16 ! SHALL FASTOPT BE LIABLE FOR ANY LOST OR ANTICIPATED PROF-
17 ! ITS, OR ANY INDIRECT, INCIDENTAL, EXEMPLARY, SPECIAL, OR
18 ! CONSEQUENTIAL DAMAGES, WHETHER OR NOT FASTOPT WAS ADVISED
19 ! OF THE POSSIBILITY OF SUCH DAMAGES.
20 !
21 ! Haftungsbeschraenkung
22 ! FastOpt gibt ausdruecklich keine Gewaehr, explizit oder indirekt,
23 ! bezueglich der Brauchbarkeit der Software fuer einen bestimmten
24 ! Zweck. Unter keinen Umstaenden ist FastOpt haftbar fuer
25 ! irgendeinen Verlust oder nicht eintretenden erwarteten Gewinn und
26 ! allen indirekten, zufaelligen, exemplarischen oder speziellen
27 ! Schaeden oder Folgeschaeden unabhaengig von einer eventuellen
28 ! Mitteilung darueber an FastOpt.
29 !
30 !******************************************************************
31 !******************************************************************
32 !** This routine was generated by Automatic differentiation. **
33 !** FastOpt: Transformation of Algorithm in Fortran, TAF 1.7.18 **
34 !******************************************************************
35 !******************************************************************
36 !WRF:MEDIATION_LAYER:SOLVER
37
38 SUBROUTINE solve_em_ad ( grid , config_flags , &
39 !
40 #include "em_dummy_args.inc"
41 !
42 )
43
44 ! Driver layer modules
45 USE module_domain
46 USE module_configure
47 USE module_driver_constants
48 USE module_machine
49 USE module_tiles
50 USE module_dm
51 ! Mediation layer modules
52 ! Model layer modules
53 USE module_model_constants
54 USE module_small_step_em
55 USE module_em
56 USE module_big_step_utilities_em
57 USE module_bc
58 USE module_bc_em
59 USE module_solvedebug_em
60 USE module_physics_addtendc
61 USE module_diffusion_em
62 ! Registry generated module
63 USE module_state_description
64 USE module_radiation_driver
65 USE module_surface_driver
66 USE module_cumulus_driver
67 USE module_microphysics_driver
68 USE module_pbl_driver
69 #ifdef WRF_CHEM
70 USE module_input_chem_data
71 #endif
72
73 USE a_module_small_step_em
74 USE a_module_em
75 USE a_module_big_step_utilities_em
76 USE a_module_bc
77 USE a_module_bc_em
78 USE a_module_diffusion_em
79
80 IMPLICIT NONE
81
82 ! Input data.
83
84 TYPE(domain) , TARGET :: grid
85
86 ! Definitions of dummy arguments to solve
87 #include <em_dummy_decl.inc>
88
89 ! WRF state bcs
90 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags
91
92 ! WRF state data
93
94 !#include "../phys/physics_drive.int"
95
96 ! Local data
97
98 INTEGER :: k_start , k_end, its, ite, jts, jte
99 INTEGER :: ids , ide , jds , jde , kds , kde , &
100 ims , ime , jms , jme , kms , kme , &
101 ips , ipe , jps , jpe , kps , kpe
102 INTEGER :: ij , iteration
103 INTEGER :: im , num_3d_m , ic , num_3d_c
104 INTEGER :: loop
105 INTEGER :: ijds, ijde
106 INTEGER :: itmpstep
107 INTEGER :: sz
108
109 ! storage for tendencies and decoupled state (generated from Registry)
110 #include <em_i1_decl.inc>
111
112 #ifdef WRFVAR
113
114 INTEGER :: rc
115 INTEGER :: number_of_small_timesteps, rk_step
116 INTEGER :: klevel,ijm,ijp,i,j,k,size1,size2 ! for prints/plots only
117 INTEGER :: idum1, idum2, dynamics_option
118
119 INTEGER :: rk_order, iwmax, jwmax, kwmax
120 REAL :: dt_rk, dts_rk, dtm, wmax
121 INTEGER :: l,kte,kk
122
123 real a_z_at_wh(1+grid%em31-grid%sm31,1+grid%em32-grid%sm32,1+grid%em33-grid%sm33)
124 real ah(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
125 real alh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
126 real ali(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
127 real alj(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
128 real alphah(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
129 real dnwh(grid%sm32:grid%em32)
130 real dt_rkh
131 real dts_rkh
132 real gammah(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
133 integer iteration1
134 real moist_h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist)
135 real moist_i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist)
136 real moist_tendh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist)
137 real moist_tendi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist)
138 real mu_1h(grid%sm31:grid%em31,grid%sm33:grid%em33)
139 real mu_1i(grid%sm31:grid%em31,grid%sm33:grid%em33)
140 real mu_2h(grid%sm31:grid%em31,grid%sm33:grid%em33)
141 real mu_2i(grid%sm31:grid%em31,grid%sm33:grid%em33)
142 real mu_2j(grid%sm31:grid%em31,grid%sm33:grid%em33)
143 integer number_of_small_timestepsh
144 real ph(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
145 real ph_1h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
146 real ph_1i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
147 real ph_2h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
148 real ph_2i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
149 real ph_2j(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
150 real ph_tendfh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
151 real ph_tendfi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
152 real pi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
153 real pj(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
154 real pm1h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
155 real pm1i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
156 integer rk_step1
157 real ru_mh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
158 real ru_tendfh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
159 real ru_tendfi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
160 real rv_mh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
161 real rv_tendfh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
162 real rv_tendfi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
163 real rw_tendfh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
164 real rw_tendfi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
165 real t_1h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
166 real t_1i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
167 real t_2h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
168 real t_2i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
169 real t_2j(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
170 real t_tendfh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
171 real t_tendfi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
172 real u_1h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
173 real u_1i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
174 real u_2h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
175 real u_2i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
176 real v_1h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
177 real v_1i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
178 real v_2h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
179 real v_2i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
180 real w_1h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
181 real w_1i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
182 real w_2h(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
183 real w_2i(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
184 real w_2j(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
185 real ww_mh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
186 real wwh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
187 real xkmhdh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
188 real xkmhdi(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
189
190 real a_pi_phyh(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
191
192 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: a_bn2h
193
194 real, DIMENSION(max(grid%ed31,grid%ed33),grid%sd32:grid%ed32,grid%spec_bdy_width,4) :: a_rqc_bth
195 real, DIMENSION(max(grid%ed31,grid%ed33),grid%sd32:grid%ed32,grid%spec_bdy_width,4) :: a_rqr_bth
196 real, DIMENSION(max(grid%ed31,grid%ed33),grid%sd32:grid%ed32,grid%spec_bdy_width,4) :: a_rqi_bth
197 real, DIMENSION(max(grid%ed31,grid%ed33),grid%sd32:grid%ed32,grid%spec_bdy_width,4) :: a_rqs_bth
198 real, DIMENSION(max(grid%ed31,grid%ed33),grid%sd32:grid%ed32,grid%spec_bdy_width,4) :: a_rqg_bth
199
200 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist) :: a_moist_tend
201
202 !----------------------------------------------------------------------------------------
203 !defined by Zaizhong Ma for saving the basic states
204 !----------------------------------------------------------------------------------------
205 !----------------------------------------------------------------------------------------
206 !*********************************************************************
207 real u_1o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
208 real u_2o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
209 real v_1o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
210 real v_2o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
211 real w_1o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
212 real w_2o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
213 real t_1o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
214 real t_2o(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
215 real mu_1o(grid%sm31:grid%em31,grid%sm33:grid%em33)
216 real mu_2o(grid%sm31:grid%em31,grid%sm33:grid%em33)
217 real po(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
218
219
220 real alp(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
221 real pp(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
222 real ph_2p(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
223 real t_2p(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
224 real mu_2p(grid%sm31:grid%em31,grid%sm33:grid%em33)
225
226 real alq(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
227 real mu_2q(grid%sm31:grid%em31,grid%sm33:grid%em33)
228 real pq(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
229 real ph_2q(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
230
231 real wwr(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
232 real u_2r(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
233 real v_2r(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
234
235 real w_2s(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
236 real wws(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
237 real u_2s(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
238 real v_2s(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
239 real t_2s(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
240 real mu_2s(grid%sm31:grid%em31,grid%sm33:grid%em33)
241 real ph_2s(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
242 real as(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
243 real t_2saves(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
244
245 real mudf_ma(grid%sm31:grid%em31,grid%sm33:grid%em33)
246 !*********************************************************************
247 !----------------------------------------------------------------------------------------
248 real ww_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
249 real xkmhd_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
250 real ru_tendf_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
251 real rv_tendf_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
252 real rw_tendf_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
253 real t_tendf_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
254 real ph_tendf_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
255 real pm1_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
256 real moist_tend_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist)
257 real dnw_keep3(grid%sm32:grid%em32)
258
259 !----------------------------------------------------------------------------------------
260 !----------------------------------------------------------------------------------------
261 real ww_m_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
262 real ru_m_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
263 real rv_m_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
264 real a_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
265 real alpha_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
266 real gamma_keep3(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
267 !----------------------------------------------------------------------------------------
268 !----------------------------------------------------------------------------------------
269
270 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: u_1_keep3
271 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: u_2_keep3
272 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: v_1_keep3
273 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: v_2_keep3
274 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: w_1_keep3
275 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: w_2_keep3
276 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: t_1_keep3
277 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: t_2_keep3
278 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: ph_1_keep3
279 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: ph_2_keep3
280 real, DIMENSION(grid%sm31:grid%em31,grid%sm33:grid%em33) :: mu_1_keep3
281 real, DIMENSION(grid%sm31:grid%em31,grid%sm33:grid%em33) :: mu_2_keep3
282 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: p_keep3
283 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: al_keep3
284 !----------------------------------------------------------------------------------------
285
286 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: u_1_keep4
287 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: u_2_keep4
288 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: v_1_keep4
289 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: v_2_keep4
290 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: w_1_keep4
291 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: w_2_keep4
292 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: t_1_keep4
293 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: t_2_keep4
294 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: ph_1_keep4
295 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: ph_2_keep4
296 real, DIMENSION(grid%sm31:grid%em31,grid%sm33:grid%em33) :: mu_1_keep4
297 real, DIMENSION(grid%sm31:grid%em31,grid%sm33:grid%em33) :: mu_2_keep4
298 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist) :: moist_keep4
299 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: p_keep4
300 real, DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: al_keep4
301 !----------------------------------------------------------------------------------------
302 integer :: iy
303
304 real ru_mz(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
305 real rv_mz(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
306 real ru_tendfz(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33)
307 real moist_tendz(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33,num_moist)
308 !----------------------------------------
309
310 #ifdef DEREF_KLUDGE
311 INTEGER :: sm31 , em31 , sm32 , em32 , sm33 , em33
312 INTEGER :: sm31x, em31x, sm32x, em32x, sm33x, em33x
313 INTEGER :: sm31y, em31y, sm32y, em32y, sm33y, em33y
314 #endif
315
316 #include <bench_solve_em_def.h>
317
318 #include "deref_kludge.h"
319
320 #define COPY_IN
321 #include <em_scalar_derefs.inc>
322 #ifdef DM_PARALLEL
323 # define REGISTER_I1
324 # include <em_data_calls.inc>
325 #endif
326
327 !<DESCRIPTION>
328 !<pre>
329 ! solve_em is the main driver for advancing a grid a single timestep.
330 ! It is a mediation-layer routine -> DM and SM calls are made where
331 ! needed for parallel processing.
332 !
333 ! solve_em can integrate the equations using 3 time-integration methods
334 !
335 ! - 3rd order Runge-Kutta time integration (recommended)
336 !
337 ! - 2nd order Runge-Kutta time integration
338 !
339 ! - Leapfrog time integration
340 ! (note: the leapfrog scheme is not correctly implemented
341 ! for most of the physics)
342 !
343 ! The main sections of solve_em are
344 !
345 ! (1) Runge-Kutta (RK) loop
346 !
347 ! (2) Non-timesplit physics (i.e., tendencies computed for updating
348 ! model state variables during the first RK sub-step (loop)
349 !
350 ! (3) Small (acoustic, sound) timestep loop - within the RK sub-steps
351 !
352 ! (4) Scalar advance for moist and chem scalar variables (and TKE)
353 ! within the RK sub-steps.
354 !
355 ! (5) time-split physics (after the RK step), currently this includes
356 ! only microphyics
357 !
358 ! A more detailed description of these sections follows.
359 !</pre>
360 !</DESCRIPTION>
361
362 #include <bench_solve_em_init.h>
363
364 ! xyh
365 ! return
366
367 !----------------------------------------------
368 ! SAVE REQUIRED INPUT VARIABLES
369 !----------------------------------------------
370 xkmhdi(:,:,:) = xkmhd(:,:,:)
371 ww_mh(:,:,:) = ww_m(:,:,:)
372 w_2j(:,:,:) = w_2(:,:,:)
373 w_1i(:,:,:) = w_1(:,:,:)
374 v_2i(:,:,:) = v_2(:,:,:)
375 v_1i(:,:,:) = v_1(:,:,:)
376 u_2i(:,:,:) = u_2(:,:,:)
377 u_1i(:,:,:) = u_1(:,:,:)
378 t_tendfi(:,:,:) = t_tendf(:,:,:)
379 t_2j(:,:,:) = t_2(:,:,:)
380 t_1i(:,:,:) = t_1(:,:,:)
381 rw_tendfi(:,:,:) = rw_tendf(:,:,:)
382 rv_tendfi(:,:,:) = rv_tendf(:,:,:)
383 rv_mh(:,:,:) = rv_m(:,:,:)
384 ru_tendfi(:,:,:) = ru_tendf(:,:,:)
385 ru_mh(:,:,:) = ru_m(:,:,:)
386 pm1i(:,:,:) = pm1(:,:,:)
387 ph_tendfi(:,:,:) = ph_tendf(:,:,:)
388 ph_2j(:,:,:) = ph_2(:,:,:)
389 ph_1i(:,:,:) = ph_1(:,:,:)
390 pj(:,:,:) = p(:,:,:)
391 number_of_small_timestepsh = number_of_small_timesteps
392 mu_2j(:,:) = mu_2(:,:)
393 mu_1i(:,:) = mu_1(:,:)
394 moist_tendi(:,:,:,:) = moist_tend(:,:,:,:)
395 moist_i(:,:,:,:) = moist(:,:,:,:)
396 gammah(:,:,:) = gamma(:,:,:)
397 dts_rkh = dts_rk
398 dt_rkh = dt_rk
399 alphah(:,:,:) = alpha(:,:,:)
400 alj(:,:,:) = al(:,:,:)
401 ah(:,:,:) = a(:,:,:)
402
403 !----------------------------------------------
404 ! RESET LOCAL ADJOINT VARIABLES
405 !----------------------------------------------
406 a_a(:,:,:) = 0.
407 a_advect_tend(:,:,:) = 0.
408 a_alpha(:,:,:) = 0.
409 a_bn2h(:,:,:) = 0.
410 a_c2a(:,:,:) = 0.
411 a_cqu(:,:,:) = 0.
412 a_cqv(:,:,:) = 0.
413 a_cqw(:,:,:) = 0.
414 a_gamma(:,:,:) = 0.
415 a_moist_tend(:,:,:,:) = 0.
416 a_mu_save(:,:) = 0.
417 a_mu_tend(:,:) = 0.
418 a_muave(:,:) = 0.
419 a_mut(:,:) = 0.
420 a_muts(:,:) = 0.
421 a_muu(:,:) = 0.
422 a_muus(:,:) = 0.
423 a_muv(:,:) = 0.
424 a_muvs(:,:) = 0.
425 a_p8w(:,:,:) = 0.
426 a_p_phy(:,:,:) = 0.
427 a_ph_save(:,:,:) = 0.
428 a_ph_tend(:,:,:) = 0.
429 a_ph_tendf(:,:,:) = 0.
430 a_pi_phyh(:,:,:) = 0.
431 a_pm1(:,:,:) = 0.
432 a_rqc_bth(:,:,:,:) = 0.
433 a_rqg_bth(:,:,:,:) = 0.
434 a_rqi_bth(:,:,:,:) = 0.
435 a_rqr_bth(:,:,:,:) = 0.
436 a_rqs_bth(:,:,:,:) = 0.
437 a_ru_m(:,:,:) = 0.
438 a_ru_tend(:,:,:) = 0.
439 a_ru_tendf(:,:,:) = 0.
440 a_rv_m(:,:,:) = 0.
441 a_rv_tend(:,:,:) = 0.
442 a_rv_tendf(:,:,:) = 0.
443 a_rw_tend(:,:,:) = 0.
444 a_rw_tendf(:,:,:) = 0.
445 a_t8w(:,:,:) = 0.
446 a_t_2save(:,:,:) = 0.
447 a_t_phy(:,:,:) = 0.
448 a_t_save(:,:,:) = 0.
449 a_t_tend(:,:,:) = 0.
450 a_t_tendf(:,:,:) = 0.
451 a_th_phy(:,:,:) = 0.
452 a_u_save(:,:,:) = 0.
453 a_v_save(:,:,:) = 0.
454 a_w_save(:,:,:) = 0.
455 a_ww1(:,:,:) = 0.
456 a_ww_m(:,:,:) = 0.
457 a_z_at_wh(:,:,:) = 0.
458
459
460 !----------------------------------------------
461 ! ROUTINE BODY
462 !----------------------------------------------
463 cfn = grid%cfn
464 ! recompute : cfn
465 cfn1 = grid%cfn1
466 ! recompute : cfn1
467 step_number = grid%step_number
468 ! recompute : step_number
469 rdx = grid%rdx
470 ! recompute : rdx
471 rdy = grid%rdy
472 ! recompute : rdy
473 cf1 = grid%cf1
474 ! recompute : cf1
475 cf2 = grid%cf2
476 ! recompute : cf2
477 cf3 = grid%cf3
478 ! recompute : cf3
479 dtbc = grid%dtbc
480 ! recompute : dtbc
481 dx = grid%dx
482 ! recompute : dx
483 dy = grid%dy
484 ! recompute : dy
485 dt = grid%dt
486 ! recompute : dt
487 rk_ord = grid%rk_ord
488 ! recompute : rk_ord
489 diff_opt = grid%diff_opt
490 ! recompute : diff_opt
491 damp_opt = grid%damp_opt
492 ! recompute : damp_opt
493 zdamp = grid%zdamp
494 ! recompute : zdamp
495 dampcoef = grid%dampcoef
496 ! recompute : dampcoef
497 khdif = grid%khdif
498 ! recompute : khdif
499 kvdif = grid%kvdif
500 ! recompute : kvdif
501 smdiv = grid%smdiv
502 ! recompute : smdiv
503 emdiv = grid%emdiv
504 ! recompute : emdiv
505 epssm = grid%epssm
506 ! recompute : epssm
507 non_hydrostatic = grid%non_hydrostatic
508 ! recompute : non_hydrostatic
509 time_step_sound = grid%time_step_sound
510 ! recompute : time_step_sound
511 kh_tke_upper_bound = grid%kh_tke_upper_bound
512 ! recompute : kh_tke_upper_bound
513 spec_bdy_width = grid%spec_bdy_width
514 ! recompute : spec_bdy_width
515 spec_zone = grid%spec_zone
516 ! recompute : spec_zone
517 relax_zone = grid%relax_zone
518 ! recompute : relax_zone
519 call get_ijk_from_grid( grid,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe )
520 ! recompute : ide,ids,ime,ims,ipe,ips,jde,jds,jme,jms,jpe,jps,kde,kds,k
521 ! me,kms,kpe,kps
522 k_start = kps
523 ! recompute : k_start
524 k_end = kpe
525 ! recompute : k_end
526 ijds = min(ids,jds)
527 ! recompute : ijds
528 ijde = max(ide,jde)
529 ! recompute : ijde
530 num_3d_m = num_moist
531 ! recompute : num_3d_m
532 call set_tiles( grid,ids,ide,jds,jde,ips,ipe,jps,jpe )
533 ! recompute : grid
534 rk_order = config_flags%rk_ord
535 ! recompute : rk_order
536 dts = dt/float(time_step_sound)
537 ! recompute : dts
538
539 !------------------------------
540 !------------------------------
541 runge_kutta_looq: do rk_step = 1, rk_order
542 if (rk_order .eq. 1) then
543 if (step_number .ne. 1) then
544 number_of_small_timesteps = 2*time_step_sound
545 dt_rk = dt
546 else
547 number_of_small_timesteps = time_step_sound
548 dt_rk = dt/2.
549 endif
550 dts_rk = dts
551 else if (rk_order .eq. 2) then
552 if (rk_step .eq. 1) then
553 dt_rk = 0.5*dt
554 dts_rk = dts
555 number_of_small_timesteps = time_step_sound/2
556 else
557 dt_rk = dt
558 dts_rk = dts
559 number_of_small_timesteps = time_step_sound
560 endif
561 else if (rk_order .eq. 3) then
562 if (rk_step .eq. 1) then
563 dt_rk = dt/3.
564 dts_rk = dt_rk
565 number_of_small_timesteps = 1
566 else if (rk_step .eq. 2) then
567 dt_rk = 0.5*dt
568 dts_rk = dts
569 number_of_small_timesteps = time_step_sound/2
570 else
571 dt_rk = dt
572 dts_rk = dts
573 number_of_small_timesteps = time_step_sound
574 endif
575 endif
576
577 do ij = 1, grid%num_tiles
578 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
579 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
580 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
581 end do
582 rk_step_is_ong: if (rk_step .eq. 1) then
583 do ij = 1, grid%num_tiles
584 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,&
585 moist_tend,chem_tend,scalar_tend,num_3d_m,num_3d_c,num_3d_s,rk_step,&
586 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
587 end do
588 do ij = 1, grid%num_tiles
589 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
590 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
591 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
592 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
593 end do
594 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
595 do ij = 1, grid%num_tiles
596 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
597 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
598 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
599 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
600 end do
601 endif
602 endif rk_step_is_ong
603 do ij = 1, grid%num_tiles
604 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
605 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
606 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
607 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
608 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
609 end do
610 do ij = 1, grid%num_tiles
611 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
612 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
613 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,&
614 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
615 endif
616 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
617 &ph_save,t_save,rk_step,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
618 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
619 if (config_flags%specified .or. config_flags%nested) then
620 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
621 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
622 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
623 endif
624 end do
625 do ij = 1, grid%num_tiles
626 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
627 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
628 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
629 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
630 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
631 if (non_hydrostatic) then
632 call calc_coef_w( a,alpha,gamma,mut,cqw,rdn,rdnw,c2a,dts,g,epssm,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
633 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
634 endif
635 end do
636 small_stepu: do iteration = 1, number_of_small_timesteps
637 do ij = 1, grid%num_tiles
638 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
639 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
640 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
641 if (config_flags%specified .or. config_flags%nested) then
642 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
643 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
644 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
645 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
646 endif
647 end do
648 do ij = 1, grid%num_tiles
649 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
650 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
651 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
652 if (config_flags%specified .or. config_flags%nested) then
653 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
654 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
655 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
656 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
657 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
658 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
659 endif
660 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration,number_of_small_timesteps,ids,&
661 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
662 if (non_hydrostatic) then
663 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,alb,&
664 &a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
665 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
666 endif
667 if (config_flags%specified .or. config_flags%nested) then
668 if (non_hydrostatic) then
669 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,&
670 &ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
671 if (config_flags%specified) then
672 call zero_grad_bdy( w_2,'w',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
673 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
674 else
675 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
676 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
677 endif
678 endif
679 endif
680 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,iteration,ids,ide,jds,jde,&
681 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
682 end do
683 end do small_stepu
684 do ij = 1, grid%num_tiles
685 call calc_mu_uv_1( config_flags,muts,muus,muvs,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
686 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
687 ! call small_step_finish( u_2,u_1,v_2,v_1,w_2,w_1,t_2,t_1,ph_2,ph_1,ww,ww1,mu_2,mu_1,mut,muts,muu,muus,muv,muvs,u_save,v_save,&
688 !&w_save,t_save,ph_save,mu_save,msfu,msfv,msft,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
689 !&grid%j_start(ij),grid%j_end(ij),k_start,k_end )
690 !-----The original
691 CALL small_step_finish( u_2, u_1, v_2, v_1, w_2, w_1, &
692 t_2, t_1, ph_2, ph_1, ww, ww1, &
693 mu_2, mu_1, &
694 mut, muts, muu, muus, muv, muvs, &
695 u_save, v_save, w_save, &
696 t_save, ph_save, mu_save, &
697 msfu, msfv, msft, &
698 h_diabatic, &
699 number_of_small_timesteps,dts_rk, &
700 ids, ide, jds, jde, kds, kde, &
701 ims, ime, jms, jme, kms, kme, &
702 grid%i_start(ij), grid%i_end(ij), &
703 grid%j_start(ij), grid%j_end(ij), &
704 k_start , k_end )
705 end do
706 moist_scalar_advancg: if (num_3d_m .ge. param_first_scalar) then
707 moist_variable_loor: do im = param_first_scalar, num_3d_m
708 moist_tile_loop_6: do ij = 1, grid%num_tiles
709 call rk_scalar_tend( im,im,config_flags,rk_step,dt_rk,ru_m,rv_m,ww_m,mut,alt,moist(ims,kms,jms,&
710 &im),moist_tend(ims,kms,jms,im),advect_tend,rqvften,qv_base, .true. ,fnm,fnp,msfu,msfv,msft,rdx,rdy,rdn,rdnw,khdif,kvdif,&
711 &xkmhd,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
712 &j_end(ij),k_start,k_end )
713 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
714 if (im .eq. p_qv) then
715 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqv_b,rqv_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
716 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
717 endif
718 endif
719 if (config_flags%nested .and. rk_step .eq. 1) then
720 if (im .eq. p_qc) then
721 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqc_b,rqc_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
722 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
723 else if (im .eq. p_qr) then
724 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqr_b,rqr_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
725 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
726 else if (im .eq. p_qi) then
727 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqi_b,rqi_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
728 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
729 else if (im .eq. p_qs) then
730 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqs_b,rqs_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
731 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
732 else if (im .eq. p_qg) then
733 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqg_b,rqg_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
734 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
735 endif
736 endif
737 end do moist_tile_loop_6
738 moist_tile_loop_5: do ij = 1, grid%num_tiles
739 call rk_update_scalar( im,im,moist_old(ims,kms,jms,im),moist(ims,kms,jms,im),moist_tend(ims,kms,jms,im),advect_tend,msft,&
740 &mu_1,mu_2,mub,rk_step,dt_rk,spec_zone,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
741 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
742 end do moist_tile_loop_5
743 end do moist_variable_loor
744 endif moist_scalar_advancg
745 do ij = 1, grid%num_tiles
746 call calc_p_rho_phi( moist,num_3d_m,al,alb,mu_2,muts,ph_2,p,pb,t_2,p0,t0,znu,dnw,rdnw,rdn,non_hydrostatic,ids,ide,jds,jde,&
747 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
748 if ( .not. non_hydrostatic) then
749 call diagnose_w( ph_tend,ph_2,ph_1,w_2,muts,dt_rk,u_2,v_2,ht,cf1,cf2,cf3,rdx,rdy,msft,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
750 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
751 endif
752 end do
753 end do runge_kutta_looq
754 ! recompute : dt_rk,moist,mu_2,muts,ph_2,ph_tend,t_2
755 if (config_flags%mp_physics .ne. 0) then
756 if (config_flags%specified .or. config_flags%nested) then
757 sz = spec_zone
758 else
759 sz = 0
760 endif
761 ! recompute : sz
762 do ij = 1, grid%num_tiles
763 its = max(grid%i_start(ij),ids+sz)
764 ! recompute : its
765 ite = min(grid%i_end(ij),ide-1-sz)
766 ! recompute : ite
767 jts = max(grid%j_start(ij),jds+sz)
768 ! recompute : jts
769 jte = min(grid%j_end(ij),jde-1-sz)
770 ! recompute : jte
771 if ( .not. non_hydrostatic) then
772 call a_diagnose_w( ph_tend,a_ph_tend,a_ph_2,a_ph_1,a_w_2,muts,a_muts,dt_rk,a_u_2,a_v_2,ht,cf1,cf2,cf3,rdx,rdy,msft,ide,jde,&
773 &ims,ime,jms,jme,kms,kme,its,ite,jts,jte,k_end )
774 endif
775 call a_calc_p_rho_phi( moist,a_moist,num_3d_m,al,a_al,alb,mu_2,a_mu_2,muts,a_muts,ph_2,a_ph_2,p,a_p,pb,t_2,a_t_2,p0,t0,dnw,&
776 &rdnw,rdn,non_hydrostatic,ide,jde,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,k_start,k_end )
777 end do
778 endif
779
780 !------------------------------
781 !------------------------------
782
783 a_runge_kutta_loop: do rk_step = rk_order, 1, -1
784
785
786 a(:,:,:) = ah(:,:,:)
787 al(:,:,:) = alj(:,:,:)
788 alpha(:,:,:) = alphah(:,:,:)
789 ! dt_rk = dt_rkh
790 ! dts_rk = dts_rkh
791 gamma(:,:,:) = gammah(:,:,:)
792 moist(:,:,:,:) = moist_i(:,:,:,:)
793 moist_tend(:,:,:,:) = moist_tendi(:,:,:,:)
794 mu_1(:,:) = mu_1i(:,:)
795 mu_2(:,:) = mu_2j(:,:)
796 ! number_of_small_timesteps = number_of_small_timestepsh
797 p(:,:,:) = pj(:,:,:)
798 ph_1(:,:,:) = ph_1i(:,:,:)
799 ph_2(:,:,:) = ph_2j(:,:,:)
800 ph_tendf(:,:,:) = ph_tendfi(:,:,:)
801 pm1(:,:,:) = pm1i(:,:,:)
802 ru_m(:,:,:) = ru_mh(:,:,:)
803 ru_tendf(:,:,:) = ru_tendfi(:,:,:)
804 rv_m(:,:,:) = rv_mh(:,:,:)
805 rv_tendf(:,:,:) = rv_tendfi(:,:,:)
806 rw_tendf(:,:,:) = rw_tendfi(:,:,:)
807 t_1(:,:,:) = t_1i(:,:,:)
808 t_2(:,:,:) = t_2j(:,:,:)
809 t_tendf(:,:,:) = t_tendfi(:,:,:)
810 u_1(:,:,:) = u_1i(:,:,:)
811 u_2(:,:,:) = u_2i(:,:,:)
812 v_1(:,:,:) = v_1i(:,:,:)
813 v_2(:,:,:) = v_2i(:,:,:)
814 w_1(:,:,:) = w_1i(:,:,:)
815 w_2(:,:,:) = w_2j(:,:,:)
816 ww_m(:,:,:) = ww_mh(:,:,:)
817 xkmhd(:,:,:) = xkmhdi(:,:,:)
818
819 runge_kutta_loor: do rk_step1 = 1, rk_step-1
820 if (rk_order .eq. 1) then
821 if (step_number .ne. 1) then
822 number_of_small_timesteps = 2*time_step_sound
823 dt_rk = dt
824 else
825 number_of_small_timesteps = time_step_sound
826 dt_rk = dt/2.
827 endif
828 dts_rk = dts
829 else if (rk_order .eq. 2) then
830 if (rk_step1 .eq. 1) then
831 dt_rk = 0.5*dt
832 dts_rk = dts
833 number_of_small_timesteps = time_step_sound/2
834 else
835 dt_rk = dt
836 dts_rk = dts
837 number_of_small_timesteps = time_step_sound
838 endif
839 else if (rk_order .eq. 3) then
840 if (rk_step1 .eq. 1) then
841 dt_rk = dt/3.
842 dts_rk = dt_rk
843 number_of_small_timesteps = 1
844 else if (rk_step1 .eq. 2) then
845 dt_rk = 0.5*dt
846 dts_rk = dts
847 number_of_small_timesteps = time_step_sound/2
848 else
849 dt_rk = dt
850 dts_rk = dts
851 number_of_small_timesteps = time_step_sound
852 endif
853 endif
854 do ij = 1, grid%num_tiles
855 call rk_step_prep( config_flags,rk_step1,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,&
856 &alb,cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij)&
857 &,grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
858 end do
859 rk_step_is_onze: if (rk_step1 .eq. 1) then
860 do ij = 1, grid%num_tiles
861 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,&
862 &rk_step1,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),&
863 &k_start,k_end )
864 end do
865 do ij = 1, grid%num_tiles
866 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,&
867 &v_phy,p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,&
868 &rqccuten,rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
869 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
870 end do
871 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
872 do ij = 1, grid%num_tiles
873 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
874 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
875 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
876 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
877 end do
878 endif
879 endif rk_step_is_onze
880
881 do ij = 1, grid%num_tiles
882 call rk_tendency( config_flags,rk_step1,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
883 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,&
884 &h_diabatic,phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,&
885 &e,sina,cosa,fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,&
886 &jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
887 end do
888
889 do ij = 1, grid%num_tiles
890 if ((config_flags%specified .or. config_flags%nested) .and. rk_step1 .eq. 1) then
891 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
892 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,&
893 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
894 endif
895 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
896 &ph_save,t_save,rk_step1,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
897 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
898 if (config_flags%specified .or. config_flags%nested) then
899 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
900 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
901 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
902 endif
903 end do
904
905
906 do ij = 1, grid%num_tiles
907 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
908 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step1,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
909 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
910 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
911 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
912 if (non_hydrostatic) then
913 call calc_coef_w( a,alpha,gamma,mut,cqw,rdn,rdnw,c2a,dts,g,epssm,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
914 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
915 endif
916 end do
917
918 small_stepzt: do iteration = 1, number_of_small_timesteps
919 do ij = 1, grid%num_tiles
920 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
921 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
922 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
923 if (config_flags%specified .or. config_flags%nested) then
924 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
925 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
926 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
927 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
928 endif
929 end do
930 do ij = 1, grid%num_tiles
931 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
932 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
933 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
934 if (config_flags%specified .or. config_flags%nested) then
935 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
936 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
937 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
938 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
939 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
940 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
941 endif
942 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration,number_of_small_timesteps,ids,&
943 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
944 if (non_hydrostatic) then
945 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,&
946 &alb,a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,&
947 &ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
948 endif
949 if (config_flags%specified .or. config_flags%nested) then
950 if (non_hydrostatic) then
951 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,&
952 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
953 if (config_flags%specified) then
954 call zero_grad_bdy( w_2,'w',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
955 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
956 else
957 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,&
958 &ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
959 endif
960 endif
961 endif
962 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,iteration,ids,ide,jds,&
963 &jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
964 end do
965 end do small_stepzt
966
967 do ij = 1, grid%num_tiles
968 call calc_mu_uv_1( config_flags,muts,muus,muvs,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
969 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
970 ! call small_step_finish( u_2,u_1,v_2,v_1,w_2,w_1,t_2,t_1,ph_2,ph_1,ww,ww1,mu_2,mu_1,mut,muts,muu,muus,muv,muvs,u_save,v_save,&
971 !&w_save,t_save,ph_save,mu_save,msfu,msfv,msft,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij)&
972 !&,grid%j_start(ij),grid%j_end(ij),k_start,k_end )
973 !-----The original
974
975
976 CALL small_step_finish( u_2, u_1, v_2, v_1, w_2, w_1, &
977 t_2, t_1, ph_2, ph_1, ww, ww1, &
978 mu_2, mu_1, &
979 mut, muts, muu, muus, muv, muvs, &
980 u_save, v_save, w_save, &
981 t_save, ph_save, mu_save, &
982 msfu, msfv, msft, &
983 h_diabatic, &
984 number_of_small_timesteps,dts_rk, &
985 ids, ide, jds, jde, kds, kde, &
986 ims, ime, jms, jme, kms, kme, &
987 grid%i_start(ij), grid%i_end(ij), &
988 grid%j_start(ij), grid%j_end(ij), &
989 k_start , k_end )
990 end do
991
992 moist_scalar_advancj: if (num_3d_m .ge. param_first_scalar) then
993 moist_variable_loou: do im = param_first_scalar, num_3d_m
994 moist_tile_loop_9d: do ij = 1, grid%num_tiles
995 call rk_scalar_tend( im,im,config_flags,rk_step1,dt_rk,ru_m,rv_m,ww_m,mut,alt,moist(ims,kms,&
996 &jms,im),moist_tend(ims,kms,jms,im),advect_tend,rqvften,qv_base, .true. ,fnm,fnp,msfu,msfv,msft,rdx,rdy,rdn,rdnw,khdif,&
997 &kvdif,xkmhd,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),&
998 &grid%j_end(ij),k_start,k_end )
999 if ((config_flags%specified .or. config_flags%nested) .and. rk_step1 .eq. 1) then
1000 if (im .eq. p_qv) then
1001 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqv_b,rqv_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1002 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1003 endif
1004 endif
1005 if (config_flags%nested .and. rk_step1 .eq. 1) then
1006 if (im .eq. p_qc) then
1007 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqc_b,rqc_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1008 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1009 else if (im .eq. p_qr) then
1010 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqr_b,rqr_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1011 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1012 else if (im .eq. p_qi) then
1013 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqi_b,rqi_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1014 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1015 else if (im .eq. p_qs) then
1016 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqs_b,rqs_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1017 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1018 else if (im .eq. p_qg) then
1019 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqg_b,rqg_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1020 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1021 endif
1022 endif
1023 end do moist_tile_loop_9d
1024 moist_tile_loop_9c: do ij = 1, grid%num_tiles
1025 call rk_update_scalar( im,im,moist_old(ims,kms,jms,im),moist(ims,kms,jms,im),moist_tend(ims,kms,jms,im),advect_tend,msft,&
1026 &mu_1,mu_2,mub,rk_step1,dt_rk,spec_zone,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
1027 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1028 end do moist_tile_loop_9c
1029 end do moist_variable_loou
1030 endif moist_scalar_advancj
1031
1032 do ij = 1, grid%num_tiles
1033 call calc_p_rho_phi( moist,num_3d_m,al,alb,mu_2,muts,ph_2,p,pb,t_2,p0,t0,znu,dnw,rdnw,rdn,non_hydrostatic,ids,ide,jds,jde,&
1034 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1035 if ( .not. non_hydrostatic) then
1036 call diagnose_w( ph_tend,ph_2,ph_1,w_2,muts,dt_rk,u_2,v_2,ht,cf1,cf2,cf3,rdx,rdy,msft,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
1037 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1038 endif
1039 end do
1040
1041 end do runge_kutta_loor
1042
1043 !------------------------------
1044 ww_m_keep3(:,:,:) = ww_m(:,:,:)
1045 rv_m_keep3(:,:,:) = rv_m(:,:,:)
1046 ru_m_keep3(:,:,:) = ru_m(:,:,:)
1047 a_keep3(:,:,:) = a(:,:,:)
1048 alpha_keep3(:,:,:) = alpha(:,:,:)
1049 gamma_keep3(:,:,:) = gamma(:,:,:)
1050 !------------------------------
1051
1052 xkmhdh(:,:,:) = xkmhd(:,:,:)
1053 w_2h(:,:,:) = w_2(:,:,:)
1054 w_1h(:,:,:) = w_1(:,:,:)
1055 v_2h(:,:,:) = v_2(:,:,:)
1056 v_1h(:,:,:) = v_1(:,:,:)
1057 u_2h(:,:,:) = u_2(:,:,:)
1058 u_1h(:,:,:) = u_1(:,:,:)
1059 t_tendfh(:,:,:) = t_tendf(:,:,:)
1060 t_2h(:,:,:) = t_2(:,:,:)
1061 t_1h(:,:,:) = t_1(:,:,:)
1062 rw_tendfh(:,:,:) = rw_tendf(:,:,:)
1063 rv_tendfh(:,:,:) = rv_tendf(:,:,:)
1064 ru_tendfh(:,:,:) = ru_tendf(:,:,:)
1065 pm1h(:,:,:) = pm1(:,:,:)
1066 ph_tendfh(:,:,:) = ph_tendf(:,:,:)
1067 ph_2h(:,:,:) = ph_2(:,:,:)
1068 ph_1h(:,:,:) = ph_1(:,:,:)
1069 ph(:,:,:) = p(:,:,:)
1070 mu_2h(:,:) = mu_2(:,:)
1071 mu_1h(:,:) = mu_1(:,:)
1072 moist_tendh(:,:,:,:) = moist_tend(:,:,:,:)
1073 moist_h(:,:,:,:) = moist(:,:,:,:)
1074 dnwh(:) = dnw(:)
1075 alh(:,:,:) = al(:,:,:)
1076 !----------------------------------------
1077 mudf_ma(:,:) = mudf(:,:)
1078 !----------------------------------------
1079
1080
1081 if (rk_order .eq. 1) then
1082 if (step_number .ne. 1) then
1083 number_of_small_timesteps = 2*time_step_sound
1084 dt_rk = dt
1085 else
1086 number_of_small_timesteps = time_step_sound
1087 dt_rk = dt/2.
1088 endif
1089 dts_rk = dts
1090 else if (rk_order .eq. 2) then
1091 if (rk_step .eq. 1) then
1092 dt_rk = 0.5*dt
1093 dts_rk = dts
1094 number_of_small_timesteps = time_step_sound/2
1095 else
1096 dt_rk = dt
1097 dts_rk = dts
1098 number_of_small_timesteps = time_step_sound
1099 endif
1100 else if (rk_order .eq. 3) then
1101 if (rk_step .eq. 1) then
1102 dt_rk = dt/3.
1103 dts_rk = dt_rk
1104 number_of_small_timesteps = 1
1105 else if (rk_step .eq. 2) then
1106 dt_rk = 0.5*dt
1107 dts_rk = dts
1108 number_of_small_timesteps = time_step_sound/2
1109 else
1110 dt_rk = dt
1111 dts_rk = dts
1112 number_of_small_timesteps = time_step_sound
1113 endif
1114 endif
1115 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
1116 do ij = 1, grid%num_tiles
1117 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
1118 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1119 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1120 end do
1121 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
1122 rk_step_is_onzc: if (rk_step .eq. 1) then
1123 do ij = 1, grid%num_tiles
1124 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,rk_step,&
1125 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1126 end do
1127 do ij = 1, grid%num_tiles
1128 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
1129 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
1130 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1131 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1132 end do
1133 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
1134 do ij = 1, grid%num_tiles
1135 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
1136 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
1137 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1138 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1139 end do
1140 endif
1141 endif rk_step_is_onzc
1142 ! recompute : moist_tend,ph_tendf,ru_tendf,rv_tendf,rw_tendf,t_tendf,xk
1143 ! mhd
1144 do ij = 1, grid%num_tiles
1145 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
1146 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
1147 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
1148 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
1149 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1150 end do
1151
1152 ! recompute : cqw,mu_tend,ph_save,ph_tend,ru_tend,ru_tendf,rv_tend,rv_t
1153 ! endf,rw_tend,rw_tendf,t_save,t_tend,t_tendf,u_save,v_save,w_save
1154 do ij = 1, grid%num_tiles
1155 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1156 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
1157 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,&
1158 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1159 endif
1160 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
1161 &ph_save,t_save,rk_step,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
1162 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1163 if (config_flags%specified .or. config_flags%nested) then
1164 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
1165 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
1166 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1167 endif
1168 end do
1169
1170 ! recompute : mu_tend,ph_tend,ru_tend,rv_tend,rw_tend,t_tend
1171 do ij = 1, grid%num_tiles
1172 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
1173 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1174 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1175 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
1176 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1177 if (non_hydrostatic) then
1178 call calc_coef_w( a,alpha,gamma,mut,cqw,rdn,rdnw,c2a,dts,g,epssm,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
1179 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1180 endif
1181 end do
1182
1183 small_stepzq: do iteration = 1, number_of_small_timesteps
1184 do ij = 1, grid%num_tiles
1185 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
1186 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1187 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1188 if (config_flags%specified .or. config_flags%nested) then
1189 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1190 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1191 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1192 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1193 endif
1194 end do
1195 do ij = 1, grid%num_tiles
1196 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
1197 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1198 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1199 if (config_flags%specified .or. config_flags%nested) then
1200 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1201 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1202 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
1203 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
1204 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
1205 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
1206 endif
1207 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration,number_of_small_timesteps,ids,&
1208 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1209 if (non_hydrostatic) then
1210 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,alb,&
1211 &a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
1212 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1213 endif
1214 if (config_flags%specified .or. config_flags%nested) then
1215 if (non_hydrostatic) then
1216 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,&
1217 &ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1218 if (config_flags%specified) then
1219 call zero_grad_bdy( w_2,'w',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
1220 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1221 else
1222 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
1223 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1224 endif
1225 endif
1226 endif
1227 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,iteration,ids,ide,jds,jde,&
1228 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1229 end do
1230 end do small_stepzq
1231 ! recompute : mu_2,mudf,muts,ph_2,pm1,ru_m,rv_m,t_2,u_2,v_2,w_2,ww_m
1232
1233 do ij = 1, grid%num_tiles
1234 call calc_mu_uv_1( config_flags,muts,muus,muvs,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1235 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1236 ! call small_step_finish( u_2,u_1,v_2,v_1,w_2,w_1,t_2,t_1,ph_2,ph_1,ww,ww1,mu_2,mu_1,mut,muts,muu,muus,muv,muvs,u_save,v_save,&
1237 !&w_save,t_save,ph_save,mu_save,msfu,msfv,msft,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1238 !&grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1239 !-----The original
1240 CALL small_step_finish( u_2, u_1, v_2, v_1, w_2, w_1, &
1241 t_2, t_1, ph_2, ph_1, ww, ww1, &
1242 mu_2, mu_1, &
1243 mut, muts, muu, muus, muv, muvs, &
1244 u_save, v_save, w_save, &
1245 t_save, ph_save, mu_save, &
1246 msfu, msfv, msft, &
1247 h_diabatic, &
1248 number_of_small_timesteps,dts_rk, &
1249 ids, ide, jds, jde, kds, kde, &
1250 ims, ime, jms, jme, kms, kme, &
1251 grid%i_start(ij), grid%i_end(ij), &
1252 grid%j_start(ij), grid%j_end(ij), &
1253 k_start , k_end )
1254 end do
1255 ! recompute : mu_2,ph_2,t_2
1256
1257
1258 moist_scalar_advanch: if (num_3d_m .ge. param_first_scalar) then
1259 moist_variable_loos: do im = param_first_scalar, num_3d_m
1260 moist_tile_loop_8: do ij = 1, grid%num_tiles
1261 call rk_scalar_tend( im,im,config_flags,rk_step,dt_rk,ru_m,rv_m,ww_m,mut,alt,moist(ims,kms,jms,&
1262 &im),moist_tend(ims,kms,jms,im),advect_tend,rqvften,qv_base, .true. ,fnm,fnp,msfu,msfv,msft,rdx,rdy,rdn,rdnw,khdif,kvdif,&
1263 &xkmhd,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
1264 &j_end(ij),k_start,k_end )
1265 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1266 if (im .eq. p_qv) then
1267 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqv_b,rqv_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1268 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1269 endif
1270 endif
1271 if (config_flags%nested .and. rk_step .eq. 1) then
1272 if (im .eq. p_qc) then
1273 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqc_b,rqc_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1274 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1275 else if (im .eq. p_qr) then
1276 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqr_b,rqr_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1277 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1278 else if (im .eq. p_qi) then
1279 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqi_b,rqi_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1280 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1281 else if (im .eq. p_qs) then
1282 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqs_b,rqs_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1283 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1284 else if (im .eq. p_qg) then
1285 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqg_b,rqg_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1286 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1287 endif
1288 endif
1289 end do moist_tile_loop_8
1290 moist_tile_loop_7: do ij = 1, grid%num_tiles
1291 call rk_update_scalar( im,im,moist_old(ims,kms,jms,im),moist(ims,kms,jms,im),moist_tend(ims,kms,jms,im),advect_tend,msft,&
1292 &mu_1,mu_2,mub,rk_step,dt_rk,spec_zone,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
1293 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1294 end do moist_tile_loop_7
1295 end do moist_variable_loos
1296 endif moist_scalar_advanch
1297
1298 ! recompute : moist
1299 do ij = 1, grid%num_tiles
1300
1301 call calc_p_rho_phi( moist,num_3d_m,al,alb,mu_2,muts,ph_2,p,pb,t_2,p0,t0,znu,dnw,rdnw,rdn,non_hydrostatic,ids,ide,jds,jde,&
1302 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1303
1304 if ( .not. non_hydrostatic) then
1305 call a_diagnose_w( ph_tend,a_ph_tend,a_ph_2,a_ph_1,a_w_2,muts,a_muts,dt_rk,a_u_2,a_v_2,ht,cf1,cf2,cf3,rdx,rdy,msft,ide,jde,&
1306 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_end )
1307 endif
1308 call a_calc_p_rho_phi( moist,a_moist,num_3d_m,al,a_al,alb,mu_2,a_mu_2,muts,a_muts,ph_2,a_ph_2,p,a_p,pb,t_2,a_t_2,p0,t0,dnw,&
1309 &rdnw,rdn,non_hydrostatic,ide,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),&
1310 &k_start,k_end )
1311 end do
1312 !------------------------------
1313 !------------------------------
1314
1315 ww_m(:,:,:) = ww_m_keep3(:,:,:)
1316 rv_m(:,:,:) = rv_m_keep3(:,:,:)
1317 ru_m(:,:,:) = ru_m_keep3(:,:,:)
1318 a(:,:,:) = a_keep3(:,:,:)
1319 alpha(:,:,:) = alpha_keep3(:,:,:)
1320 gamma(:,:,:) = gamma_keep3(:,:,:)
1321 !------------------------------
1322
1323 xkmhd(:,:,:) = xkmhdh(:,:,:)
1324 w_2(:,:,:) = w_2h(:,:,:)
1325 w_1(:,:,:) = w_1h(:,:,:)
1326 v_2(:,:,:) = v_2h(:,:,:)
1327 v_1(:,:,:) = v_1h(:,:,:)
1328 u_2(:,:,:) = u_2h(:,:,:)
1329 u_1(:,:,:) = u_1h(:,:,:)
1330 t_tendf(:,:,:) = t_tendfh(:,:,:)
1331 t_2(:,:,:) = t_2h(:,:,:)
1332 t_1(:,:,:) = t_1h(:,:,:)
1333 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
1334 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
1335 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
1336 pm1(:,:,:) = pm1h(:,:,:)
1337 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
1338 ph_2(:,:,:) = ph_2h(:,:,:)
1339 ph_1(:,:,:) = ph_1h(:,:,:)
1340 p(:,:,:) = ph(:,:,:)
1341 mu_2(:,:) = mu_2h(:,:)
1342 mu_1(:,:) = mu_1h(:,:)
1343 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
1344 moist(:,:,:,:) = moist_h(:,:,:,:)
1345 dnw(:) = dnwh(:)
1346 al(:,:,:) = alh(:,:,:)
1347 !----------------------------------------
1348 mudf(:,:) = mudf_ma(:,:)
1349 !----------------------------------------
1350
1351
1352 if (rk_order .eq. 1) then
1353 if (step_number .ne. 1) then
1354 number_of_small_timesteps = 2*time_step_sound
1355 dt_rk = dt
1356 else
1357 number_of_small_timesteps = time_step_sound
1358 dt_rk = dt/2.
1359 endif
1360 dts_rk = dts
1361 else if (rk_order .eq. 2) then
1362 if (rk_step .eq. 1) then
1363 dt_rk = 0.5*dt
1364 dts_rk = dts
1365 number_of_small_timesteps = time_step_sound/2
1366 else
1367 dt_rk = dt
1368 dts_rk = dts
1369 number_of_small_timesteps = time_step_sound
1370 endif
1371 else if (rk_order .eq. 3) then
1372 if (rk_step .eq. 1) then
1373 dt_rk = dt/3.
1374 dts_rk = dt_rk
1375 number_of_small_timesteps = 1
1376 else if (rk_step .eq. 2) then
1377 dt_rk = 0.5*dt
1378 dts_rk = dts
1379 number_of_small_timesteps = time_step_sound/2
1380 else
1381 dt_rk = dt
1382 dts_rk = dts
1383 number_of_small_timesteps = time_step_sound
1384 endif
1385 endif
1386 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
1387 do ij = 1, grid%num_tiles
1388 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
1389 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1390 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1391 end do
1392 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
1393 rk_step_is_onzg: if (rk_step .eq. 1) then
1394 do ij = 1, grid%num_tiles
1395 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,&
1396 moist_tend,chem_tend,scalar_tend,num_3d_m,num_3d_c,num_3d_s,rk_step,&
1397 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1398 end do
1399 do ij = 1, grid%num_tiles
1400 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
1401 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
1402 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1403 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1404 end do
1405 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
1406 do ij = 1, grid%num_tiles
1407 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
1408 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
1409 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1410 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1411 end do
1412 endif
1413 endif rk_step_is_onzg
1414 ! recompute : moist_tend,ph_tendf,ru_tendf,rv_tendf,rw_tendf,t_tendf,xk
1415 ! mhd
1416 do ij = 1, grid%num_tiles
1417 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
1418 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
1419 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
1420 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
1421 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1422 end do
1423
1424 ! recompute : cqw,mu_tend,ph_save,ph_tend,ru_tend,ru_tendf,rv_tend,rv_t
1425 ! endf,rw_tend,rw_tendf,t_save,t_tend,t_tendf,u_save,v_save,w_save
1426 do ij = 1, grid%num_tiles
1427 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1428 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
1429 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,&
1430 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1431 endif
1432 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
1433 &ph_save,t_save,rk_step,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
1434 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1435 if (config_flags%specified .or. config_flags%nested) then
1436 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
1437 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
1438 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1439 endif
1440 end do
1441
1442 ! recompute : mu_tend,ph_tend,ru_tend,rv_tend,rw_tend,t_tend
1443 do ij = 1, grid%num_tiles
1444 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
1445 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1446 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1447 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
1448 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1449 if (non_hydrostatic) then
1450 call calc_coef_w( a,alpha,gamma,mut,cqw,rdn,rdnw,c2a,dts,g,epssm,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
1451 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1452 endif
1453 end do
1454
1455 !------------------------------
1456 !------------------------------
1457
1458 small_stepzw: do iteration = 1, number_of_small_timesteps
1459 do ij = 1, grid%num_tiles
1460 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
1461 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1462 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1463 if (config_flags%specified .or. config_flags%nested) then
1464 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1465 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1466 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1467 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1468 endif
1469 end do
1470 do ij = 1, grid%num_tiles
1471 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
1472 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1473 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1474 if (config_flags%specified .or. config_flags%nested) then
1475 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1476 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1477 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
1478 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
1479 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
1480 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
1481 endif
1482 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration,number_of_small_timesteps,ids,&
1483 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1484 if (non_hydrostatic) then
1485 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,alb,&
1486 &a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
1487 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1488 endif
1489 if (config_flags%specified .or. config_flags%nested) then
1490 if (non_hydrostatic) then
1491 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,&
1492 &ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1493 if (config_flags%specified) then
1494 call zero_grad_bdy( w_2,'w',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
1495 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1496 else
1497 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
1498 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1499 endif
1500 endif
1501 endif
1502 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,iteration,ids,ide,jds,jde,&
1503 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1504 end do
1505 end do small_stepzw
1506 ! recompute : mu_2,mudf,muts,ph_2,pm1,ru_m,rv_m,t_2,u_2,v_2,w_2,ww_m
1507
1508 do ij = 1, grid%num_tiles
1509 call calc_mu_uv_1( config_flags,muts,muus,muvs,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1510 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1511 ! call small_step_finish( u_2,u_1,v_2,v_1,w_2,w_1,t_2,t_1,ph_2,ph_1,ww,ww1,mu_2,mu_1,mut,muts,muu,muus,muv,muvs,u_save,v_save,&
1512 !&w_save,t_save,ph_save,mu_save,msfu,msfv,msft,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1513 !&grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1514 !-----The original
1515 CALL small_step_finish( u_2, u_1, v_2, v_1, w_2, w_1, &
1516 t_2, t_1, ph_2, ph_1, ww, ww1, &
1517 mu_2, mu_1, &
1518 mut, muts, muu, muus, muv, muvs, &
1519 u_save, v_save, w_save, &
1520 t_save, ph_save, mu_save, &
1521 msfu, msfv, msft, &
1522 h_diabatic, &
1523 number_of_small_timesteps,dts_rk, &
1524 ids, ide, jds, jde, kds, kde, &
1525 ims, ime, jms, jme, kms, kme, &
1526 grid%i_start(ij), grid%i_end(ij), &
1527 grid%j_start(ij), grid%j_end(ij), &
1528 k_start , k_end )
1529 end do
1530 ! recompute : mu_2,ph_2,t_2
1531 !---------------------------------------------------------------------------------
1532
1533 a_moist_scalar_advance: if (num_3d_m .ge. param_first_scalar) then
1534
1535 !----------------------------------------
1536 !1111
1537 !saved by Zaizhong Ma
1538 !----------------------------------------
1539 u_1_keep4(:,:,:) = u_1(:,:,:)
1540 u_2_keep4(:,:,:) = u_2(:,:,:)
1541 v_1_keep4(:,:,:) = v_1(:,:,:)
1542 v_2_keep4(:,:,:) = v_2(:,:,:)
1543 w_1_keep4(:,:,:) = w_1(:,:,:)
1544 w_2_keep4(:,:,:) = w_2(:,:,:)
1545 ph_1_keep4(:,:,:) = ph_1(:,:,:)
1546 ph_2_keep4(:,:,:) = ph_2(:,:,:)
1547 t_1_keep4(:,:,:) = t_1(:,:,:)
1548 t_2_keep4(:,:,:) = t_2(:,:,:)
1549 mu_1_keep4(:,:) = mu_1(:,:)
1550 mu_2_keep4(:,:) = mu_2(:,:)
1551 p_keep4(:,:,:) = p(:,:,:)
1552 al_keep4(:,:,:) = al(:,:,:)
1553 moist_keep4(:,:,:,:) = moist(:,:,:,:)
1554 !----------------------------------------
1555 rv_mz(:,:,:) = rv_m(:,:,:)
1556 ru_tendfz(:,:,:) = ru_tendf(:,:,:)
1557 ru_mz(:,:,:) = ru_m(:,:,:)
1558 moist_tendz(:,:,:,:) = moist_tend(:,:,:,:)
1559 !----------------------------------------
1560
1561 do im = num_3d_m, param_first_scalar, -1
1562 !----------------------------------------
1563 ! 1111
1564 !recovered by Zaizhong Ma
1565 !----------------------------------------
1566 u_1(:,:,:) = u_1_keep4(:,:,:)
1567 u_2(:,:,:) = u_2_keep4(:,:,:)
1568 v_1(:,:,:) = v_1_keep4(:,:,:)
1569 v_2(:,:,:) = v_2_keep4(:,:,:)
1570 w_1(:,:,:) = w_1_keep4(:,:,:)
1571 w_2(:,:,:) = w_2_keep4(:,:,:)
1572 ph_1(:,:,:) = ph_1_keep4(:,:,:)
1573 ph_2(:,:,:) = ph_2_keep4(:,:,:)
1574 t_1(:,:,:) = t_1_keep4(:,:,:)
1575 t_2(:,:,:) = t_2_keep4(:,:,:)
1576 mu_1(:,:) = mu_1_keep4(:,:)
1577 mu_2(:,:) = mu_2_keep4(:,:)
1578 p(:,:,:) = p_keep4(:,:,:)
1579 al(:,:,:) = al_keep4(:,:,:)
1580 moist(:,:,:,:) = moist_keep4(:,:,:,:)
1581 !----------------------------------------
1582 !----------------------------------------
1583 rv_m(:,:,:) = rv_mz(:,:,:)
1584 ru_tendf(:,:,:) = ru_tendfz(:,:,:)
1585 ru_m(:,:,:) = ru_mz(:,:,:)
1586 moist_tend(:,:,:,:) = moist_tendz(:,:,:,:)
1587 !----------------------------------------
1588
1589 moist_variable_loor2: do iy = param_first_scalar, num_3d_m - 1
1590
1591 moist_tile_loop_10: do ij = 1, grid%num_tiles
1592 call rk_scalar_tend( iy,iy,config_flags,rk_step,dt_rk,ru_m,rv_m,ww_m,mut,alt,moist(ims,kms,jms,&
1593 &iy),moist_tend(ims,kms,jms,iy),advect_tend,rqvften,qv_base, .true. ,fnm,fnp,msfu,msfv,msft,rdx,rdy,rdn,rdnw,khdif,kvdif,&
1594 &xkmhd,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
1595 &j_end(ij),k_start,k_end )
1596 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1597 if (iy .eq. p_qv) then
1598 call spec_bdy_scalar( moist_tend(ims,kms,jms,iy),rqv_b,rqv_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1599 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1600 endif
1601 endif
1602 if (config_flags%nested .and. rk_step .eq. 1) then
1603 if (iy .eq. p_qc) then
1604 call spec_bdy_scalar( moist_tend(ims,kms,jms,iy),rqc_b,rqc_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1605 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1606 else if (iy .eq. p_qr) then
1607 call spec_bdy_scalar( moist_tend(ims,kms,jms,iy),rqr_b,rqr_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1608 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1609 else if (iy .eq. p_qi) then
1610 call spec_bdy_scalar( moist_tend(ims,kms,jms,iy),rqi_b,rqi_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1611 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1612 else if (iy .eq. p_qs) then
1613 call spec_bdy_scalar( moist_tend(ims,kms,jms,iy),rqs_b,rqs_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1614 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1615 else if (iy .eq. p_qg) then
1616 call spec_bdy_scalar( moist_tend(ims,kms,jms,iy),rqg_b,rqg_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1617 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1618 endif
1619 endif
1620 end do moist_tile_loop_10
1621
1622 moist_tile_loop_11: do ij = 1, grid%num_tiles
1623 call rk_update_scalar( iy,iy,moist_old(ims,kms,jms,iy),moist(ims,kms,jms,iy),moist_tend(ims,kms,jms,iy),advect_tend,msft,&
1624 &mu_1,mu_2,mub,rk_step,dt_rk,spec_zone,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
1625 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1626 end do moist_tile_loop_11
1627
1628 end do moist_variable_loor2
1629
1630 moist_tile_loop_12: do ij = 1, grid%num_tiles
1631 call rk_scalar_tend( im,im,config_flags,rk_step,dt_rk,ru_m,rv_m,ww_m,mut,alt,moist(ims,kms,jms,&
1632 &im),moist_tend(ims,kms,jms,im),advect_tend,rqvften,qv_base, .true. ,fnm,fnp,msfu,msfv,msft,rdx,rdy,rdn,rdnw,khdif,kvdif,&
1633 &xkmhd,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
1634 &j_end(ij),k_start,k_end )
1635 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1636 if (im .eq. p_qv) then
1637 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqv_b,rqv_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1638 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1639 endif
1640 endif
1641 if (config_flags%nested .and. rk_step .eq. 1) then
1642 if (im .eq. p_qc) then
1643 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqc_b,rqc_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1644 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1645 else if (im .eq. p_qr) then
1646 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqr_b,rqr_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1647 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1648 else if (im .eq. p_qi) then
1649 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqi_b,rqi_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1650 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1651 else if (im .eq. p_qs) then
1652 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqs_b,rqs_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1653 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1654 else if (im .eq. p_qg) then
1655 call spec_bdy_scalar( moist_tend(ims,kms,jms,im),rqg_b,rqg_bt,spec_bdy_width,spec_zone,config_flags,ijds,ijde,ids,ide,jds,jde,kds,&
1656 &kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1657 endif
1658 endif
1659 end do moist_tile_loop_12
1660
1661 !-------------------------------------------------------
1662
1663 ! recompute : advect_tend,moist_tend
1664 do ij = 1, grid%num_tiles
1665 call a_rk_update_scalar( im,im,moist_old(ims,kms,jms,im),a_moist_old(ims,kms,jms,im),moist(ims,kms,jms,im),a_moist(ims,kms,&
1666 &jms,im),moist_tend(ims,kms,jms,im),a_moist_tend(ims,kms,jms,im),advect_tend,a_advect_tend,msft,mu_1,a_mu_1,mu_2,a_mu_2,&
1667 &mub,rk_step,dt_rk,spec_zone,config_flags,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1668 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1669 end do
1670
1671 do ij = 1, grid%num_tiles
1672 if (config_flags%nested .and. rk_step .eq. 1) then
1673 if (im .eq. p_qc) then
1674 call a_spec_bdy_scalar( a_moist_tend(ims,kms,jms,im),a_rqc_bth,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,&
1675 &kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1676 else if (im .eq. p_qr) then
1677 call a_spec_bdy_scalar( a_moist_tend(ims,kms,jms,im),a_rqr_bth,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,&
1678 &kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1679 else if (im .eq. p_qi) then
1680 call a_spec_bdy_scalar( a_moist_tend(ims,kms,jms,im),a_rqi_bth,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,&
1681 &kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1682 else if (im .eq. p_qs) then
1683 call a_spec_bdy_scalar( a_moist_tend(ims,kms,jms,im),a_rqs_bth,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,&
1684 &kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1685 else if (im .eq. p_qg) then
1686 call a_spec_bdy_scalar( a_moist_tend(ims,kms,jms,im),a_rqg_bth,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,&
1687 &kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1688 endif
1689 endif
1690 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1691 if (im .eq. p_qv) then
1692 call a_spec_bdy_scalar( a_moist_tend(ims,kms,jms,im),a_rqv_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,&
1693 &kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1694 endif
1695 endif
1696 call a_rk_scalar_tend( im,im,config_flags,rk_step,ru_m,a_ru_m,rv_m,a_rv_m,ww_m,a_ww_m,mut,a_mut,alt,a_alt,moist(ims,kms,jms,im),a_moist(ims,kms,jms,im),a_moist_tend(ims,kms,jms,im),&
1697 &a_advect_tend,qv_base, .true. ,fnm,fnp,msfu,msfv,msft,rdx,rdy,rdn,rdnw,kvdif,xkmhd,a_xkmhd,ids,ide,jds,jde,kde,&
1698 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1699 end do
1700 end do
1701 endif a_moist_scalar_advance
1702
1703 !------------------------------
1704 !------------------------------
1705
1706 !------------------------------
1707 ww_m(:,:,:) = ww_m_keep3(:,:,:)
1708 rv_m(:,:,:) = rv_m_keep3(:,:,:)
1709 ru_m(:,:,:) = ru_m_keep3(:,:,:)
1710 a(:,:,:) = a_keep3(:,:,:)
1711 alpha(:,:,:) = alpha_keep3(:,:,:)
1712 gamma(:,:,:) = gamma_keep3(:,:,:)
1713 !------------------------------
1714
1715 xkmhd(:,:,:) = xkmhdh(:,:,:)
1716 w_2(:,:,:) = w_2h(:,:,:)
1717 w_1(:,:,:) = w_1h(:,:,:)
1718 v_2(:,:,:) = v_2h(:,:,:)
1719 v_1(:,:,:) = v_1h(:,:,:)
1720 u_2(:,:,:) = u_2h(:,:,:)
1721 u_1(:,:,:) = u_1h(:,:,:)
1722 t_tendf(:,:,:) = t_tendfh(:,:,:)
1723 t_2(:,:,:) = t_2h(:,:,:)
1724 t_1(:,:,:) = t_1h(:,:,:)
1725 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
1726 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
1727 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
1728 pm1(:,:,:) = pm1h(:,:,:)
1729 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
1730 ph_2(:,:,:) = ph_2h(:,:,:)
1731 ph_1(:,:,:) = ph_1h(:,:,:)
1732 p(:,:,:) = ph(:,:,:)
1733 mu_2(:,:) = mu_2h(:,:)
1734 mu_1(:,:) = mu_1h(:,:)
1735 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
1736 moist(:,:,:,:) = moist_h(:,:,:,:)
1737 dnw(:) = dnwh(:)
1738 al(:,:,:) = alh(:,:,:)
1739 !----------------------------------------
1740 mudf(:,:) = mudf_ma(:,:)
1741 !----------------------------------------
1742 if (rk_order .eq. 1) then
1743 if (step_number .ne. 1) then
1744 number_of_small_timesteps = 2*time_step_sound
1745 dt_rk = dt
1746 else
1747 number_of_small_timesteps = time_step_sound
1748 dt_rk = dt/2.
1749 endif
1750 dts_rk = dts
1751 else if (rk_order .eq. 2) then
1752 if (rk_step .eq. 1) then
1753 dt_rk = 0.5*dt
1754 dts_rk = dts
1755 number_of_small_timesteps = time_step_sound/2
1756 else
1757 dt_rk = dt
1758 dts_rk = dts
1759 number_of_small_timesteps = time_step_sound
1760 endif
1761 else if (rk_order .eq. 3) then
1762 if (rk_step .eq. 1) then
1763 dt_rk = dt/3.
1764 dts_rk = dt_rk
1765 number_of_small_timesteps = 1
1766 else if (rk_step .eq. 2) then
1767 dt_rk = 0.5*dt
1768 dts_rk = dts
1769 number_of_small_timesteps = time_step_sound/2
1770 else
1771 dt_rk = dt
1772 dts_rk = dts
1773 number_of_small_timesteps = time_step_sound
1774 endif
1775 endif
1776 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
1777 do ij = 1, grid%num_tiles
1778 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
1779 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1780 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1781 end do
1782 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
1783 rk_step_is_onzz: if (rk_step .eq. 1) then
1784 do ij = 1, grid%num_tiles
1785 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,rk_step,&
1786 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1787 end do
1788 do ij = 1, grid%num_tiles
1789 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
1790 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
1791 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1792 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1793 end do
1794 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
1795 do ij = 1, grid%num_tiles
1796 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
1797 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
1798 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1799 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1800 end do
1801 endif
1802 endif rk_step_is_onzz
1803 ! recompute : moist_tend,ph_tendf,ru_tendf,rv_tendf,rw_tendf,t_tendf,xk
1804 ! mhd
1805 do ij = 1, grid%num_tiles
1806 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
1807 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
1808 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
1809 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
1810 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1811 end do
1812
1813 ! recompute : cqw,mu_tend,ph_save,ph_tend,ru_tend,ru_tendf,rv_tend,rv_t
1814 ! endf,rw_tend,rw_tendf,t_save,t_tend,t_tendf,u_save,v_save,w_save
1815 do ij = 1, grid%num_tiles
1816 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
1817 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
1818 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,&
1819 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1820 endif
1821 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
1822 &ph_save,t_save,rk_step,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
1823 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1824 if (config_flags%specified .or. config_flags%nested) then
1825 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
1826 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
1827 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1828 endif
1829 end do
1830
1831 ! recompute : mu_tend,ph_tend,ru_tend,rv_tend,rw_tend,t_tend
1832 do ij = 1, grid%num_tiles
1833 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
1834 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1835 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1836 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
1837 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1838 if (non_hydrostatic) then
1839 call calc_coef_w( a,alpha,gamma,mut,cqw,rdn,rdnw,c2a,dts,g,epssm,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
1840 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1841 endif
1842 end do
1843 ! recompute : a,al,alpha,c2a,gamma,mu_1,mu_2,mu_save,mudf,p,ph_1,ph_2,p
1844 ! h_save,pm1,t_1,t_2,t_save,u_1,u_2,u_save,v_1,v_2,v_save,w_1,w_2,w_save
1845 ! ,ww1
1846
1847 !------------------------------
1848 !------------------------------
1849
1850 small_stepzz: do iteration = 1, number_of_small_timesteps
1851 do ij = 1, grid%num_tiles
1852 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
1853 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
1854 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1855 if (config_flags%specified .or. config_flags%nested) then
1856 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1857 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1858 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1859 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1860 endif
1861 end do
1862 do ij = 1, grid%num_tiles
1863 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
1864 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
1865 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1866 if (config_flags%specified .or. config_flags%nested) then
1867 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
1868 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1869 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
1870 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
1871 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
1872 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
1873 endif
1874 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration,number_of_small_timesteps,ids,&
1875 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1876 if (non_hydrostatic) then
1877 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,alb,&
1878 &a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
1879 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1880 endif
1881 if (config_flags%specified .or. config_flags%nested) then
1882 if (non_hydrostatic) then
1883 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,&
1884 &ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1885 if (config_flags%specified) then
1886 call zero_grad_bdy( w_2,'w',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,&
1887 &kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1888 else
1889 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
1890 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1891 endif
1892 endif
1893 endif
1894 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,iteration,ids,ide,jds,jde,&
1895 &kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1896 end do
1897 end do small_stepzz
1898
1899 !------------------------------------------------------
1900
1901 do ij = 1, grid%num_tiles
1902
1903 call calc_mu_uv_1( config_flags,muts,muus,muvs,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
1904 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
1905
1906 ! recompute : muus,muvs
1907 call 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,a_muu,muus,a_muus,&
1908 &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,a_mu_save,msfu,msfv,msft,ide,&
1909 &jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij) )
1910
1911
1912 ! call 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_ww,a_ww1,a_mu_2,mut,a_mut,muts,a_muts,muu,a_muu,muus,a_muus,&
1913 !&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,a_mu_save,msfu,msfv,msft,h_diabatic,number_of_small_timesteps,dts_rk,ide,&
1914 !&jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij) )
1915
1916 call a_calc_mu_uv_1( config_flags,a_muts,a_muus,a_muvs,ids,ide,jds,jde,ims,ime,jms,jme,grid%i_start(ij),grid%i_end(ij),grid%&
1917 &j_start(ij),grid%j_end(ij) )
1918 end do
1919
1920
1921 !------------------------------
1922 !------------------------------
1923
1924 a_small_steps: do iteration = number_of_small_timesteps, 1, -1
1925
1926
1927 !------------------------------
1928 ww_m(:,:,:) = ww_m_keep3(:,:,:)
1929 rv_m(:,:,:) = rv_m_keep3(:,:,:)
1930 ru_m(:,:,:) = ru_m_keep3(:,:,:)
1931 a(:,:,:) = a_keep3(:,:,:)
1932 alpha(:,:,:) = alpha_keep3(:,:,:)
1933 gamma(:,:,:) = gamma_keep3(:,:,:)
1934 !------------------------------
1935
1936 xkmhd(:,:,:) = xkmhdh(:,:,:)
1937 w_2(:,:,:) = w_2h(:,:,:)
1938 w_1(:,:,:) = w_1h(:,:,:)
1939 v_2(:,:,:) = v_2h(:,:,:)
1940 v_1(:,:,:) = v_1h(:,:,:)
1941 u_2(:,:,:) = u_2h(:,:,:)
1942 u_1(:,:,:) = u_1h(:,:,:)
1943 t_tendf(:,:,:) = t_tendfh(:,:,:)
1944 t_2(:,:,:) = t_2h(:,:,:)
1945 t_1(:,:,:) = t_1h(:,:,:)
1946 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
1947 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
1948 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
1949 pm1(:,:,:) = pm1h(:,:,:)
1950 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
1951 ph_2(:,:,:) = ph_2h(:,:,:)
1952 ph_1(:,:,:) = ph_1h(:,:,:)
1953 p(:,:,:) = ph(:,:,:)
1954 mu_2(:,:) = mu_2h(:,:)
1955 mu_1(:,:) = mu_1h(:,:)
1956 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
1957 moist(:,:,:,:) = moist_h(:,:,:,:)
1958 dnw(:) = dnwh(:)
1959 al(:,:,:) = alh(:,:,:)
1960 !----------------------------------------
1961 mudf(:,:) = mudf_ma(:,:)
1962 !----------------------------------------
1963 if (rk_order .eq. 1) then
1964 if (step_number .ne. 1) then
1965 number_of_small_timesteps = 2*time_step_sound
1966 dt_rk = dt
1967 else
1968 number_of_small_timesteps = time_step_sound
1969 dt_rk = dt/2.
1970 endif
1971 dts_rk = dts
1972 else if (rk_order .eq. 2) then
1973 if (rk_step .eq. 1) then
1974 dt_rk = 0.5*dt
1975 dts_rk = dts
1976 number_of_small_timesteps = time_step_sound/2
1977 else
1978 dt_rk = dt
1979 dts_rk = dts
1980 number_of_small_timesteps = time_step_sound
1981 endif
1982 else if (rk_order .eq. 3) then
1983 if (rk_step .eq. 1) then
1984 dt_rk = dt/3.
1985 dts_rk = dt_rk
1986 number_of_small_timesteps = 1
1987 else if (rk_step .eq. 2) then
1988 dt_rk = 0.5*dt
1989 dts_rk = dts
1990 number_of_small_timesteps = time_step_sound/2
1991 else
1992 dt_rk = dt
1993 dts_rk = dts
1994 number_of_small_timesteps = time_step_sound
1995 endif
1996 endif
1997 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
1998 do ij = 1, grid%num_tiles
1999 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
2000 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2001 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2002 end do
2003
2004 rk_step_is_onzd: if (rk_step .eq. 1) then
2005 do ij = 1, grid%num_tiles
2006 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,rk_step,&
2007 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2008 end do
2009 do ij = 1, grid%num_tiles
2010 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
2011 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
2012 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
2013 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2014 end do
2015 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
2016 do ij = 1, grid%num_tiles
2017 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
2018 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
2019 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2020 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2021 end do
2022 endif
2023 endif rk_step_is_onzd
2024
2025 do ij = 1, grid%num_tiles
2026 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
2027 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
2028 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
2029 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
2030 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2031 end do
2032 do ij = 1, grid%num_tiles
2033 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
2034 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
2035 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,&
2036 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2037 endif
2038 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
2039 &ph_save,t_save,rk_step,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
2040 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2041 if (config_flags%specified .or. config_flags%nested) then
2042 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
2043 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
2044 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2045 endif
2046 end do
2047
2048 do ij = 1, grid%num_tiles
2049 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
2050 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2051 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2052 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
2053 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2054 if (non_hydrostatic) then
2055 call calc_coef_w( a,alpha,gamma,mut,cqw,rdn,rdnw,c2a,dts,g,epssm,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%&
2056 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2057 endif
2058 end do
2059
2060 small_stepzr: do iteration1 = 1, iteration-1
2061 do ij = 1, grid%num_tiles
2062 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
2063 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2064 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2065 if (config_flags%specified .or. config_flags%nested) then
2066 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
2067 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2068 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
2069 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2070 endif
2071 end do
2072 do ij = 1, grid%num_tiles
2073 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
2074 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration1,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2075 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2076 if (config_flags%specified .or. config_flags%nested) then
2077 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,&
2078 &ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2079 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
2080 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
2081 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
2082 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
2083 endif
2084
2085 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration1,number_of_small_timesteps,ids,&
2086 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2087
2088 if (non_hydrostatic) then
2089 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,&
2090 &alb,a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,&
2091 &ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2092 endif
2093 if (config_flags%specified .or. config_flags%nested) then
2094 if (non_hydrostatic) then
2095 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,&
2096 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2097 if (config_flags%specified) then
2098
2099 CALL zero_grad_bdy ( w_2, &
2100 'w' , config_flags, &
2101 spec_zone, &
2102 ids,ide, jds,jde, kds,kde, & ! domain dims
2103 ims,ime, jms,jme, kms,kme, & ! memory dims
2104 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2105 grid%i_start(ij), grid%i_end(ij), &
2106 grid%j_start(ij), grid%j_end(ij), &
2107 k_start , k_end )
2108
2109 else
2110 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,&
2111 &ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2112 endif
2113 endif
2114 endif
2115 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,iteration1,ids,ide,jds,&
2116 &jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2117 end do
2118 end do small_stepzr
2119
2120 !----------------------------------------
2121 alq(:,:,:) = al(:,:,:)
2122 mu_2q(:,:) = mu_2(:,:)
2123 pq(:,:,:) = p(:,:,:)
2124 ph_2q(:,:,:) = ph_2(:,:,:)
2125
2126 do ij = 1, grid%num_tiles
2127 call advance_uv( u_2,ru_tend,v_2,rv_tend,p,pb,ph_2,php,alt,al,mu_2,muu,cqu,muv,cqv,mudf,rdx,rdy,dts,cf1,cf2,cf3,fnm,fnp,&
2128 &emdiv,rdnw,config_flags,spec_zone,non_hydrostatic,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2129 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2130 if (config_flags%specified .or. config_flags%nested) then
2131 call spec_bdyupdate( u_2,ru_tend,dts_rk,'u',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
2132 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2133 call spec_bdyupdate( v_2,rv_tend,dts_rk,'v',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
2134 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2135 endif
2136 end do
2137
2138 do ij = 1, grid%num_tiles
2139
2140 !----------------------------------------
2141 wwr(:,:,:) = ww(:,:,:)
2142 u_2r(:,:,:) = u_2(:,:,:)
2143 v_2r(:,:,:) = v_2(:,:,:)
2144
2145
2146 call advance_mu_t( ww,ww1,u_2,u_save,v_2,v_save,mu_2,mut,muave,muts,muu,muv,mudf,ru_m,rv_m,ww_m,t_2,t_save,t_2save,t_tend,&
2147 &mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,msfv,msft,iteration,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2148 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2149
2150 if (config_flags%specified .or. config_flags%nested) then
2151 call spec_bdyupdate( t_2,t_tend,dts_rk,'t',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,&
2152 &jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2153 call spec_bdyupdate( mu_2,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
2154 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
2155 call spec_bdyupdate( muts,mu_tend,dts_rk,'m',config_flags,spec_zone,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,ips,ipe,jps,&
2156 &jpe,1,1,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
2157 endif
2158 call sumflux( u_2,v_2,ww,u_save,v_save,ww1,muu,muv,ru_m,rv_m,ww_m,epssm,msfu,msfv,iteration,number_of_small_timesteps,ids,&
2159 &ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2160
2161
2162 !----------------------------------------
2163 w_2s(:,:,:) = w_2(:,:,:)
2164 wws(:,:,:) = ww(:,:,:)
2165 u_2s(:,:,:) = u_2(:,:,:)
2166 v_2s(:,:,:) = v_2(:,:,:)
2167 mu_2s(:,:) = mu_2(:,:)
2168 t_2s(:,:,:) = t_2(:,:,:)
2169 ph_2s(:,:,:) = ph_2(:,:,:)
2170 as(:,:,:) = a(:,:,:)
2171 t_2saves(:,:,:) = t_2save(:,:,:)
2172
2173 if (non_hydrostatic) then
2174 call advance_w( w_2,rw_tend,ww,u_2,v_2,mu_2,mut,muave,muts,t_2save,t_2,t_save,ph_2,ph_save,phb,ph_tend,ht,c2a,cqw,alt,alb,&
2175 &a,alpha,gamma,rdx,rdy,dts,t0,epssm,dnw,fnm,fnp,rdnw,rdn,cf1,cf2,cf3,msft,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,&
2176 &jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2177 endif
2178 if (config_flags%specified .or. config_flags%nested) then
2179 if (non_hydrostatic) then
2180
2181 call spec_bdyupdate_ph( ph_save,ph_2,ph_tend,mu_tend,muts,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,&
2182 &ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2183 if (config_flags%specified) then
2184 call zero_grad_bdy( w_2,'w',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,&
2185 &kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2186 else
2187 call spec_bdyupdate( w_2,rw_tend,dts_rk,'h',config_flags,spec_zone,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,&
2188 &ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2189 endif
2190 endif
2191 endif
2192 call a_calc_p_rho( al,a_al,p,a_p,ph_2,a_ph_2,alt,a_alt,t_2,a_t_2,t_save,a_t_save,c2a,a_c2a,a_pm1,mu_2,a_mu_2,muts,a_muts,znu,&
2193 &t0,rdnw,dnw,smdiv,non_hydrostatic,iteration,ide,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%&
2194 &j_start(ij),grid%j_end(ij),k_start,k_end )
2195
2196 if (config_flags%specified .or. config_flags%nested) then
2197 if (non_hydrostatic) then
2198 if (config_flags%specified) then
2199 call a_zero_grad_bdy( a_w_2,'w',spec_zone,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
2200 &grid%j_start(ij),grid%j_end(ij),k_start )
2201 else
2202 call a_spec_bdyupdate( a_w_2,a_rw_tend,dts_rk,'h',spec_zone,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%&
2203 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2204 endif
2205
2206 call a_spec_bdyupdate_ph( ph_save,a_ph_save,ph_2,a_ph_2,ph_tend,a_ph_tend,mu_tend,a_mu_tend,muts,a_muts,dts_rk,'h',&
2207 &spec_zone,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),&
2208 &k_start,k_end )
2209 endif
2210 endif
2211
2212 !----------------------------------------
2213 w_2(:,:,:) = w_2s(:,:,:)
2214 ww(:,:,:) = wws(:,:,:)
2215 u_2(:,:,:) = u_2s(:,:,:)
2216 v_2(:,:,:) = v_2s(:,:,:)
2217 mu_2(:,:) = mu_2s(:,:)
2218 t_2(:,:,:) = t_2s(:,:,:)
2219 ph_2(:,:,:) = ph_2s(:,:,:)
2220 a(:,:,:) = as(:,:,:)
2221 t_2save(:,:,:) = t_2saves(:,:,:)
2222
2223
2224 if (non_hydrostatic) then
2225 call a_advance_w( w_2,a_w_2,rw_tend,a_rw_tend,ww,a_ww,u_2,a_u_2,v_2,a_v_2,mu_2,a_mu_2,mut,a_mut,muave,a_muave,muts,a_muts,&
2226 &t_2save,a_t_2save,t_2,a_t_2,t_save,a_t_save,ph_2,a_ph_2,ph_save,a_ph_save,phb,ph_tend,a_ph_tend,ht,c2a,a_c2a,cqw,a_cqw,&
2227 &alt,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,config_flags,ids,&
2228 &ide,jds,jde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2229 endif
2230
2231 call a_sumflux( a_u_2,a_v_2,a_ww,u_save,a_u_save,v_save,a_v_save,a_ww1,muu,a_muu,muv,a_muv,a_ru_m,a_rv_m,a_ww_m,msfu,msfv,&
2232 &iteration,number_of_small_timesteps,ide,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),&
2233 &grid%j_end(ij),k_start,k_end )
2234
2235 if (config_flags%specified .or. config_flags%nested) then
2236 call a_spec_bdyupdate( a_muts,a_mu_tend,dts_rk,'m',spec_zone,ids,ide,jds,jde,1,ims,ime,jms,jme,1,1,grid%i_start(ij),grid%&
2237 &i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
2238 call a_spec_bdyupdate( a_mu_2,a_mu_tend,dts_rk,'m',spec_zone,ids,ide,jds,jde,1,ims,ime,jms,jme,1,1,grid%i_start(ij),grid%&
2239 &i_end(ij),grid%j_start(ij),grid%j_end(ij),1,1 )
2240 call a_spec_bdyupdate( a_t_2,a_t_tend,dts_rk,'t',spec_zone,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),&
2241 &grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2242 endif
2243
2244 !----------------------------------------
2245 u_2(:,:,:) = u_2r(:,:,:)
2246 v_2(:,:,:) = v_2r(:,:,:)
2247 ww(:,:,:) = wwr(:,:,:)
2248
2249 call a_advance_mu_t( ww,a_ww,ww1,a_ww1,u_2,a_u_2,u_save,a_u_save,v_2,a_v_2,v_save,a_v_save,a_mu_2,a_mut,a_muave,a_muts,muu,&
2250 &a_muu,muv,a_muv,a_mudf,a_t_2,t_save,a_t_save,a_t_2save,a_t_tend,mu_tend,a_mu_tend,rdx,rdy,dts,epssm,dnw,fnm,fnp,rdnw,msfu,&
2251 &msfv,msft,config_flags,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
2252 &j_end(ij),k_start,k_end )
2253
2254 end do
2255 !----------------------------------------
2256 al(:,:,:) = alq(:,:,:)
2257 mu_2(:,:) = mu_2q(:,:)
2258 p(:,:,:) = pq(:,:,:)
2259 ph_2(:,:,:) = ph_2q(:,:,:)
2260
2261 do ij = 1, grid%num_tiles
2262 if (config_flags%specified .or. config_flags%nested) then
2263 call a_spec_bdyupdate( a_v_2,a_rv_tend,dts_rk,'v',spec_zone,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),&
2264 &grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2265 call a_spec_bdyupdate( a_u_2,a_ru_tend,dts_rk,'u',spec_zone,ids,ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),&
2266 &grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2267 endif
2268 call a_advance_uv( a_u_2,a_ru_tend,a_v_2,a_rv_tend,p,a_p,pb,ph_2,a_ph_2,php,a_php,alt,a_alt,al,a_al,mu_2,a_mu_2,muu,a_muu,&
2269 &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,non_hydrostatic,ids,&
2270 &ide,jds,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2271
2272 end do
2273 end do a_small_steps
2274
2275 !------------------------------
2276 !------------------------------
2277
2278 ww_m(:,:,:) = ww_m_keep3(:,:,:)
2279 rv_m(:,:,:) = rv_m_keep3(:,:,:)
2280 ru_m(:,:,:) = ru_m_keep3(:,:,:)
2281 a(:,:,:) = a_keep3(:,:,:)
2282 alpha(:,:,:) = alpha_keep3(:,:,:)
2283 gamma(:,:,:) = gamma_keep3(:,:,:)
2284 !------------------------------
2285
2286 xkmhd(:,:,:) = xkmhdh(:,:,:)
2287 w_2(:,:,:) = w_2h(:,:,:)
2288 w_1(:,:,:) = w_1h(:,:,:)
2289 v_2(:,:,:) = v_2h(:,:,:)
2290 v_1(:,:,:) = v_1h(:,:,:)
2291 u_2(:,:,:) = u_2h(:,:,:)
2292 u_1(:,:,:) = u_1h(:,:,:)
2293 t_tendf(:,:,:) = t_tendfh(:,:,:)
2294 t_2(:,:,:) = t_2h(:,:,:)
2295 t_1(:,:,:) = t_1h(:,:,:)
2296 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
2297 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
2298 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
2299 pm1(:,:,:) = pm1h(:,:,:)
2300 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
2301 ph_2(:,:,:) = ph_2h(:,:,:)
2302 ph_1(:,:,:) = ph_1h(:,:,:)
2303 p(:,:,:) = ph(:,:,:)
2304 mu_2(:,:) = mu_2h(:,:)
2305 mu_1(:,:) = mu_1h(:,:)
2306 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
2307 moist(:,:,:,:) = moist_h(:,:,:,:)
2308 dnw(:) = dnwh(:)
2309 al(:,:,:) = alh(:,:,:)
2310 !----------------------------------------
2311 mudf(:,:) = mudf_ma(:,:)
2312 !----------------------------------------
2313 if (rk_order .eq. 1) then
2314 if (step_number .ne. 1) then
2315 number_of_small_timesteps = 2*time_step_sound
2316 dt_rk = dt
2317 else
2318 number_of_small_timesteps = time_step_sound
2319 dt_rk = dt/2.
2320 endif
2321 dts_rk = dts
2322 else if (rk_order .eq. 2) then
2323 if (rk_step .eq. 1) then
2324 dt_rk = 0.5*dt
2325 dts_rk = dts
2326 number_of_small_timesteps = time_step_sound/2
2327 else
2328 dt_rk = dt
2329 dts_rk = dts
2330 number_of_small_timesteps = time_step_sound
2331 endif
2332 else if (rk_order .eq. 3) then
2333 if (rk_step .eq. 1) then
2334 dt_rk = dt/3.
2335 dts_rk = dt_rk
2336 number_of_small_timesteps = 1
2337 else if (rk_step .eq. 2) then
2338 dt_rk = 0.5*dt
2339 dts_rk = dts
2340 number_of_small_timesteps = time_step_sound/2
2341 else
2342 dt_rk = dt
2343 dts_rk = dts
2344 number_of_small_timesteps = time_step_sound
2345 endif
2346 endif
2347 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
2348
2349 do ij = 1, grid%num_tiles
2350 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
2351 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2352 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2353 end do
2354 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
2355
2356 rk_step_is_onzf: if (rk_step .eq. 1) then
2357 do ij = 1, grid%num_tiles
2358 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,rk_step,&
2359 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2360 end do
2361 do ij = 1, grid%num_tiles
2362 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
2363 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
2364 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
2365 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2366 end do
2367 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
2368 do ij = 1, grid%num_tiles
2369 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
2370 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
2371 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2372 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2373 end do
2374 endif
2375 endif rk_step_is_onzf
2376
2377 ! recompute : moist_tend,ph_tendf,ru_tendf,rv_tendf,rw_tendf,t_tendf,xk
2378 ! mhd
2379
2380 do ij = 1, grid%num_tiles
2381 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
2382 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
2383 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
2384 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
2385 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2386 end do
2387
2388 ! recompute : cqw,mu_tend,ph_save,ph_tend,ru_tend,ru_tendf,rv_tend,rv_t
2389 ! endf,rw_tend,rw_tendf,t_save,t_tend,t_tendf,u_save,v_save,w_save
2390
2391
2392 do ij = 1, grid%num_tiles
2393 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
2394 call relax_bdy_dry( config_flags,u_save,v_save,ph_save,t_save,w_save,mu_tend,ru,rv,ph_2,t_2,w_2,mu_2,mut,u_b,v_b,ph_b,t_b,&
2395 &w_b,mu_b,u_bt,v_bt,ph_bt,t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,&
2396 &ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2397 endif
2398 call rk_addtend_dry( ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,u_save,v_save,w_save,&
2399 &ph_save,t_save,rk_step,h_diabatic,mut,msft,msfu,msfv,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
2400 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2401 if (config_flags%specified .or. config_flags%nested) then
2402 call spec_bdy_dry( config_flags,ru_tend,rv_tend,ph_tend,t_tend,rw_tend,mu_tend,u_b,v_b,ph_b,t_b,w_b,mu_b,u_bt,v_bt,ph_bt,&
2403 &t_bt,w_bt,mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,&
2404 &grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2405 endif
2406 end do
2407
2408
2409 !------------------------------
2410 !------------------------------
2411
2412
2413 do ij = 1, grid%num_tiles
2414
2415 !----------------------------------------
2416 u_1o(:,:,:) = u_1(:,:,:)
2417 u_2o(:,:,:) = u_2(:,:,:)
2418 v_1o(:,:,:) = v_1(:,:,:)
2419 v_2o(:,:,:) = v_2(:,:,:)
2420 w_1o(:,:,:) = w_1(:,:,:)
2421 w_2o(:,:,:) = w_2(:,:,:)
2422 t_1o(:,:,:) = t_1(:,:,:)
2423 t_2o(:,:,:) = t_2(:,:,:)
2424 mu_1o(:,:) = mu_1(:,:)
2425 mu_2o(:,:) = mu_2(:,:)
2426 po(:,:,:) = p(:,:,:)
2427
2428 call small_step_prep( u_1,u_2,v_1,v_2,w_1,w_2,t_1,t_2,ph_1,ph_2,mub,mu_1,mu_2,muu,muus,muv,muvs,mut,muts,mudf,u_save,v_save,&
2429 &w_save,t_save,ph_save,mu_save,ww,ww1,dnw,c2a,pb,p,alt,msfu,msfv,msft,rk_step,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2430 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2431
2432 !----------------------------------------
2433 alp(:,:,:) = al(:,:,:)
2434 pp(:,:,:) = p(:,:,:)
2435 ph_2p(:,:,:) = ph_2(:,:,:)
2436 t_2p(:,:,:) = t_2(:,:,:)
2437 mu_2p(:,:) = mu_2(:,:)
2438
2439 call calc_p_rho( al,p,ph_2,alt,t_2,t_save,c2a,pm1,mu_2,muts,znu,t0,rdnw,dnw,smdiv,non_hydrostatic,0,ids,ide,jds,jde,kds,kde,&
2440 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2441
2442 if (non_hydrostatic) then
2443 call 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,ide,jde,kde,ims,ime,&
2444 &jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij) )
2445 endif
2446
2447 !----------------------------------------
2448 al(:,:,:) = alp(:,:,:)
2449 p(:,:,:) = pp(:,:,:)
2450 ph_2(:,:,:) = ph_2p(:,:,:)
2451 t_2(:,:,:) = t_2p(:,:,:)
2452 mu_2(:,:) = mu_2p(:,:)
2453
2454 call a_calc_p_rho( al,a_al,p,a_p,ph_2,a_ph_2,alt,a_alt,t_2,a_t_2,t_save,a_t_save,c2a,a_c2a,a_pm1,mu_2,a_mu_2,muts,a_muts,znu,&
2455 &t0,rdnw,dnw,smdiv,non_hydrostatic,0,ide,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
2456 &j_end(ij),k_start,k_end )
2457
2458 !----------------------------------------
2459 u_1(:,:,:) = u_1o(:,:,:)
2460 u_2(:,:,:) = u_2o(:,:,:)
2461 v_1(:,:,:) = v_1o(:,:,:)
2462 v_2(:,:,:) = v_2o(:,:,:)
2463 w_1(:,:,:) = w_1o(:,:,:)
2464 w_2(:,:,:) = w_2o(:,:,:)
2465 t_1(:,:,:) = t_1o(:,:,:)
2466 t_2(:,:,:) = t_2o(:,:,:)
2467 mu_1(:,:) = mu_1o(:,:)
2468 mu_2(:,:) = mu_2o(:,:)
2469 p(:,:,:) = po(:,:,:)
2470
2471 call 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,a_ph_1,a_ph_2,mub,mu_1,&
2472 &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,a_mudf,a_u_save,a_v_save,a_w_save,&
2473 &a_t_save,a_ph_save,a_mu_save,a_ww,a_ww1,a_c2a,pb,p,a_p,alt,a_alt,msfu,msfv,msft,rk_step,ide,jde,kde,ims,ime,jms,jme,&
2474 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2475
2476 end do
2477
2478 !------------------------------
2479 !------------------------------
2480
2481 xkmhd(:,:,:) = xkmhdh(:,:,:)
2482 w_2(:,:,:) = w_2h(:,:,:)
2483 w_1(:,:,:) = w_1h(:,:,:)
2484 v_2(:,:,:) = v_2h(:,:,:)
2485 v_1(:,:,:) = v_1h(:,:,:)
2486 u_2(:,:,:) = u_2h(:,:,:)
2487 u_1(:,:,:) = u_1h(:,:,:)
2488 t_tendf(:,:,:) = t_tendfh(:,:,:)
2489 t_2(:,:,:) = t_2h(:,:,:)
2490 t_1(:,:,:) = t_1h(:,:,:)
2491 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
2492 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
2493 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
2494 pm1(:,:,:) = pm1h(:,:,:)
2495 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
2496 ph_2(:,:,:) = ph_2h(:,:,:)
2497 ph_1(:,:,:) = ph_1h(:,:,:)
2498 p(:,:,:) = ph(:,:,:)
2499 mu_2(:,:) = mu_2h(:,:)
2500 mu_1(:,:) = mu_1h(:,:)
2501 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
2502 moist(:,:,:,:) = moist_h(:,:,:,:)
2503 dnw(:) = dnwh(:)
2504 al(:,:,:) = alh(:,:,:)
2505 !----------------------------------------
2506 !---------- added by zzma for recomputing basic states -------------
2507 if (rk_order .eq. 1) then
2508 if (step_number .ne. 1) then
2509 number_of_small_timesteps = 2*time_step_sound
2510 dt_rk = dt
2511 else
2512 number_of_small_timesteps = time_step_sound
2513 dt_rk = dt/2.
2514 endif
2515 dts_rk = dts
2516 else if (rk_order .eq. 2) then
2517 if (rk_step .eq. 1) then
2518 dt_rk = 0.5*dt
2519 dts_rk = dts
2520 number_of_small_timesteps = time_step_sound/2
2521 else
2522 dt_rk = dt
2523 dts_rk = dts
2524 number_of_small_timesteps = time_step_sound
2525 endif
2526 else if (rk_order .eq. 3) then
2527 if (rk_step .eq. 1) then
2528 dt_rk = dt/3.
2529 dts_rk = dt_rk
2530 number_of_small_timesteps = 1
2531 else if (rk_step .eq. 2) then
2532 dt_rk = 0.5*dt
2533 dts_rk = dts
2534 number_of_small_timesteps = time_step_sound/2
2535 else
2536 dt_rk = dt
2537 dts_rk = dts
2538 number_of_small_timesteps = time_step_sound
2539 endif
2540 endif
2541 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
2542
2543 do ij = 1, grid%num_tiles
2544 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
2545 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2546 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2547 end do
2548 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
2549
2550 rk_step_is_onzs: if (rk_step .eq. 1) then
2551 do ij = 1, grid%num_tiles
2552 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,rk_step,&
2553 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2554 end do
2555 do ij = 1, grid%num_tiles
2556 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
2557 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
2558 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
2559 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2560 end do
2561 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
2562 do ij = 1, grid%num_tiles
2563 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
2564 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
2565 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2566 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2567 end do
2568 endif
2569 endif rk_step_is_onzs
2570 ! recompute : moist_tend,ph_tendf,ru_tendf,rv_tendf,rw_tendf,t_tendf,xk
2571 ! mhd
2572
2573 do ij = 1, grid%num_tiles
2574 call rk_tendency( config_flags,rk_step,ru_tend,rv_tend,rw_tend,ph_tend,t_tend,ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,&
2575 &mu_tend,u_save,v_save,w_save,ph_save,t_save,mu_save,rthften,ru,rv,rw,ww,u_2,v_2,w_2,t_2,ph_2,u_1,v_1,w_1,t_1,ph_1,h_diabatic,&
2576 &phb,t_init,mu_2,mut,muu,muv,mub,al,alt,p,pb,php,cqu,cqv,cqw,u_base,v_base,t_base,qv_base,z_base,msfu,msfv,msft,f,e,sina,cosa,&
2577 &fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,cf1,cf2,cf3,cfn,cfn1,num_3d_m,non_hydrostatic,ids,ide,jds,jde,kds,kde,&
2578 &ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2579 end do
2580
2581 !------------------------------
2582 !------------------------------
2583
2584 do ij = 1, grid%num_tiles
2585 if (config_flags%specified .or. config_flags%nested) then
2586 call a_spec_bdy_dry( config_flags,a_ru_tend,a_rv_tend,a_ph_tend,a_t_tend,a_rw_tend,a_mu_tend,a_u_bt,a_v_bt,a_ph_bt,a_t_bt,&
2587 &a_w_bt,a_mu_bt,spec_bdy_width,spec_zone,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2588 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2589 endif
2590 call a_rk_addtend_dry( a_ru_tend,a_rv_tend,a_rw_tend,a_ph_tend,a_t_tend,a_ru_tendf,a_rv_tendf,a_rw_tendf,a_ph_tendf,a_t_tendf,&
2591 &a_u_save,a_v_save,a_w_save,a_ph_save,a_t_save,rk_step,h_diabatic,a_mut,msft,msfu,msfv,ide,jde,ims,ime,jms,jme,kms,kme,grid%&
2592 &i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2593 if ((config_flags%specified .or. config_flags%nested) .and. rk_step .eq. 1) then
2594 call a_relax_bdy_dry( config_flags,a_u_save,a_v_save,a_ph_save,a_t_save,a_w_save,a_mu_tend,a_ru,a_rv,ph_2,a_ph_2,t_2,a_t_2,&
2595 &w_2,a_w_2,a_mu_2,mut,a_mut,a_u_b,a_v_b,a_ph_b,a_t_b,a_w_b,a_mu_b,a_u_bt,a_v_bt,a_ph_bt,a_t_bt,a_w_bt,a_mu_bt,spec_bdy_width,&
2596 &spec_zone,relax_zone,dtbc,fcx,gcx,ijds,ijde,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
2597 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2598 endif
2599 end do
2600
2601 !------------------------------
2602 !------------------------------
2603
2604 xkmhd(:,:,:) = xkmhdh(:,:,:)
2605 w_2(:,:,:) = w_2h(:,:,:)
2606 w_1(:,:,:) = w_1h(:,:,:)
2607 v_2(:,:,:) = v_2h(:,:,:)
2608 v_1(:,:,:) = v_1h(:,:,:)
2609 u_2(:,:,:) = u_2h(:,:,:)
2610 u_1(:,:,:) = u_1h(:,:,:)
2611 t_tendf(:,:,:) = t_tendfh(:,:,:)
2612 t_2(:,:,:) = t_2h(:,:,:)
2613 t_1(:,:,:) = t_1h(:,:,:)
2614 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
2615 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
2616 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
2617 pm1(:,:,:) = pm1h(:,:,:)
2618 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
2619 ph_2(:,:,:) = ph_2h(:,:,:)
2620 ph_1(:,:,:) = ph_1h(:,:,:)
2621 p(:,:,:) = ph(:,:,:)
2622 mu_2(:,:) = mu_2h(:,:)
2623 mu_1(:,:) = mu_1h(:,:)
2624 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
2625 moist(:,:,:,:) = moist_h(:,:,:,:)
2626 dnw(:) = dnwh(:)
2627 al(:,:,:) = alh(:,:,:)
2628 !----------------------------------------
2629
2630 if (rk_order .eq. 1) then
2631 if (step_number .ne. 1) then
2632 number_of_small_timesteps = 2*time_step_sound
2633 dt_rk = dt
2634 else
2635 number_of_small_timesteps = time_step_sound
2636 dt_rk = dt/2.
2637 endif
2638 dts_rk = dts
2639 else if (rk_order .eq. 2) then
2640 if (rk_step .eq. 1) then
2641 dt_rk = 0.5*dt
2642 dts_rk = dts
2643 number_of_small_timesteps = time_step_sound/2
2644 else
2645 dt_rk = dt
2646 dts_rk = dts
2647 number_of_small_timesteps = time_step_sound
2648 endif
2649 else if (rk_order .eq. 3) then
2650 if (rk_step .eq. 1) then
2651 dt_rk = dt/3.
2652 dts_rk = dt_rk
2653 number_of_small_timesteps = 1
2654 else if (rk_step .eq. 2) then
2655 dt_rk = 0.5*dt
2656 dts_rk = dts
2657 number_of_small_timesteps = time_step_sound/2
2658 else
2659 dt_rk = dt
2660 dts_rk = dts
2661 number_of_small_timesteps = time_step_sound
2662 endif
2663 endif
2664 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
2665 do ij = 1, grid%num_tiles
2666 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
2667 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2668 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2669 end do
2670 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
2671 rk_step_is_onzr: if (rk_step .eq. 1) then
2672 do ij = 1, grid%num_tiles
2673 call init_zero_tendency( ru_tendf,rv_tendf,rw_tendf,ph_tendf,t_tendf,tke_tend,moist_tend,chem_tend,num_3d_m,num_3d_c,rk_step,&
2674 &ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2675 end do
2676 do ij = 1, grid%num_tiles
2677 call phy_prep( config_flags,mut,u_2,v_2,p,pb,alt,ph_2,phb,t_2,tsk,moist,num_3d_m,mu_3d,rho,th_phy,p_phy,pi_phy,u_phy,v_phy,&
2678 &p8w,t_phy,t8w,z,z_at_w,dz8w,fnm,fnp,rthraten,rthblten,rublten,rvblten,rqvblten,rqcblten,rqiblten,rthcuten,rqvcuten,rqccuten,&
2679 &rqrcuten,rqicuten,rqscuten,rthften,rqvften,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),&
2680 &grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2681 end do
2682 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
2683 do ij = 1, grid%num_tiles
2684 call calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,xkmv,xkhh,xkhv,bn2,khdif,kvdif,div,defor11,&
2685 &defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,t8w,th_phy,t_phy,p_phy,moist,dn,dnw,dx,dy,rdz,rdzw,&
2686 &mix_cr_len,num_3d_m,cf1,cf2,cf3,warm_rain,kh_tke_upper_bound,kv_tke_upper_bound,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,&
2687 &kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2688 end do
2689 endif
2690 endif rk_step_is_onzr
2691 ! recompute : moist_tend,ph_tendf,ru_tendf,rv_tendf,rw_tendf,t_tendf,xk
2692 ! mhd
2693
2694
2695 do ij = 1, grid%num_tiles
2696 call a_rk_tendency( config_flags,rk_step,a_ru_tend,a_rv_tend,a_rw_tend,a_ph_tend,a_t_tend,a_ru_tendf,a_rv_tendf,a_rw_tendf,&
2697 &a_t_tendf,a_mu_tend,a_u_save,a_v_save,a_w_save,a_ph_save,a_t_save,ru,a_ru,rv,a_rv,rw,a_rw,ww,a_ww,u_2,a_u_2,v_2,a_v_2,w_2,&
2698 &a_w_2,t_2,a_t_2,ph_2,a_ph_2,u_1,a_u_1,v_1,a_v_1,w_1,a_w_1,t_1,a_t_1,ph_1,a_ph_1,phb,t_init,mu_2,a_mu_2,mut,a_mut,muu,a_muu,&
2699 &muv,a_muv,mub,al,a_al,alt,a_alt,p,a_p,pb,php,a_php,cqu,a_cqu,cqv,a_cqv,cqw,a_cqw,u_base,v_base,z_base,msfu,msfv,msft,f,e,sina,&
2700 &cosa,fnm,fnp,rdn,rdnw,dt,rdx,rdy,kvdif,xkmhd,a_xkmhd,cf1,cf2,cf3,cfn,cfn1,non_hydrostatic,ids,ide,jds,jde,kde,ims,&
2701 &ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2702 end do
2703
2704
2705 !------------------------------
2706 !------------------------------
2707
2708 xkmhd(:,:,:) = xkmhdh(:,:,:)
2709 w_2(:,:,:) = w_2h(:,:,:)
2710 w_1(:,:,:) = w_1h(:,:,:)
2711 v_2(:,:,:) = v_2h(:,:,:)
2712 v_1(:,:,:) = v_1h(:,:,:)
2713 u_2(:,:,:) = u_2h(:,:,:)
2714 u_1(:,:,:) = u_1h(:,:,:)
2715 t_tendf(:,:,:) = t_tendfh(:,:,:)
2716 t_2(:,:,:) = t_2h(:,:,:)
2717 t_1(:,:,:) = t_1h(:,:,:)
2718 rw_tendf(:,:,:) = rw_tendfh(:,:,:)
2719 rv_tendf(:,:,:) = rv_tendfh(:,:,:)
2720 ru_tendf(:,:,:) = ru_tendfh(:,:,:)
2721 pm1(:,:,:) = pm1h(:,:,:)
2722 ph_tendf(:,:,:) = ph_tendfh(:,:,:)
2723 ph_2(:,:,:) = ph_2h(:,:,:)
2724 ph_1(:,:,:) = ph_1h(:,:,:)
2725 p(:,:,:) = ph(:,:,:)
2726 mu_2(:,:) = mu_2h(:,:)
2727 mu_1(:,:) = mu_1h(:,:)
2728 moist_tend(:,:,:,:) = moist_tendh(:,:,:,:)
2729 moist(:,:,:,:) = moist_h(:,:,:,:)
2730 dnw(:) = dnwh(:)
2731 al(:,:,:) = alh(:,:,:)
2732 !----------------------------------------
2733 if (rk_order .eq. 1) then
2734 if (step_number .ne. 1) then
2735 number_of_small_timesteps = 2*time_step_sound
2736 dt_rk = dt
2737 else
2738 number_of_small_timesteps = time_step_sound
2739 dt_rk = dt/2.
2740 endif
2741 dts_rk = dts
2742 else if (rk_order .eq. 2) then
2743 if (rk_step .eq. 1) then
2744 dt_rk = 0.5*dt
2745 dts_rk = dts
2746 number_of_small_timesteps = time_step_sound/2
2747 else
2748 dt_rk = dt
2749 dts_rk = dts
2750 number_of_small_timesteps = time_step_sound
2751 endif
2752 else if (rk_order .eq. 3) then
2753 if (rk_step .eq. 1) then
2754 dt_rk = dt/3.
2755 dts_rk = dt_rk
2756 number_of_small_timesteps = 1
2757 else if (rk_step .eq. 2) then
2758 dt_rk = 0.5*dt
2759 dts_rk = dts
2760 number_of_small_timesteps = time_step_sound/2
2761 else
2762 dt_rk = dt
2763 dts_rk = dts
2764 number_of_small_timesteps = time_step_sound
2765 endif
2766 endif
2767 ! recompute : dt_rk,dts_rk,number_of_small_timesteps
2768 do ij = 1, grid%num_tiles
2769 call rk_step_prep( config_flags,rk_step,u_2,v_2,w_2,t_2,ph_2,mu_2,moist,ru,rv,rw,ww,php,alt,muu,muv,mub,mut,phb,pb,p,al,alb,&
2770 &cqu,cqv,cqw,msfu,msfv,msft,fnm,fnp,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%&
2771 &i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2772 end do
2773 ! recompute : alt,cqu,cqv,cqw,mut,muu,muv,php,ru,rv,rw,ww
2774
2775 a_rk_step_is_one: if (rk_step .eq. 1) then
2776 ! recompute : p8w,p_phy,t8w,t_phy,th_phy
2777 if (diff_opt .eq. 2 .or. diff_opt .eq. 1) then
2778 do ij = 1, grid%num_tiles
2779 call a_calculate_km_kh( config_flags,dt,dampcoef,zdamp,damp_opt,xkmh,xkmhd,a_xkmhd,xkmv,xkhh,xkhv,bn2,a_bn2h,khdif,div,&
2780 &defor11,defor22,defor33,defor12,defor13,defor23,tke_2(ims,kms,jms),p8w,a_p8w,t8w,a_t8w,th_phy,a_th_phy,t_phy,a_t_phy,&
2781 &p_phy,a_p_phy,moist,a_moist,dn,dnw,dx,dy,rdz,rdzw,num_3d_m,cf1,cf2,cf3,kh_tke_upper_bound,ids,ide,jds,jde,kde,ims,ime,&
2782 &jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2783 end do
2784 endif
2785 do ij = 1, grid%num_tiles
2786 call a_phy_prep( p,a_p,pb,ph_2,a_ph_2,phb,t_2,a_t_2,th_phy,a_th_phy,p_phy,a_p_phy,pi_phy,a_pi_phyh,a_p8w,t_phy,a_t_phy,a_t8w,&
2787 &z,a_z,z_at_w,a_z_at_wh,fnm,fnp,ide,jde,kde,ims,ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%&
2788 &j_end(ij),k_start,k_end )
2789 end do
2790 do ij = 1, grid%num_tiles
2791 call a_init_zero_tendency( a_ru_tendf,a_rv_tendf,a_rw_tendf,a_ph_tendf,a_t_tendf,a_moist_tend,num_3d_m,ims,ime,jms,jme,kms,&
2792 &kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2793 end do
2794 endif a_rk_step_is_one
2795
2796
2797 dnw(:) = dnwh(:)
2798 moist(:,:,:,:) = moist_h(:,:,:,:)
2799 mu_2(:,:) = mu_2h(:,:)
2800 u_2(:,:,:) = u_2h(:,:,:)
2801 v_2(:,:,:) = v_2h(:,:,:)
2802 w_2(:,:,:) = w_2h(:,:,:)
2803
2804 do ij = 1, grid%num_tiles
2805 call a_rk_step_prep( config_flags,u_2,a_u_2,v_2,a_v_2,w_2,a_w_2,a_ph_2,mu_2,a_mu_2,moist,a_moist,a_ru,a_rv,a_rw,a_ww,a_php,&
2806 &a_alt,muu,a_muu,muv,a_muv,mub,mut,a_mut,a_al,a_cqu,a_cqv,a_cqw,msfu,msfv,msft,dnw,rdx,rdy,num_3d_m,ids,ide,jds,jde,kde,ims,&
2807 &ime,jms,jme,kms,kme,grid%i_start(ij),grid%i_end(ij),grid%j_start(ij),grid%j_end(ij),k_start,k_end )
2808 end do
2809 end do a_runge_kutta_loop
2810
2811 !----------------------------------------------
2812 ! FREE DYNAMIC MEMORY
2813 !----------------------------------------------
2814 #endif
2815
2816 end subroutine solve_em_ad
2817