Line data Source code
1 : module dycore_budget
2 : use shr_kind_mod, only: r8=>shr_kind_r8
3 : implicit none
4 :
5 : public :: print_budget
6 : real(r8), parameter :: eps = 1.0E-7_r8
7 : real(r8), parameter :: eps_mass = 1.0E-12_r8
8 :
9 : real(r8), save :: previous_dEdt_adiabatic_dycore = 0.0_r8
10 : real(r8), save :: previous_dEdt_dry_mass_adjust = 0.0_r8
11 : real(r8), save :: previous_dEdt_phys_dyn_coupl_err = 0.0_r8
12 : !=========================================================================================
13 : contains
14 : !=========================================================================================
15 :
16 369408 : subroutine print_budget(hstwr)
17 :
18 : use spmd_utils, only: masterproc
19 : use cam_abortutils, only: endrun
20 : use cam_logfile, only: iulog
21 : use cam_budget, only: cam_budget_get_global, is_cam_budget, thermo_budget_histfile_num, thermo_budget_history
22 : use cam_thermo, only: thermo_budget_vars_descriptor, thermo_budget_num_vars, thermo_budget_vars_massv, &
23 : teidx, seidx, keidx, poidx
24 : use dimensions_mod, only: use_cslam
25 : use control_mod, only: ftype
26 :
27 : ! arguments
28 : logical, intent(in) :: hstwr(:)
29 :
30 : ! Local variables
31 : character(len=*), parameter :: subname = 'dycore_budget:print_budgets:'
32 : !
33 : ! physics energy tendencies
34 : !
35 : integer :: idx(4)
36 : real(r8) :: dEdt_param_physE(4) ! dE/dt CAM physics using physics E formula (phAP-phBP)
37 : real(r8) :: dEdt_param_dynE(4) ! dE/dt CAM physics using dycore E (dyAP-dyBP)
38 :
39 : real(r8) :: dEdt_efix_physE(4) ! dE/dt energy fixer using physics E formula (phBP-phBF)
40 : real(r8) :: dEdt_efix_dynE(4) ! dE/dt energy fixer using dycore E formula (dyBP-dyBF)
41 :
42 : real(r8) :: dEdt_dme_adjust_physE(4) ! dE/dt dry mass adjustment using physics E formula (phAM-phAP)
43 : real(r8) :: dEdt_dme_adjust_dynE(4) ! dE/dt dry mass adjustment using dycore E (dyAM-dyAP)
44 :
45 : real(r8) :: dEdt_param_efix_physE(4) ! dE/dt CAM physics + energy fixer using physics E formula (phAP-phBF)
46 : real(r8) :: dEdt_param_efix_dynE(4) ! dE/dt CAM physics + energy fixer using dycore E formula (dyAP-dyBF)
47 :
48 : real(r8) :: dEdt_phys_total_dynE(4) ! dE/dt physics total using dycore E (dyAM-dyBF)
49 : ! physics total = parameterizations + efix + dry-mass adjustment
50 : !
51 : ! SE dycore specific energy tendencies
52 : !
53 : real(r8) :: dEdt_phys_total_in_dyn(4) ! dEdt of physics total in dynamical core
54 : ! physics total = parameterizations + efix + dry-mass adjustment
55 : real(r8) :: dEdt_dycore_phys ! dEdt dycore (estimated in physics)
56 : !
57 : ! mass budgets physics
58 : !
59 : real(r8) :: dMdt_efix ! mass tendency energy fixer
60 : real(r8) :: dMdt_parameterizations ! mass tendency physics paramterizations
61 : real(r8) :: dMdt_dme_adjust ! mass tendency dry-mass adjustment
62 : real(r8) :: dMdt_phys_total ! mass tendency physics total (energy fixer + parameterizations + dry-mass adjustment)
63 : !
64 : ! mass budgets dynamics
65 : !
66 : real(r8) :: dMdt_floating_dyn ! mass tendency floating dynamics (dAL-dBL)
67 : real(r8) :: dMdt_vert_remap ! mass tendency vertical remapping (dAR-dAD)
68 : real(r8) :: dMdt_del4_fric_heat ! mass tendency del4 frictional heating (dAH-dCH)
69 : real(r8) :: dMdt_del4_tot ! mass tendency del4 + del4 frictional heating (dAH-dBH)
70 : real(r8) :: dMdt_residual ! mass tendency residual (time truncation errors)
71 : real(r8) :: dMdt_phys_total_in_dyn ! mass tendency physics total in dycore
72 : real(r8) :: dMdt_PDC ! mass tendency physics-dynamics coupling
73 : !
74 : ! energy budgets dynamics
75 : !
76 : real(r8) :: dEdt_floating_dyn ! dE/dt floating dynamics (dAL-dBL)
77 : real(r8) :: dEdt_vert_remap ! dE/dt vertical remapping (dAR-dAD)
78 : real(r8) :: dEdt_del4 ! dE/dt del4 (dCH-dBH)
79 : real(r8) :: dEdt_del4_fric_heat ! dE/dt del4 frictional heating (dAH-dCH)
80 : real(r8) :: dEdt_del4_tot ! dE/dt del4 + del4 fricitional heating (dAH-dBH)
81 : real(r8) :: dEdt_del2_sponge ! dE/dt del2 sponge (dAS-dBS)
82 : real(r8) :: dEdt_del2_del4_tot ! dE/dt explicit diffusion total
83 : real(r8) :: dEdt_residual ! dE/dt residual (dEdt_floating_dyn-dEdt_del2_del4_tot)
84 : real(r8) :: dEdt_dycore_dyn ! dE/dt adiabatic dynamical core (calculated in dycore)
85 : !
86 : ! physics-dynamics coupling variables
87 : !
88 : real(r8) :: E_dBF(4) ! E of dynamics state at the end of dycore integration (on dycore deomposition)
89 : real(r8) :: E_dyBF(4) ! E of physics state using dycore E
90 :
91 :
92 : real(r8) :: diff, tmp ! dummy variables
93 : integer :: m_cnst, i
94 : character(LEN=*), parameter :: fmt = "(a40,a15,a1,F6.2,a1,F6.2,a1,E10.2,a5)"
95 : character(LEN=*), parameter :: fmtf = "(a48,F8.4,a6)"
96 : character(LEN=*), parameter :: fmtm = "(a48,E8.2,a9)"
97 : character(LEN=15) :: str(4)
98 : character(LEN=5) :: pf ! pass or fail identifier
99 : !--------------------------------------------------------------------------------------
100 :
101 369408 : if (masterproc .and. thermo_budget_history .and. hstwr(thermo_budget_histfile_num)) then
102 0 : idx(1) = teidx !total energy index
103 0 : idx(2) = seidx !enthaly index
104 0 : idx(3) = keidx !kinetic energy index
105 0 : idx(4) = poidx !surface potential energy index
106 0 : str(1) = "(total )"
107 0 : str(2) = "(enthalpy )"
108 0 : str(3) = "(kinetic )"
109 0 : str(4) = "(srf potential)"
110 0 : do i=1,4
111 : !
112 : ! CAM physics energy tendencies
113 : !
114 0 : call cam_budget_get_global('phAP-phBP',idx(i),dEdt_param_physE(i))
115 : call cam_budget_get_global('phBP-phBF',idx(i),dEdt_efix_physE(i))
116 : call cam_budget_get_global('phAM-phAP',idx(i),dEdt_dme_adjust_physE(i))
117 : call cam_budget_get_global('phAP-phBF',idx(i),dEdt_param_efix_physE(i))
118 : !
119 : ! CAM physics energy tendencies using dycore energy formula scaling
120 : ! temperature tendencies for consistency with CAM physics
121 : !
122 : call cam_budget_get_global('dyAP-dyBP',idx(i),dEdt_param_dynE(i))
123 : call cam_budget_get_global('dyBP-dyBF',idx(i),dEdt_efix_dynE(i))
124 : call cam_budget_get_global('dyAM-dyAP',idx(i),dEdt_dme_adjust_dynE(i))
125 : call cam_budget_get_global('dyAP-dyBF',idx(i),dEdt_param_efix_dynE(i))
126 : call cam_budget_get_global('dyAM-dyBF',idx(i),dEdt_phys_total_dynE(i))
127 : call cam_budget_get_global('dyBF' ,idx(i),E_dyBF(i))!state beginning physics
128 : !
129 : ! CAM physics energy tendencies in dynamical core
130 : !
131 : call cam_budget_get_global('dBD-dAF',idx(i),dEdt_phys_total_in_dyn(i))
132 0 : call cam_budget_get_global('dBF' ,idx(i),E_dBF(i)) !state passed to physics
133 : end do
134 :
135 0 : call cam_budget_get_global('dAL-dBL',teidx,dEdt_floating_dyn)
136 0 : call cam_budget_get_global('dAR-dAD',teidx,dEdt_vert_remap)
137 0 : dEdt_dycore_dyn = dEdt_floating_dyn+dEdt_vert_remap
138 :
139 0 : call cam_budget_get_global('dCH-dBH',teidx,dEdt_del4)
140 0 : call cam_budget_get_global('dAH-dCH',teidx,dEdt_del4_fric_heat)
141 0 : call cam_budget_get_global('dAH-dBH',teidx,dEdt_del4_tot)
142 0 : call cam_budget_get_global('dAS-dBS',teidx,dEdt_del2_sponge)
143 0 : dEdt_del2_del4_tot = dEdt_del4_tot+dEdt_del2_sponge
144 0 : dEdt_residual = dEdt_floating_dyn-dEdt_del2_del4_tot
145 :
146 0 : write(iulog,*)" "
147 0 : write(iulog,*)"======================================================================"
148 0 : write(iulog,*)"Total energy diagnostics introduced in Lauritzen and Williamson (2019)"
149 0 : write(iulog,*)"(DOI:10.1029/2018MS001549)"
150 0 : write(iulog,*)"======================================================================"
151 0 : write(iulog,*)" "
152 0 : write(iulog,*)"Globally and vertically integrated total energy (E) diagnostics are"
153 0 : write(iulog,*)"computed at various points in the physics and dynamics loops to compute"
154 0 : write(iulog,*)"energy tendencies (dE/dt) and check for consistency (e.g., is E of"
155 0 : write(iulog,*)"state passed to physics computed using dycore state variables the same"
156 0 : write(iulog,*)"E of the state in the beginning of physics computed using the physics"
157 0 : write(iulog,*)"representation of the state)"
158 0 : write(iulog,*)" "
159 0 : write(iulog,*)"Energy stages in physics:"
160 0 : write(iulog,*)"-------------------------"
161 0 : write(iulog,*)" "
162 0 : write(iulog,*)" xxBF: state passed to parameterizations, before energy fixer"
163 0 : write(iulog,*)" xxBP: after energy fixer, before parameterizations"
164 0 : write(iulog,*)" xxAP: after last phys_update in parameterizations and state "
165 0 : write(iulog,*)" saved for energy fixer"
166 0 : write(iulog,*)" xxAM: after dry mass adjustment"
167 0 : write(iulog,*)" history files saved off here"
168 0 : write(iulog,*)" "
169 0 : write(iulog,*)"where xx='ph','dy' "
170 0 : write(iulog,*)" "
171 0 : write(iulog,*)"Suffix ph is CAM physics total energy"
172 0 : write(iulog,*)"(eq. 111 in Lauritzen et al. 2022; 10.1029/2022MS003117)"
173 0 : write(iulog,*)" "
174 0 : write(iulog,*)"Suffix dy is dycore energy computed in CAM physics using"
175 0 : write(iulog,*)"CAM physics state variables"
176 0 : write(iulog,*)" "
177 0 : write(iulog,*)" "
178 0 : write(iulog,*)"Energy stages in dynamics (specific to the SE dycore)"
179 0 : write(iulog,*)"-----------------------------------------------------"
180 0 : write(iulog,*)" "
181 0 : write(iulog,*)"suffix (d)"
182 0 : write(iulog,*)"dED: state from end of previous dynamics (= pBF + time sampling)"
183 0 : write(iulog,*)" loop over vertical remapping and physics dribbling -------- (nsplit) -------"
184 0 : write(iulog,*)" (dribbling and remapping always done together) |"
185 0 : write(iulog,*)" dAF: state from previous remapping |"
186 0 : write(iulog,*)" dBD: state after physics dribble, before dynamics |"
187 0 : write(iulog,*)" loop over vertical Lagrangian dynamics --------rsplit------------- |"
188 0 : write(iulog,*)" dynamics here | |"
189 0 : write(iulog,*)" loop over hyperviscosity ----------hypervis_sub------------ | |"
190 0 : write(iulog,*)" dBH state before hyperviscosity | | |"
191 0 : write(iulog,*)" dCH state after hyperviscosity | | |"
192 0 : write(iulog,*)" dAH state after hyperviscosity momentum heating | | |"
193 0 : write(iulog,*)" end hyperviscosity loop ----------------------------------- | |"
194 0 : write(iulog,*)" dBS state before del2 sponge | | |"
195 0 : write(iulog,*)" dAS state after del2+mom heating sponge | | |"
196 0 : write(iulog,*)" end of vertical Lagrangian dynamics loop ------------------------- |"
197 0 : write(iulog,*)" dAD state after dynamics, before vertical remapping |"
198 0 : write(iulog,*)" dAR state after vertical remapping |"
199 0 : write(iulog,*)" end of remapping loop ------------------------------------------------------"
200 0 : write(iulog,*)"dBF state passed to parameterizations = state after last remapping "
201 0 : write(iulog,*)" "
202 0 : write(iulog,*)" "
203 0 : write(iulog,*)"FYI: all difference (diff) below are absolute normalized differences"
204 0 : write(iulog,*)" "
205 0 : write(iulog,*)"Consistency check 0:"
206 0 : write(iulog,*)"--------------------"
207 0 : write(iulog,*)" "
208 0 : write(iulog,*)"For energetic consistency we require that dE/dt [W/m^2] from energy "
209 0 : write(iulog,*)"fixer and all parameterizations computed using physics E and"
210 0 : write(iulog,*)"dycore in physics E are the same! Checking:"
211 0 : write(iulog,*)" "
212 0 : write(iulog,*) " xx=ph xx=dy norm. diff."
213 0 : write(iulog,*) " ----- ----- -----------"
214 0 : do i=1,4
215 0 : diff = abs_diff(dEdt_efix_physE(i),dEdt_efix_dynE(i),pf=pf)
216 0 : write(iulog,fmt)"dE/dt energy fixer (xxBP-xxBF) ",str(i)," ",dEdt_efix_physE(i), " ", &
217 0 : dEdt_efix_dynE(i)," ",diff,pf
218 0 : diff = abs_diff(dEdt_param_physE(i),dEdt_param_dynE(i),pf=pf)
219 0 : write(iulog,fmt)"dE/dt all parameterizations (xxAP-xxBP) ",str(i)," ",dEdt_param_physE(i)," ", &
220 0 : dEdt_param_dynE(i)," ",diff,pf
221 0 : write(iulog,*) " "
222 0 : if (diff>eps) then
223 0 : write(iulog,*)"FAIL"
224 0 : call endrun(subname//"dE/dt's in physics inconsistent")
225 : end if
226 : end do
227 0 : write(iulog,*)" "
228 0 : write(iulog,*)" "
229 0 : write(iulog,*)"dE/dt from dry-mass adjustment will differ if dynamics and physics use"
230 0 : write(iulog,*)"different energy definitions! Checking:"
231 0 : write(iulog,*)" "
232 0 : write(iulog,*) " xx=ph xx=dy diff"
233 0 : write(iulog,*) " ----- ----- ----"
234 0 : do i=1,4
235 0 : diff = dEdt_dme_adjust_physE(i)-dEdt_dme_adjust_dynE(i)
236 0 : write(iulog,fmt)"dE/dt dry mass adjustment (xxAM-xxAP) ",str(i)," ",dEdt_dme_adjust_physE(i)," ", &
237 0 : dEdt_dme_adjust_dynE(i)," ",diff
238 : end do
239 0 : write(iulog,*)" "
240 0 : write(iulog,*)" "
241 : !
242 : ! these diagnostics only make sense time-step to time-step
243 : !
244 0 : write(iulog,*)" "
245 0 : write(iulog,*)"Some energy budget observations:"
246 0 : write(iulog,*)"--------------------------------"
247 0 : write(iulog,*)" "
248 0 : write(iulog,*)"Note that total energy fixer fixes:"
249 0 : write(iulog,*) " "
250 0 : write(iulog,*) "-dE/dt energy fixer(t=n) = dE/dt dry mass adjustment (t=n-1) +"
251 0 : write(iulog,*) " dE/dt adiabatic dycore (t=n-1) +"
252 0 : write(iulog,*) " dE/dt physics-dynamics coupling errors (t=n-1)"
253 0 : write(iulog,*) " "
254 0 : write(iulog,*) "(equation 23 in Lauritzen and Williamson (2019))"
255 0 : write(iulog,*) " "
256 :
257 0 : tmp = previous_dEdt_phys_dyn_coupl_err+previous_dEdt_adiabatic_dycore+previous_dEdt_dry_mass_adjust
258 0 : diff = abs_diff(-dEdt_efix_dynE(1),tmp,pf)
259 0 : if (.not.use_cslam) then
260 0 : write(iulog,*) "Check if that is the case:", pf, diff
261 0 : write(iulog,*) " "
262 0 : if (abs(diff)>eps) then
263 0 : write(iulog,*) "dE/dt energy fixer(t=n) = ",dEdt_efix_dynE(1)
264 0 : write(iulog,*) "dE/dt dry mass adjustment (t=n-1) = ",previous_dEdt_dry_mass_adjust
265 0 : write(iulog,*) "dE/dt adiabatic dycore (t=n-1) = ",previous_dEdt_adiabatic_dycore
266 0 : write(iulog,*) "dE/dt physics-dynamics coupling errors (t=n-1) = ",previous_dEdt_phys_dyn_coupl_err
267 : end if
268 : else
269 0 : previous_dEdt_phys_dyn_coupl_err = dEdt_efix_dynE(1)+previous_dEdt_dry_mass_adjust+previous_dEdt_adiabatic_dycore
270 0 : write(iulog,*) "dE/dt energy fixer(t=n) = ",dEdt_efix_dynE(1)
271 0 : write(iulog,*) "dE/dt dry mass adjustment (t=n-1) = ",previous_dEdt_dry_mass_adjust
272 0 : write(iulog,*) "dE/dt adiabatic dycore (t=n-1) = ",previous_dEdt_adiabatic_dycore
273 0 : write(iulog,*) "dE/dt physics-dynamics coupling errors (t=n-1) = ",previous_dEdt_phys_dyn_coupl_err
274 0 : write(iulog,*) " "
275 0 : write(iulog,*) "Note: when running CSLAM the physics-dynamics coupling error is diagnosed"
276 0 : write(iulog,*) " (using equation above) rather than explicitly computed"
277 0 : write(iulog,*) " "
278 0 : write(iulog,*) " "
279 0 : write(iulog,*) "Physics-dynamics coupling errors include: "
280 0 : write(iulog,*) " "
281 0 : write(iulog,*) " -dE/dt adiabatic dycore is computed on GLL grid;"
282 0 : write(iulog,*) " error in mapping to physics grid"
283 0 : write(iulog,*) " -dE/dt physics tendencies mapped to GLL grid"
284 0 : write(iulog,*) " (tracer tendencies mapped non-conservatively!)"
285 0 : write(iulog,*) " -dE/dt dynamics state mapped to GLL grid"
286 : end if
287 0 : write(iulog,*) ""
288 0 : if (.not.use_cslam) then
289 0 : dEdt_dycore_phys = -dEdt_efix_dynE(1)-previous_dEdt_phys_dyn_coupl_err-previous_dEdt_dry_mass_adjust
290 0 : write(iulog,*) "Hence the dycore E dissipation estimated from energy fixer "
291 0 : write(iulog,'(A39,F6.2,A6)') "based on previous time-step values is ",dEdt_dycore_phys," W/M^2"
292 0 : write(iulog,*) " "
293 : end if
294 0 : write(iulog,*) " "
295 0 : write(iulog,*) "-------------------------------------------------------------------"
296 0 : write(iulog,*) " Consistency check 1: state passed to physics same as end dynamics?"
297 0 : write(iulog,*) "-------------------------------------------------------------------"
298 0 : write(iulog,*) " "
299 0 : write(iulog,*) "Is globally integrated total energy of state at the end of dynamics (dBF)"
300 0 : write(iulog,*) "and beginning of physics (using dynamics in physics energy; dyBF) the same?"
301 0 : write(iulog,*) ""
302 0 : if (.not.use_cslam) then
303 0 : if (abs(E_dyBF(1))>eps) then
304 0 : diff = abs_diff(E_dBF(1),E_dyBF(1))
305 0 : if (abs(diff)<eps) then
306 0 : write(iulog,'(A23,E8.3)')"yes. (dBF-dyBF)/dyBF =",diff
307 : else
308 0 : write(iulog,*)"Error in physics dynamics coupling!"
309 0 : write(iulog,*)" "
310 0 : do i=1,4
311 0 : write(iulog,*) str(i),":"
312 0 : write(iulog,*) "======"
313 0 : diff = abs_diff(E_dBF(i),E_dyBF(i),pf=pf)
314 0 : write(iulog,*) "diff, E_dBF, E_dyBF ",diff,E_dBF(i),E_dyBF(i)
315 0 : write(iulog,*) " "
316 : end do
317 0 : call endrun(subname//"Error in physics dynamics coupling")
318 : end if
319 : end if
320 : else
321 0 : write(iulog,*)" "
322 0 : write(iulog,*)"Since you are using a separate physics grid, the state in dynamics"
323 0 : write(iulog,*)"will not be the same on the physics grid since it is"
324 0 : write(iulog,*)"interpolated from the dynamics to the physics grid"
325 0 : write(iulog,*)" "
326 0 : do i=1,4
327 0 : write(iulog,*) str(i),":"
328 0 : write(iulog,*) "======"
329 0 : diff = abs_diff(E_dBF(i),E_dyBF(i),pf=pf)
330 0 : write(iulog,*) "diff, E_dBF, E_dyBF ",diff,E_dBF(i),E_dyBF(i)
331 0 : write(iulog,*) " "
332 : end do
333 : end if
334 :
335 0 : write(iulog,*)" "
336 0 : write(iulog,*)"-------------------------------------------------------------------------"
337 0 : write(iulog,*)" Consistency check 2: total energy increment in dynamics same as physics?"
338 0 : write(iulog,*)"-------------------------------------------------------------------------"
339 0 : write(iulog,*)" "
340 0 : if (.not.use_cslam) then
341 0 : previous_dEdt_phys_dyn_coupl_err = dEdt_phys_total_in_dyn(1)-dEdt_phys_total_dynE(1)
342 0 : diff = abs_diff(dEdt_phys_total_dynE(1),dEdt_phys_total_in_dyn(1),pf=pf)
343 0 : write(iulog,'(A40,E8.2,A7,A5)')" dE/dt physics-dynamics coupling errors ",diff," W/M^2 ",pf
344 0 : if (abs(diff)>eps) then
345 : !
346 : ! if errors print details
347 : !
348 0 : if (ftype==1) then
349 0 : write(iulog,*) ""
350 0 : write(iulog,*) " You are using ftype==1 so physics-dynamics coupling errors should be round-off!"
351 0 : write(iulog,*) ""
352 0 : write(iulog,*) " Because of failure provide detailed diagnostics below:"
353 0 : write(iulog,*) ""
354 : else
355 0 : write(iulog,*) ""
356 0 : write(iulog,*) " Since ftype<>1 there are physics dynamics coupling errors"
357 0 : write(iulog,*) ""
358 0 : write(iulog,*) " Break-down below:"
359 0 : write(iulog,*) ""
360 : end if
361 :
362 0 : do i=1,4
363 0 : write(iulog,*) str(i),":"
364 0 : write(iulog,*) "======"
365 0 : diff = abs_diff(dEdt_phys_total_dynE(i),dEdt_phys_total_in_dyn(i),pf=pf)
366 0 : write(iulog,*) "dE/dt physics-dynamics coupling errors (diff) ",diff
367 0 : write(iulog,*) "dE/dt physics total in dynamics (dBD-dAF) ",dEdt_phys_total_in_dyn(i)
368 0 : write(iulog,*) "dE/dt physics total in physics (dyAM-dyBF) ",dEdt_phys_total_dynE(i)
369 0 : write(iulog,*) " "
370 0 : write(iulog,*) " physics total = parameterizations + efix + dry-mass adjustment"
371 0 : write(iulog,*) " "
372 : end do
373 : ! Temporarily disable endrun until energy bias for consistancy check 2 is better understood.
374 : ! if (ftype==1) then
375 : ! call endrun(subname//"Physics-dynamics coupling error. See atm.log")
376 : ! end if
377 : end if
378 : else
379 0 : write(iulog,'(a47,F6.2,a6)')" dE/dt physics tendency in dynamics (dBD-dAF) ",dEdt_phys_total_in_dyn(1)," W/M^2"
380 0 : write(iulog,'(a47,F6.2,a6)')" dE/dt physics tendency in physics (dyAM-dyBF) ",dEdt_phys_total_dynE(1)," W/M^2"
381 0 : write(iulog,*)" "
382 0 : write(iulog,*) " When runnig with a physics grid this consistency check does not make sense"
383 0 : write(iulog,*) " since it is computed on the GLL grid whereas we enforce energy conservation"
384 0 : write(iulog,*) " on the physics grid. To assess the errors of running dynamics on GLL"
385 0 : write(iulog,*) " grid, tracers on CSLAM grid and physics on physics grid we use the energy"
386 0 : write(iulog,*) " fixer check from above:"
387 0 : write(iulog,*) " "
388 0 : write(iulog,*) " dE/dt physics-dynamics coupling errors (t=n-1) =",previous_dEdt_phys_dyn_coupl_err
389 0 : write(iulog,*) ""
390 : end if
391 0 : write(iulog,*)" "
392 0 : write(iulog,*)"------------------------------------------------------------"
393 0 : write(iulog,*)" SE dycore energy tendencies"
394 0 : write(iulog,*)"------------------------------------------------------------"
395 0 : write(iulog,*)" "
396 0 : write(iulog,fmtf)" dE/dt dycore ",dEdt_dycore_dyn," W/M^2"
397 0 : write(iulog,*)" "
398 0 : write(iulog,*)"Adiabatic dynamics can be divided into quasi-horizontal and vertical remapping: "
399 0 : write(iulog,*)" "
400 0 : write(iulog,fmtf)" dE/dt floating dynamics (dAD-dBD) ",dEdt_floating_dyn," W/M^2"
401 0 : write(iulog,fmtf)" dE/dt vertical remapping (dAR-dAD) ",dEdt_vert_remap," W/M^2"
402 :
403 0 : write(iulog,*) " "
404 0 : write(iulog,*) "Breakdown of floating dynamics:"
405 0 : write(iulog,*) " "
406 0 : write(iulog,fmtf)" dE/dt hypervis del4 (dCH-dBH) ",dEdt_del4, " W/M^2"
407 0 : write(iulog,fmtf)" dE/dt hypervis frictional heating (dAH-dCH) ",dEdt_del4_fric_heat," W/M^2"
408 0 : write(iulog,fmtf)" dE/dt hypervis del4 total (dAH-dBH) ",dEdt_del4_tot, " W/M^2"
409 0 : write(iulog,fmtf)" dE/dt hypervis sponge del2 (dAS-dBS) ",dEdt_del2_sponge, " W/M^2"
410 0 : write(iulog,fmtf)" dE/dt explicit diffusion total ",dEdt_del2_del4_tot, " W/M^2"
411 0 : write(iulog,*) " "
412 0 : write(iulog,fmtf)" dE/dt residual (time-truncation errors,...) ",dEdt_residual, " W/M^2"
413 0 : write(iulog,*)" "
414 0 : write(iulog,*)" "
415 0 : write(iulog,*)"------------------------------------------------------------"
416 0 : write(iulog,*)"Tracer mass budgets"
417 0 : write(iulog,*)"------------------------------------------------------------"
418 0 : write(iulog,*)" "
419 0 : write(iulog,*)"Below the physics-dynamics coupling error is computed as "
420 0 : write(iulog,*)"dMASS/dt physics tendency in dycore (dBD-dAF) minus"
421 0 : write(iulog,*)"dMASS/dt total physics (pAM-pBF)"
422 0 : write(iulog,*)" "
423 0 : write(iulog,*)" "
424 0 : do m_cnst=1,thermo_budget_num_vars
425 0 : if (thermo_budget_vars_massv(m_cnst)) then
426 0 : write(iulog,*)thermo_budget_vars_descriptor(m_cnst)
427 0 : write(iulog,*)"------------------------------"
428 0 : call cam_budget_get_global('phBP-phBF',m_cnst,dMdt_efix)
429 0 : call cam_budget_get_global('phAM-phAP',m_cnst,dMdt_dme_adjust)
430 0 : call cam_budget_get_global('phAP-phBP',m_cnst,dMdt_parameterizations)
431 0 : call cam_budget_get_global('phAM-phBF',m_cnst,dMdt_phys_total)
432 : !
433 : ! total energy fixer should not affect mass - checking
434 : !
435 0 : if (abs(dMdt_efix)>eps_mass) then
436 0 : write(iulog,*) "dMASS/dt energy fixer (pBP-pBF) ",dMdt_efix," Pa/m^2/s"
437 0 : write(iulog,*) "ERROR: Mass not conserved in energy fixer. ABORT"
438 0 : call endrun(subname//"Mass not conserved in energy fixer. See atm.log")
439 : endif
440 : !
441 : ! dry-mass adjustmnt should not affect mass - checking
442 : !
443 0 : if (abs(dMdt_dme_adjust)>eps_mass) then
444 0 : write(iulog,*)"dMASS/dt dry mass adjustment (pAM-pAP) ",dMdt_dme_adjust," Pa/m^2/s"
445 0 : write(iulog,*) "ERROR: Mass not conserved in dry mass adjustment. ABORT"
446 0 : call endrun(subname//"Mass not conserved in dry mass adjustment. See atm.log")
447 : end if
448 : !
449 : ! all of the mass-tendency should come from parameterization - checking
450 : !
451 0 : if (abs(dMdt_parameterizations-dMdt_phys_total)>eps_mass) then
452 0 : write(iulog,*) "Error: dMASS/dt parameterizations (pAP-pBP) .ne. dMASS/dt physics total (pAM-pBF)"
453 0 : write(iulog,*) "dMASS/dt parameterizations (pAP-pBP) ",dMdt_parameterizations," Pa/m^2/s"
454 0 : write(iulog,*) "dMASS/dt physics total (pAM-pBF) ",dMdt_phys_total," Pa/m^2/s"
455 0 : call endrun(subname//"mass change not only due to parameterizations. See atm.log")
456 : end if
457 0 : write(iulog,*)" "
458 : !
459 : ! detailed mass budget in dynamical core
460 : !
461 0 : if (is_cam_budget('dAD').and.is_cam_budget('dBD').and.is_cam_budget('dAR').and.is_cam_budget('dCH')) then
462 0 : call cam_budget_get_global('dAL-dBL',m_cnst,dMdt_floating_dyn)
463 0 : call cam_budget_get_global('dAR-dAD',m_cnst,dMdt_vert_remap)
464 0 : tmp = dMdt_floating_dyn+dMdt_vert_remap
465 0 : diff = abs_diff(tmp,0.0_r8,pf=pf)
466 0 : write(iulog,fmtm)" dMASS/dt total adiabatic dynamics ",diff,pf
467 : !
468 : ! check for mass-conservation in the adiabatic dynamical core -
469 : ! if not conserved provide detailed break-down
470 : !
471 0 : if (abs(diff)>eps_mass) then
472 0 : write(iulog,*) "Error: mass non-conservation in dynamical core"
473 0 : write(iulog,*) "(detailed budget below)"
474 0 : write(iulog,*) " "
475 0 : write(iulog,*)"dMASS/dt 2D dynamics (dAL-dBL) ",dMdt_floating_dyn," Pa/m^2/s"
476 0 : write(iulog,*)"dE/dt vertical remapping (dAR-dAD) ",dMdt_vert_remap
477 0 : write(iulog,*)" "
478 0 : write(iulog,*)"Breakdown of 2D dynamics:"
479 0 : write(iulog,*)" "
480 0 : call cam_budget_get_global('dAH-dCH',m_cnst,dMdt_del4_fric_heat)
481 0 : call cam_budget_get_global('dAH-dBH',m_cnst,dMdt_del4_tot)
482 0 : write(iulog,*)"dMASS/dt hypervis (dAH-dBH) ",dMdt_del4_tot," Pa/m^2/s"
483 0 : write(iulog,*)"dMASS/dt frictional heating (dAH-dCH) ",dMdt_del4_fric_heat," Pa/m^2/s"
484 0 : dMdt_residual = dMdt_floating_dyn-dMdt_del4_tot
485 0 : write(iulog,*)"dMASS/dt residual (time truncation errors)",dMdt_residual," Pa/m^2/s"
486 : end if
487 : end if
488 0 : if (is_cam_budget('dBD').and.is_cam_budget('dAF')) then
489 : !
490 : ! check if mass change in physics is the same as dynamical core
491 : !
492 0 : call cam_budget_get_global('dBD-dAF',m_cnst,dMdt_phys_total_in_dyn)
493 0 : dMdt_PDC = dMdt_phys_total-dMdt_phys_total_in_dyn
494 0 : write(iulog,fmtm)" Mass physics-dynamics coupling error ",dMdt_PDC," Pa/m^2/s"
495 0 : write(iulog,*)" "
496 0 : if (abs(dMdt_PDC)>eps_mass) then
497 0 : write(iulog,fmtm)" dMASS/dt physics tendency in dycore (dBD-dAF) ",dMdt_phys_total_in_dyn," Pa/m^2/s"
498 0 : write(iulog,fmtm)" dMASS/dt total physics ",dMdt_phys_total," Pa/m^2/s"
499 : end if
500 : end if
501 : end if
502 : end do
503 : !
504 : ! save adiabatic dycore dE/dt and dry-mass adjustment to avoid samping error
505 : !
506 0 : previous_dEdt_adiabatic_dycore = dEdt_dycore_dyn
507 0 : previous_dEdt_dry_mass_adjust = dEdt_dme_adjust_dynE(1)
508 : end if
509 369408 : end subroutine print_budget
510 : !=========================================================================================
511 0 : function abs_diff(a,b,pf)
512 : real(r8), intent(in) :: a,b
513 : character(LEN=5), optional, intent(out):: pf
514 : real(r8) :: abs_diff
515 0 : if (abs(b)>eps) then
516 0 : abs_diff = abs((b-a)/b)
517 : else
518 0 : abs_diff = abs(b-a)
519 : end if
520 0 : If (present(pf)) then
521 0 : if (abs_diff>eps) then
522 0 : pf = ' FAIL'
523 : else
524 0 : pf = ' PASS'
525 : end if
526 : end if
527 369408 : end function abs_diff
528 : end module dycore_budget
|