Line data Source code
1 : !-------------------------------------------------------------------------------
2 : !physics data types module
3 : !-------------------------------------------------------------------------------
4 : module physics_types
5 :
6 : use shr_kind_mod, only: r8 => shr_kind_r8
7 : use ppgrid, only: pcols, pver
8 : use constituents, only: pcnst, qmin, cnst_name, cnst_get_ind
9 : use geopotential, only: geopotential_t
10 : use physconst, only: zvir, gravit, cpair, rair
11 : use air_composition, only: cpairv, rairv
12 : use phys_grid, only: get_ncols_p, get_rlon_all_p, get_rlat_all_p, get_gcol_all_p
13 : use cam_logfile, only: iulog
14 : use cam_abortutils, only: endrun
15 : use phys_control, only: waccmx_is
16 : use shr_const_mod, only: shr_const_rwv
17 :
18 : implicit none
19 : private ! Make default type private to the module
20 :
21 : ! Public types:
22 :
23 : public physics_state
24 : public physics_tend
25 : public physics_ptend
26 :
27 : ! Public interfaces
28 :
29 : public physics_update
30 : public physics_state_check ! Check state object for invalid data.
31 : public physics_ptend_reset
32 : public physics_ptend_init
33 : public physics_state_set_grid
34 : public physics_dme_adjust ! adjust dry mass and energy for change in water
35 : public physics_state_copy ! copy a physics_state object
36 : public physics_ptend_copy ! copy a physics_ptend object
37 : public physics_ptend_sum ! accumulate physics_ptend objects
38 : public physics_ptend_scale ! Multiply physics_ptend objects by a constant factor.
39 : public physics_tend_init ! initialize a physics_tend object
40 :
41 : public set_state_pdry ! calculate dry air masses in state variable
42 : public set_wet_to_dry
43 : public set_dry_to_wet
44 : public physics_type_alloc
45 :
46 : public physics_state_alloc ! allocate individual components within state
47 : public physics_state_dealloc ! deallocate individual components within state
48 : public physics_tend_alloc ! allocate individual components within tend
49 : public physics_tend_dealloc ! deallocate individual components within tend
50 : public physics_ptend_alloc ! allocate individual components within tend
51 : public physics_ptend_dealloc ! deallocate individual components within tend
52 :
53 : public physics_cnst_limit ! apply limiters to constituents (waccmx)
54 : !-------------------------------------------------------------------------------
55 : integer, parameter, public :: phys_te_idx = 1
56 : integer ,parameter, public :: dyn_te_idx = 2
57 :
58 : type physics_state
59 : integer :: &
60 : lchnk, &! chunk index
61 : ngrdcol, &! -- Grid -- number of active columns (on the grid)
62 : psetcols=0, &! -- -- max number of columns set - if subcols = pcols*psubcols, else = pcols
63 : ncol=0 ! -- -- sum of nsubcol for all ngrdcols - number of active columns
64 : real(r8), dimension(:), allocatable :: &
65 : lat, &! latitude (radians)
66 : lon, &! longitude (radians)
67 : ps, &! surface pressure
68 : psdry, &! dry surface pressure
69 : phis, &! surface geopotential
70 : ulat, &! unique latitudes (radians)
71 : ulon ! unique longitudes (radians)
72 : real(r8), dimension(:,:),allocatable :: &
73 : t, &! temperature (K)
74 : u, &! zonal wind (m/s)
75 : v, &! meridional wind (m/s)
76 : s, &! dry static energy
77 : omega, &! vertical pressure velocity (Pa/s)
78 : pmid, &! midpoint pressure (Pa)
79 : pmiddry, &! midpoint pressure dry (Pa)
80 : pdel, &! layer thickness (Pa)
81 : pdeldry, &! layer thickness dry (Pa)
82 : rpdel, &! reciprocal of layer thickness (Pa)
83 : rpdeldry,&! recipricol layer thickness dry (Pa)
84 : lnpmid, &! ln(pmid)
85 : lnpmiddry,&! log midpoint pressure dry (Pa)
86 : exner, &! inverse exner function w.r.t. surface pressure (ps/p)^(R/cp)
87 : zm ! geopotential height above surface at midpoints (m)
88 :
89 : real(r8), dimension(:,:,:),allocatable :: &
90 : q ! constituent mixing ratio (kg/kg moist or dry air depending on type)
91 :
92 : real(r8), dimension(:,:),allocatable :: &
93 : pint, &! interface pressure (Pa)
94 : pintdry, &! interface pressure dry (Pa)
95 : lnpint, &! ln(pint)
96 : lnpintdry,&! log interface pressure dry (Pa)
97 : zi ! geopotential height above surface at interfaces (m)
98 :
99 : real(r8), dimension(:,:),allocatable :: &
100 : ! Second dimension is (phys_te_idx) CAM physics total energy and
101 : ! (dyn_te_idx) dycore total energy computed in physics
102 : te_ini, &! vertically integrated total (kinetic + static) energy of initial state
103 : te_cur ! vertically integrated total (kinetic + static) energy of current state
104 : real(r8), dimension(:), allocatable :: &
105 : tw_ini, &! vertically integrated total water of initial state
106 : tw_cur ! vertically integrated total water of new state
107 : real(r8), dimension(:,:),allocatable :: &
108 : temp_ini, &! Temperature of initial state (used for energy computations)
109 : z_ini ! Height of initial state (used for energy computations)
110 : integer :: count ! count of values with significant energy or water imbalances
111 : integer, dimension(:),allocatable :: &
112 : latmapback, &! map from column to unique lat for that column
113 : lonmapback, &! map from column to unique lon for that column
114 : cid ! unique column id
115 : integer :: ulatcnt, &! number of unique lats in chunk
116 : uloncnt ! number of unique lons in chunk
117 :
118 : end type physics_state
119 :
120 : !-------------------------------------------------------------------------------
121 : type physics_tend
122 :
123 : integer :: psetcols=0 ! max number of columns set- if subcols = pcols*psubcols, else = pcols
124 :
125 : real(r8), dimension(:,:),allocatable :: dtdt, dudt, dvdt
126 : real(r8), dimension(:), allocatable :: flx_net
127 : real(r8), dimension(:), allocatable :: &
128 : te_tnd, &! cumulative boundary flux of total energy
129 : tw_tnd ! cumulative boundary flux of total water
130 : end type physics_tend
131 :
132 : !-------------------------------------------------------------------------------
133 : ! This is for tendencies returned from individual parameterizations
134 : type physics_ptend
135 :
136 : integer :: psetcols=0 ! max number of columns set- if subcols = pcols*psubcols, else = pcols
137 :
138 : character*24 :: name ! name of parameterization which produced tendencies.
139 :
140 : logical :: &
141 : ls = .false., &! true if dsdt is returned
142 : lu = .false., &! true if dudt is returned
143 : lv = .false. ! true if dvdt is returned
144 :
145 : logical,dimension(pcnst) :: lq = .false. ! true if dqdt() is returned
146 :
147 : integer :: &
148 : top_level, &! top level index for which nonzero tendencies have been set
149 : bot_level ! bottom level index for which nonzero tendencies have been set
150 :
151 : real(r8), dimension(:,:),allocatable :: &
152 : s, &! heating rate (J/kg/s)
153 : u, &! u momentum tendency (m/s/s)
154 : v ! v momentum tendency (m/s/s)
155 : real(r8), dimension(:,:,:),allocatable :: &
156 : q ! consituent tendencies (kg/kg/s)
157 :
158 : ! boundary fluxes
159 : real(r8), dimension(:),allocatable ::&
160 : hflux_srf, &! net heat flux at surface (W/m2)
161 : hflux_top, &! net heat flux at top of model (W/m2)
162 : taux_srf, &! net zonal stress at surface (Pa)
163 : taux_top, &! net zonal stress at top of model (Pa)
164 : tauy_srf, &! net meridional stress at surface (Pa)
165 : tauy_top ! net meridional stress at top of model (Pa)
166 : real(r8), dimension(:,:),allocatable ::&
167 : cflx_srf, &! constituent flux at surface (kg/m2/s)
168 : cflx_top ! constituent flux top of model (kg/m2/s)
169 :
170 : end type physics_ptend
171 :
172 :
173 : !===============================================================================
174 : contains
175 : !===============================================================================
176 1536 : subroutine physics_type_alloc(phys_state, phys_tend, begchunk, endchunk, psetcols)
177 : implicit none
178 : type(physics_state), pointer :: phys_state(:)
179 : type(physics_tend), pointer :: phys_tend(:)
180 : integer, intent(in) :: begchunk, endchunk
181 : integer, intent(in) :: psetcols
182 :
183 : integer :: ierr=0, lchnk
184 :
185 12288 : allocate(phys_state(begchunk:endchunk), stat=ierr)
186 1536 : if( ierr /= 0 ) then
187 0 : write(iulog,*) 'physics_types: phys_state allocation error = ',ierr
188 0 : call endrun('physics_types: failed to allocate physics_state array')
189 : end if
190 :
191 9216 : do lchnk=begchunk,endchunk
192 9216 : call physics_state_alloc(phys_state(lchnk),lchnk,pcols)
193 : end do
194 :
195 12288 : allocate(phys_tend(begchunk:endchunk), stat=ierr)
196 1536 : if( ierr /= 0 ) then
197 0 : write(iulog,*) 'physics_types: phys_tend allocation error = ',ierr
198 0 : call endrun('physics_types: failed to allocate physics_tend array')
199 : end if
200 :
201 9216 : do lchnk=begchunk,endchunk
202 9216 : call physics_tend_alloc(phys_tend(lchnk),phys_state(lchnk)%psetcols)
203 : end do
204 :
205 1536 : end subroutine physics_type_alloc
206 : !===============================================================================
207 3624960 : subroutine physics_update(state, ptend, dt, tend)
208 : !-----------------------------------------------------------------------
209 : ! Update the state and or tendency structure with the parameterization tendencies
210 : !-----------------------------------------------------------------------
211 : use scamMod, only: scm_crm_mode, single_column
212 : use phys_control, only: phys_getopts
213 : use cam_thermo, only: cam_thermo_dry_air_update ! Routine which updates physconst variables (WACCM-X)
214 : use air_composition, only: dry_air_species_num
215 : use qneg_module , only: qneg3
216 :
217 : !------------------------------Arguments--------------------------------
218 : type(physics_ptend), intent(inout) :: ptend ! Parameterization tendencies
219 :
220 : type(physics_state), intent(inout) :: state ! Physics state variables
221 :
222 : real(r8), intent(in) :: dt ! time step
223 :
224 : type(physics_tend ), intent(inout), optional :: tend ! Physics tendencies over timestep
225 : ! tend is usually only needed by calls from physpkg.
226 : !
227 : !---------------------------Local storage-------------------------------
228 : integer :: k,m ! column,level,constituent indices
229 : integer :: ixcldice, ixcldliq ! indices for CLDICE and CLDLIQ
230 : integer :: ixnumice, ixnumliq
231 : integer :: ixnumsnow, ixnumrain
232 : integer :: ncol ! number of columns
233 : integer :: ixh, ixh2 ! constituent indices for H, H2
234 :
235 7249920 : real(r8) :: zvirv(state%psetcols,pver) ! Local zvir array pointer
236 :
237 3624960 : real(r8),allocatable :: cpairv_loc(:,:)
238 3624960 : real(r8),allocatable :: rairv_loc(:,:)
239 :
240 : ! PERGRO limits cldliq/ice for macro/microphysics:
241 : character(len=24), parameter :: pergro_cldlim_names(4) = &
242 : (/ "stratiform", "cldwat ", "micro_mg ", "macro_park" /)
243 :
244 : ! cldliq/ice limits that are always on.
245 : character(len=24), parameter :: cldlim_names(2) = &
246 : (/ "convect_deep", "zm_conv_tend" /)
247 :
248 : ! Whether to do validation of state on each call.
249 : logical :: state_debug_checks
250 :
251 : !-----------------------------------------------------------------------
252 :
253 : ! The column radiation model does not update the state
254 3624960 : if(single_column.and.scm_crm_mode) return
255 :
256 :
257 : !-----------------------------------------------------------------------
258 : ! If no fields are set, then return
259 197280000 : if (.not. (any(ptend%lq(:)) .or. ptend%ls .or. ptend%lu .or. ptend%lv)) then
260 299520 : ptend%name = "none"
261 299520 : ptend%psetcols = 0
262 299520 : return
263 : end if
264 :
265 : !-----------------------------------------------------------------------
266 : ! Check that the state/tend/ptend are all dimensioned with the same number of columns
267 3325440 : if (state%psetcols /= ptend%psetcols) then
268 : call endrun('ERROR in physics_update with ptend%name='//trim(ptend%name) &
269 0 : //': state and ptend must have the same number of psetcols.')
270 : end if
271 :
272 3325440 : if (present(tend)) then
273 1793280 : if (state%psetcols /= tend%psetcols) then
274 : call endrun('ERROR in physics_update with ptend%name='//trim(ptend%name) &
275 0 : //': state and tend must have the same number of psetcols.')
276 : end if
277 : end if
278 :
279 :
280 : !-----------------------------------------------------------------------
281 3325440 : call phys_getopts(state_debug_checks_out=state_debug_checks)
282 :
283 3325440 : ncol = state%ncol
284 :
285 : ! Update u,v fields
286 3325440 : if(ptend%lu) then
287 71704320 : do k = ptend%top_level, ptend%bot_level
288 1088693760 : state%u (:ncol,k) = state%u (:ncol,k) + ptend%u(:ncol,k) * dt
289 70694400 : if (present(tend)) &
290 741984000 : tend%dudt(:ncol,k) = tend%dudt(:ncol,k) + ptend%u(:ncol,k)
291 : end do
292 : end if
293 :
294 3325440 : if(ptend%lv) then
295 66524160 : do k = ptend%top_level, ptend%bot_level
296 1010042880 : state%v (:ncol,k) = state%v (:ncol,k) + ptend%v(:ncol,k) * dt
297 65587200 : if (present(tend)) &
298 663260160 : tend%dvdt(:ncol,k) = tend%dvdt(:ncol,k) + ptend%v(:ncol,k)
299 : end do
300 : end if
301 :
302 : ! Update constituents, all schemes use time split q: no tendency kept
303 3325440 : call cnst_get_ind('CLDICE', ixcldice, abort=.false.)
304 3325440 : call cnst_get_ind('CLDLIQ', ixcldliq, abort=.false.)
305 : ! Check for number concentration of cloud liquid and cloud ice (if not present
306 : ! the indices will be set to -1)
307 3325440 : call cnst_get_ind('NUMICE', ixnumice, abort=.false.)
308 3325440 : call cnst_get_ind('NUMLIQ', ixnumliq, abort=.false.)
309 3325440 : call cnst_get_ind('NUMRAI', ixnumrain, abort=.false.)
310 3325440 : call cnst_get_ind('NUMSNO', ixnumsnow, abort=.false.)
311 :
312 675064320 : do m = 1, pcnst
313 675064320 : if(ptend%lq(m)) then
314 14252256000 : do k = ptend%top_level, ptend%bot_level
315 >21659*10^7 : state%q(:ncol,k,m) = state%q(:ncol,k,m) + ptend%q(:ncol,k,m) * dt
316 : end do
317 :
318 : ! now test for mixing ratios which are too small
319 : ! don't call qneg3 for number concentration variables
320 : if (m /= ixnumice .and. m /= ixnumliq .and. &
321 200736000 : m /= ixnumrain .and. m /= ixnumsnow ) then
322 194991360 : call qneg3(trim(ptend%name), state%lchnk, ncol, state%psetcols, pver, m, m, qmin(m:m), state%q(:,1:pver,m:m))
323 : else
324 407869440 : do k = ptend%top_level, ptend%bot_level
325 : ! checks for number concentration
326 6192721920 : state%q(:ncol,k,m) = max(1.e-12_r8,state%q(:ncol,k,m))
327 6198466560 : state%q(:ncol,k,m) = min(1.e10_r8,state%q(:ncol,k,m))
328 : end do
329 : end if
330 :
331 : end if
332 :
333 : end do
334 :
335 : !------------------------------------------------------------------------
336 : ! This is a temporary fix for the large H, H2 in WACCM-X
337 : ! Well, it was supposed to be temporary, but it has been here
338 : ! for a while now.
339 : !------------------------------------------------------------------------
340 3325440 : if ( waccmx_is('ionosphere') .or. waccmx_is('neutral') ) then
341 0 : call cnst_get_ind('H', ixh)
342 0 : do k = ptend%top_level, ptend%bot_level
343 0 : state%q(:ncol,k,ixh) = min(state%q(:ncol,k,ixh), 0.01_r8)
344 : end do
345 :
346 0 : call cnst_get_ind('H2', ixh2)
347 0 : do k = ptend%top_level, ptend%bot_level
348 0 : state%q(:ncol,k,ixh2) = min(state%q(:ncol,k,ixh2), 6.e-5_r8)
349 : end do
350 : endif
351 :
352 : ! Special tests for cloud liquid and ice:
353 : ! Enforce a minimum non-zero value.
354 3325440 : if (ixcldliq > 1) then
355 3325440 : if(ptend%lq(ixcldliq)) then
356 : #ifdef PERGRO
357 : if ( any(ptend%name == pergro_cldlim_names) ) &
358 : call state_cnst_min_nz(1.e-12_r8, ixcldliq, ixnumliq)
359 : #endif
360 4711680 : if ( any(ptend%name == cldlim_names) ) &
361 80640 : call state_cnst_min_nz(1.e-36_r8, ixcldliq, ixnumliq)
362 : end if
363 : end if
364 :
365 3325440 : if (ixcldice > 1) then
366 3325440 : if(ptend%lq(ixcldice)) then
367 : #ifdef PERGRO
368 : if ( any(ptend%name == pergro_cldlim_names) ) &
369 : call state_cnst_min_nz(1.e-12_r8, ixcldice, ixnumice)
370 : #endif
371 4953600 : if ( any(ptend%name == cldlim_names) ) &
372 80640 : call state_cnst_min_nz(1.e-36_r8, ixcldice, ixnumice)
373 : end if
374 : end if
375 :
376 : !------------------------------------------------------------------------
377 : ! Get indices for molecular weights and call WACCM-X cam_thermo_update
378 : !------------------------------------------------------------------------
379 3325440 : if (dry_air_species_num>0) then
380 0 : call cam_thermo_dry_air_update(state%q, state%t, state%lchnk, state%ncol)
381 : endif
382 :
383 : !-----------------------------------------------------------------------
384 : ! cpairv_loc and rairv_loc need to be allocated to a size which matches state and ptend
385 : ! If psetcols == pcols, the cpairv is the correct size and just copy
386 : ! If psetcols > pcols and all cpairv match cpair, then assign the constant cpair
387 9976320 : allocate(cpairv_loc(state%psetcols,pver))
388 3325440 : if (state%psetcols == pcols) then
389 3960599040 : cpairv_loc(:,:) = cpairv(:,:,state%lchnk)
390 0 : else if (state%psetcols > pcols .and. all(cpairv(:,:,:) == cpair)) then
391 0 : cpairv_loc(:,:) = cpair
392 : else
393 0 : call endrun('physics_update: cpairv is not allowed to vary when subcolumns are turned on')
394 : end if
395 9976320 : allocate(rairv_loc(state%psetcols,pver))
396 3325440 : if (state%psetcols == pcols) then
397 3960599040 : rairv_loc(:,:) = rairv(:,:,state%lchnk)
398 0 : else if (state%psetcols > pcols .and. all(rairv(:,:,:) == rair)) then
399 0 : rairv_loc(:,:) = rair
400 : else
401 0 : call endrun('physics_update: rairv_loc is not allowed to vary when subcolumns are turned on')
402 : end if
403 :
404 3325440 : if ( waccmx_is('ionosphere') .or. waccmx_is('neutral') ) then
405 0 : zvirv(:,:) = shr_const_rwv / rairv_loc(:,:) - 1._r8
406 : else
407 3960599040 : zvirv(:,:) = zvir
408 : endif
409 :
410 : !-------------------------------------------------------------------------------------------------------------
411 : ! Update temperature from dry static energy (moved from above for WACCM-X so updating after cpairv_loc update)
412 : !-------------------------------------------------------------------------------------------------------------
413 :
414 3325440 : if(ptend%ls) then
415 163856640 : do k = ptend%top_level, ptend%bot_level
416 2487851520 : state%t(:ncol,k) = state%t(:ncol,k) + ptend%s(:ncol,k)*dt/cpairv_loc(:ncol,k)
417 161548800 : if (present(tend)) &
418 1360070400 : tend%dtdt(:ncol,k) = tend%dtdt(:ncol,k) + ptend%s(:ncol,k)/cpairv_loc(:ncol,k)
419 : end do
420 : end if
421 :
422 : ! Derive new geopotential fields if heating or water tendency not 0.
423 :
424 3325440 : if (ptend%ls .or. ptend%lq(1)) then
425 : call geopotential_t ( &
426 : state%lnpint, state%lnpmid, state%pint , state%pmid , state%pdel , state%rpdel , &
427 : state%t , state%q(:,:,:), rairv_loc(:,:), gravit , zvirv , &
428 2380800 : state%zi , state%zm , ncol )
429 : ! update dry static energy for use in next process
430 169036800 : do k = ptend%top_level, ptend%bot_level
431 166656000 : state%s(:ncol,k) = state%t(:ncol,k)*cpairv_loc(:ncol,k) &
432 2735539200 : + gravit*state%zm(:ncol,k) + state%phis(:ncol)
433 : end do
434 : end if
435 :
436 3325440 : if (state_debug_checks) call physics_state_check(state, ptend%name)
437 :
438 3325440 : deallocate(cpairv_loc, rairv_loc)
439 :
440 : ! Deallocate ptend
441 3325440 : call physics_ptend_dealloc(ptend)
442 :
443 3325440 : ptend%name = "none"
444 675064320 : ptend%lq(:) = .false.
445 3325440 : ptend%ls = .false.
446 3325440 : ptend%lu = .false.
447 3325440 : ptend%lv = .false.
448 6950400 : ptend%psetcols = 0
449 :
450 : contains
451 :
452 161280 : subroutine state_cnst_min_nz(lim, qix, numix)
453 : ! Small utility function for setting minimum nonzero
454 : ! constituent concentrations.
455 :
456 : ! Lower limit and constituent index
457 : real(r8), intent(in) :: lim
458 : integer, intent(in) :: qix
459 : ! Number concentration that goes with qix.
460 : ! Ignored if <= 0 (and therefore constituent is not present).
461 : integer, intent(in) :: numix
462 :
463 161280 : if (numix > 0) then
464 : ! Where q is too small, zero mass and number
465 : ! concentration.
466 522063360 : where (state%q(:ncol,:,qix) < lim)
467 161280 : state%q(:ncol,:,qix) = 0._r8
468 161280 : state%q(:ncol,:,numix) = 0._r8
469 : end where
470 : else
471 : ! If no number index, just do mass.
472 0 : where (state%q(:ncol,:,qix) < lim)
473 0 : state%q(:ncol,:,qix) = 0._r8
474 : end where
475 : end if
476 :
477 3624960 : end subroutine state_cnst_min_nz
478 :
479 :
480 : end subroutine physics_update
481 :
482 : !===============================================================================
483 :
484 3559680 : subroutine physics_state_check(state, name)
485 : !-----------------------------------------------------------------------
486 : ! Check a physics_state object for invalid data (e.g NaNs, negative
487 : ! temperatures).
488 : !-----------------------------------------------------------------------
489 : use shr_infnan_mod, only: assignment(=), &
490 : shr_infnan_posinf, shr_infnan_neginf
491 : use shr_assert_mod, only: shr_assert_in_domain
492 : use constituents, only: pcnst
493 :
494 : !------------------------------Arguments--------------------------------
495 : ! State to check.
496 : type(physics_state), intent(in) :: state
497 : ! Name of the package responsible for this state.
498 : character(len=*), intent(in), optional :: name
499 :
500 : !---------------------------Local storage-------------------------------
501 : ! Shortened name for ncol.
502 : integer :: ncol
503 : ! Double precision positive/negative infinity.
504 : real(r8) :: posinf_r8, neginf_r8
505 : ! Canned message.
506 : character(len=64) :: msg
507 : ! Constituent index
508 : integer :: m
509 :
510 : !-----------------------------------------------------------------------
511 :
512 3559680 : ncol = state%ncol
513 :
514 3559680 : posinf_r8 = shr_infnan_posinf
515 3559680 : neginf_r8 = shr_infnan_neginf
516 :
517 : ! It may be reasonable to check some of the integer components of the
518 : ! state as well, but this is not yet implemented.
519 :
520 : ! Check for NaN first to avoid any IEEE exceptions.
521 :
522 3559680 : if (present(name)) then
523 : msg = "NaN produced in physics_state by package "// &
524 3559680 : trim(name)//"."
525 : else
526 0 : msg = "NaN found in physics_state."
527 : end if
528 :
529 : ! 1-D variables
530 0 : call shr_assert_in_domain(state%ps(:ncol), is_nan=.false., &
531 3559680 : varname="state%ps", msg=msg)
532 0 : call shr_assert_in_domain(state%psdry(:ncol), is_nan=.false., &
533 3559680 : varname="state%psdry", msg=msg)
534 0 : call shr_assert_in_domain(state%phis(:ncol), is_nan=.false., &
535 3559680 : varname="state%phis", msg=msg)
536 0 : call shr_assert_in_domain(state%te_ini(:ncol,:), is_nan=.false., &
537 3559680 : varname="state%te_ini", msg=msg)
538 0 : call shr_assert_in_domain(state%te_cur(:ncol,:), is_nan=.false., &
539 3559680 : varname="state%te_cur", msg=msg)
540 0 : call shr_assert_in_domain(state%tw_ini(:ncol), is_nan=.false., &
541 3559680 : varname="state%tw_ini", msg=msg)
542 0 : call shr_assert_in_domain(state%tw_cur(:ncol), is_nan=.false., &
543 3559680 : varname="state%tw_cur", msg=msg)
544 0 : call shr_assert_in_domain(state%temp_ini(:ncol,:), is_nan=.false., &
545 3559680 : varname="state%temp_ini", msg=msg)
546 0 : call shr_assert_in_domain(state%z_ini(:ncol,:), is_nan=.false., &
547 3559680 : varname="state%z_ini", msg=msg)
548 :
549 : ! 2-D variables (at midpoints)
550 0 : call shr_assert_in_domain(state%t(:ncol,:), is_nan=.false., &
551 3559680 : varname="state%t", msg=msg)
552 0 : call shr_assert_in_domain(state%u(:ncol,:), is_nan=.false., &
553 3559680 : varname="state%u", msg=msg)
554 0 : call shr_assert_in_domain(state%v(:ncol,:), is_nan=.false., &
555 3559680 : varname="state%v", msg=msg)
556 0 : call shr_assert_in_domain(state%s(:ncol,:), is_nan=.false., &
557 3559680 : varname="state%s", msg=msg)
558 0 : call shr_assert_in_domain(state%omega(:ncol,:), is_nan=.false., &
559 3559680 : varname="state%omega", msg=msg)
560 0 : call shr_assert_in_domain(state%pmid(:ncol,:), is_nan=.false., &
561 3559680 : varname="state%pmid", msg=msg)
562 0 : call shr_assert_in_domain(state%pmiddry(:ncol,:), is_nan=.false., &
563 3559680 : varname="state%pmiddry", msg=msg)
564 0 : call shr_assert_in_domain(state%pdel(:ncol,:), is_nan=.false., &
565 3559680 : varname="state%pdel", msg=msg)
566 0 : call shr_assert_in_domain(state%pdeldry(:ncol,:), is_nan=.false., &
567 3559680 : varname="state%pdeldry", msg=msg)
568 0 : call shr_assert_in_domain(state%rpdel(:ncol,:), is_nan=.false., &
569 3559680 : varname="state%rpdel", msg=msg)
570 0 : call shr_assert_in_domain(state%rpdeldry(:ncol,:), is_nan=.false., &
571 3559680 : varname="state%rpdeldry", msg=msg)
572 0 : call shr_assert_in_domain(state%lnpmid(:ncol,:), is_nan=.false., &
573 3559680 : varname="state%lnpmid", msg=msg)
574 0 : call shr_assert_in_domain(state%lnpmiddry(:ncol,:), is_nan=.false., &
575 3559680 : varname="state%lnpmiddry", msg=msg)
576 0 : call shr_assert_in_domain(state%exner(:ncol,:), is_nan=.false., &
577 3559680 : varname="state%exner", msg=msg)
578 0 : call shr_assert_in_domain(state%zm(:ncol,:), is_nan=.false., &
579 3559680 : varname="state%zm", msg=msg)
580 :
581 : ! 2-D variables (at interfaces)
582 0 : call shr_assert_in_domain(state%pint(:ncol,:), is_nan=.false., &
583 3559680 : varname="state%pint", msg=msg)
584 0 : call shr_assert_in_domain(state%pintdry(:ncol,:), is_nan=.false., &
585 3559680 : varname="state%pintdry", msg=msg)
586 0 : call shr_assert_in_domain(state%lnpint(:ncol,:), is_nan=.false., &
587 3559680 : varname="state%lnpint", msg=msg)
588 0 : call shr_assert_in_domain(state%lnpintdry(:ncol,:), is_nan=.false., &
589 3559680 : varname="state%lnpintdry", msg=msg)
590 0 : call shr_assert_in_domain(state%zi(:ncol,:), is_nan=.false., &
591 3559680 : varname="state%zi", msg=msg)
592 :
593 : ! 3-D variables
594 0 : call shr_assert_in_domain(state%q(:ncol,:,:), is_nan=.false., &
595 3559680 : varname="state%q", msg=msg)
596 :
597 : ! Now run other checks (i.e. values are finite and within a range that
598 : ! is physically meaningful).
599 :
600 3559680 : if (present(name)) then
601 : msg = "Invalid value produced in physics_state by package "// &
602 3559680 : trim(name)//"."
603 : else
604 0 : msg = "Invalid value found in physics_state."
605 : end if
606 :
607 : ! 1-D variables
608 : call shr_assert_in_domain(state%ps(:ncol), lt=posinf_r8, gt=0._r8, &
609 3559680 : varname="state%ps", msg=msg)
610 : call shr_assert_in_domain(state%psdry(:ncol), lt=posinf_r8, gt=0._r8, &
611 3559680 : varname="state%psdry", msg=msg)
612 : call shr_assert_in_domain(state%phis(:ncol), lt=posinf_r8, gt=neginf_r8, &
613 3559680 : varname="state%phis", msg=msg)
614 : call shr_assert_in_domain(state%te_ini(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
615 3559680 : varname="state%te_ini", msg=msg)
616 : call shr_assert_in_domain(state%te_cur(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
617 3559680 : varname="state%te_cur", msg=msg)
618 : call shr_assert_in_domain(state%tw_ini(:ncol), lt=posinf_r8, gt=neginf_r8, &
619 3559680 : varname="state%tw_ini", msg=msg)
620 : call shr_assert_in_domain(state%tw_cur(:ncol), lt=posinf_r8, gt=neginf_r8, &
621 3559680 : varname="state%tw_cur", msg=msg)
622 : call shr_assert_in_domain(state%temp_ini(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
623 3559680 : varname="state%temp_ini", msg=msg)
624 : call shr_assert_in_domain(state%z_ini(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
625 3559680 : varname="state%z_ini", msg=msg)
626 :
627 : ! 2-D variables (at midpoints)
628 : call shr_assert_in_domain(state%t(:ncol,:), lt=posinf_r8, gt=0._r8, &
629 3559680 : varname="state%t", msg=msg)
630 : call shr_assert_in_domain(state%u(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
631 3559680 : varname="state%u", msg=msg)
632 : call shr_assert_in_domain(state%v(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
633 3559680 : varname="state%v", msg=msg)
634 : call shr_assert_in_domain(state%s(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
635 3559680 : varname="state%s", msg=msg)
636 : call shr_assert_in_domain(state%omega(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
637 3559680 : varname="state%omega", msg=msg)
638 : call shr_assert_in_domain(state%pmid(:ncol,:), lt=posinf_r8, gt=0._r8, &
639 3559680 : varname="state%pmid", msg=msg)
640 : call shr_assert_in_domain(state%pmiddry(:ncol,:), lt=posinf_r8, gt=0._r8, &
641 3559680 : varname="state%pmiddry", msg=msg)
642 : call shr_assert_in_domain(state%pdel(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
643 3559680 : varname="state%pdel", msg=msg)
644 : call shr_assert_in_domain(state%pdeldry(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
645 3559680 : varname="state%pdeldry", msg=msg)
646 : call shr_assert_in_domain(state%rpdel(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
647 3559680 : varname="state%rpdel", msg=msg)
648 : call shr_assert_in_domain(state%rpdeldry(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
649 3559680 : varname="state%rpdeldry", msg=msg)
650 : call shr_assert_in_domain(state%lnpmid(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
651 3559680 : varname="state%lnpmid", msg=msg)
652 : call shr_assert_in_domain(state%lnpmiddry(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
653 3559680 : varname="state%lnpmiddry", msg=msg)
654 : call shr_assert_in_domain(state%exner(:ncol,:), lt=posinf_r8, gt=0._r8, &
655 3559680 : varname="state%exner", msg=msg)
656 : call shr_assert_in_domain(state%zm(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
657 3559680 : varname="state%zm", msg=msg)
658 :
659 : ! 2-D variables (at interfaces)
660 : call shr_assert_in_domain(state%pint(:ncol,:), lt=posinf_r8, gt=0._r8, &
661 3559680 : varname="state%pint", msg=msg)
662 : call shr_assert_in_domain(state%pintdry(:ncol,:), lt=posinf_r8, gt=0._r8, &
663 3559680 : varname="state%pintdry", msg=msg)
664 : call shr_assert_in_domain(state%lnpint(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
665 3559680 : varname="state%lnpint", msg=msg)
666 : call shr_assert_in_domain(state%lnpintdry(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
667 3559680 : varname="state%lnpintdry", msg=msg)
668 : call shr_assert_in_domain(state%zi(:ncol,:), lt=posinf_r8, gt=neginf_r8, &
669 3559680 : varname="state%zi", msg=msg)
670 :
671 : ! 3-D variables
672 722615040 : do m = 1,pcnst
673 : call shr_assert_in_domain(state%q(:ncol,:,m), lt=posinf_r8, gt=neginf_r8, &
674 722615040 : varname="state%q ("//trim(cnst_name(m))//")", msg=msg)
675 : end do
676 :
677 3559680 : end subroutine physics_state_check
678 :
679 : !===============================================================================
680 :
681 2096640 : subroutine physics_ptend_sum(ptend, ptend_sum, ncol)
682 : !-----------------------------------------------------------------------
683 : ! Add ptend fields to ptend_sum for ptend logical flags = .true.
684 : ! Where ptend logical flags = .false, don't change ptend_sum
685 : !-----------------------------------------------------------------------
686 :
687 : !------------------------------Arguments--------------------------------
688 : type(physics_ptend), intent(in) :: ptend ! New parameterization tendencies
689 : type(physics_ptend), intent(inout) :: ptend_sum ! Sum of incoming ptend_sum and ptend
690 : integer, intent(in) :: ncol ! number of columns
691 :
692 : !---------------------------Local storage-------------------------------
693 : integer :: i,k,m ! column,level,constituent indices
694 : integer :: psetcols ! maximum number of columns
695 : integer :: ierr = 0
696 :
697 : !-----------------------------------------------------------------------
698 2096640 : if (ptend%psetcols /= ptend_sum%psetcols) then
699 0 : call endrun('physics_ptend_sum error: ptend and ptend_sum must have the same value for psetcols')
700 : end if
701 :
702 2096640 : if (ncol > ptend_sum%psetcols) then
703 0 : call endrun('physics_ptend_sum error: ncol must be less than or equal to psetcols')
704 : end if
705 :
706 2096640 : psetcols = ptend_sum%psetcols
707 :
708 2096640 : ptend_sum%top_level = ptend%top_level
709 2096640 : ptend_sum%bot_level = ptend%bot_level
710 :
711 : ! Update u,v fields
712 2096640 : if(ptend%lu) then
713 564480 : if (.not. allocated(ptend_sum%u)) then
714 967680 : allocate(ptend_sum%u(psetcols,pver), stat=ierr)
715 322560 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%u')
716 384168960 : ptend_sum%u=0.0_r8
717 :
718 967680 : allocate(ptend_sum%taux_srf(psetcols), stat=ierr)
719 322560 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%taux_srf')
720 5483520 : ptend_sum%taux_srf=0.0_r8
721 :
722 645120 : allocate(ptend_sum%taux_top(psetcols), stat=ierr)
723 322560 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%taux_top')
724 5483520 : ptend_sum%taux_top=0.0_r8
725 : end if
726 564480 : ptend_sum%lu = .true.
727 :
728 40078080 : do k = ptend%top_level, ptend%bot_level
729 609073920 : do i = 1, ncol
730 608509440 : ptend_sum%u(i,k) = ptend_sum%u(i,k) + ptend%u(i,k)
731 : end do
732 : end do
733 8692992 : do i = 1, ncol
734 8128512 : ptend_sum%taux_srf(i) = ptend_sum%taux_srf(i) + ptend%taux_srf(i)
735 8692992 : ptend_sum%taux_top(i) = ptend_sum%taux_top(i) + ptend%taux_top(i)
736 : end do
737 : end if
738 :
739 2096640 : if(ptend%lv) then
740 564480 : if (.not. allocated(ptend_sum%v)) then
741 967680 : allocate(ptend_sum%v(psetcols,pver), stat=ierr)
742 322560 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%v')
743 384168960 : ptend_sum%v=0.0_r8
744 :
745 967680 : allocate(ptend_sum%tauy_srf(psetcols), stat=ierr)
746 322560 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%tauy_srf')
747 5483520 : ptend_sum%tauy_srf=0.0_r8
748 :
749 645120 : allocate(ptend_sum%tauy_top(psetcols), stat=ierr)
750 322560 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%tauy_top')
751 5483520 : ptend_sum%tauy_top=0.0_r8
752 : end if
753 564480 : ptend_sum%lv = .true.
754 :
755 40078080 : do k = ptend%top_level, ptend%bot_level
756 609073920 : do i = 1, ncol
757 608509440 : ptend_sum%v(i,k) = ptend_sum%v(i,k) + ptend%v(i,k)
758 : end do
759 : end do
760 8692992 : do i = 1, ncol
761 8128512 : ptend_sum%tauy_srf(i) = ptend_sum%tauy_srf(i) + ptend%tauy_srf(i)
762 8692992 : ptend_sum%tauy_top(i) = ptend_sum%tauy_top(i) + ptend%tauy_top(i)
763 : end do
764 : end if
765 :
766 :
767 2096640 : if(ptend%ls) then
768 1290240 : if (.not. allocated(ptend_sum%s)) then
769 1451520 : allocate(ptend_sum%s(psetcols,pver), stat=ierr)
770 483840 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%s')
771 576253440 : ptend_sum%s=0.0_r8
772 :
773 1451520 : allocate(ptend_sum%hflux_srf(psetcols), stat=ierr)
774 483840 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%hflux_srf')
775 8225280 : ptend_sum%hflux_srf=0.0_r8
776 :
777 967680 : allocate(ptend_sum%hflux_top(psetcols), stat=ierr)
778 483840 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%hflux_top')
779 8225280 : ptend_sum%hflux_top=0.0_r8
780 : end if
781 1290240 : ptend_sum%ls = .true.
782 :
783 91607040 : do k = ptend%top_level, ptend%bot_level
784 1392168960 : do i = 1, ncol
785 1390878720 : ptend_sum%s(i,k) = ptend_sum%s(i,k) + ptend%s(i,k)
786 : end do
787 : end do
788 19869696 : do i = 1, ncol
789 18579456 : ptend_sum%hflux_srf(i) = ptend_sum%hflux_srf(i) + ptend%hflux_srf(i)
790 19869696 : ptend_sum%hflux_top(i) = ptend_sum%hflux_top(i) + ptend%hflux_top(i)
791 : end do
792 : end if
793 :
794 43384320 : if (any(ptend%lq(:))) then
795 :
796 2016000 : if (.not. allocated(ptend_sum%q)) then
797 2903040 : allocate(ptend_sum%q(psetcols,pver,pcnst), stat=ierr)
798 725760 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%q')
799 >17460*10^7 : ptend_sum%q=0.0_r8
800 :
801 2177280 : allocate(ptend_sum%cflx_srf(psetcols,pcnst), stat=ierr)
802 725760 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%cflx_srf')
803 2492985600 : ptend_sum%cflx_srf=0.0_r8
804 :
805 1451520 : allocate(ptend_sum%cflx_top(psetcols,pcnst), stat=ierr)
806 725760 : if ( ierr /= 0 ) call endrun('physics_ptend_sum error: allocation error for ptend_sum%cflx_top')
807 2492985600 : ptend_sum%cflx_top=0.0_r8
808 : end if
809 :
810 409248000 : do m = 1, pcnst
811 409248000 : if(ptend%lq(m)) then
812 118218240 : ptend_sum%lq(m) = .true.
813 8393495040 : do k = ptend%top_level, ptend%bot_level
814 >12755*10^7 : do i = 1,ncol
815 >12743*10^7 : ptend_sum%q(i,k,m) = ptend_sum%q(i,k,m) + ptend%q(i,k,m)
816 : end do
817 : end do
818 1820560896 : do i = 1,ncol
819 1702342656 : ptend_sum%cflx_srf(i,m) = ptend_sum%cflx_srf(i,m) + ptend%cflx_srf(i,m)
820 1820560896 : ptend_sum%cflx_top(i,m) = ptend_sum%cflx_top(i,m) + ptend%cflx_top(i,m)
821 : end do
822 : end if
823 : end do
824 :
825 : end if
826 :
827 2096640 : end subroutine physics_ptend_sum
828 :
829 : !===============================================================================
830 :
831 725760 : subroutine physics_ptend_scale(ptend, fac, ncol)
832 : !-----------------------------------------------------------------------
833 : ! Scale ptend fields for ptend logical flags = .true.
834 : ! Where ptend logical flags = .false, don't change ptend.
835 : !
836 : ! Assumes that input ptend is valid (e.g. that
837 : ! ptend%lu .eqv. allocated(ptend%u)), and therefore
838 : ! does not check allocation status of individual arrays.
839 : !-----------------------------------------------------------------------
840 :
841 : !------------------------------Arguments--------------------------------
842 : type(physics_ptend), intent(inout) :: ptend ! Incoming ptend
843 : real(r8), intent(in) :: fac ! Factor to multiply ptend by.
844 : integer, intent(in) :: ncol ! number of columns
845 :
846 : !---------------------------Local storage-------------------------------
847 : integer :: m ! constituent index
848 :
849 : !-----------------------------------------------------------------------
850 :
851 : ! Update u,v fields
852 725760 : if (ptend%lu) &
853 : call multiply_tendency(ptend%u, &
854 241920 : ptend%taux_srf, ptend%taux_top)
855 :
856 725760 : if (ptend%lv) &
857 : call multiply_tendency(ptend%v, &
858 241920 : ptend%tauy_srf, ptend%tauy_top)
859 :
860 : ! Heat
861 725760 : if (ptend%ls) &
862 : call multiply_tendency(ptend%s, &
863 725760 : ptend%hflux_srf, ptend%hflux_top)
864 :
865 : ! Update constituents
866 147329280 : do m = 1, pcnst
867 146603520 : if (ptend%lq(m)) &
868 : call multiply_tendency(ptend%q(:,:,m), &
869 53948160 : ptend%cflx_srf(:,m), ptend%cflx_top(:,m))
870 : end do
871 :
872 :
873 : contains
874 :
875 54432000 : subroutine multiply_tendency(tend_arr, flx_srf, flx_top)
876 : real(r8), intent(inout) :: tend_arr(:,:) ! Tendency array (pcols, plev)
877 : real(r8), intent(inout) :: flx_srf(:) ! Surface flux (or stress)
878 : real(r8), intent(inout) :: flx_top(:) ! Top-of-model flux (or stress)
879 :
880 : integer :: k
881 :
882 3864672000 : do k = ptend%top_level, ptend%bot_level
883 58732128000 : tend_arr(:ncol,k) = tend_arr(:ncol,k) * fac
884 : end do
885 838252800 : flx_srf(:ncol) = flx_srf(:ncol) * fac
886 838252800 : flx_top(:ncol) = flx_top(:ncol) * fac
887 :
888 54432000 : end subroutine multiply_tendency
889 :
890 : end subroutine physics_ptend_scale
891 :
892 : !===============================================================================
893 :
894 0 : subroutine physics_ptend_copy(ptend, ptend_cp)
895 :
896 : !-----------------------------------------------------------------------
897 : ! Copy a physics_ptend object. Allocate ptend_cp internally before copy.
898 : !-----------------------------------------------------------------------
899 :
900 : type(physics_ptend), intent(in) :: ptend ! ptend source
901 : type(physics_ptend), intent(out) :: ptend_cp ! copy of ptend
902 :
903 : !-----------------------------------------------------------------------
904 :
905 0 : ptend_cp%name = ptend%name
906 :
907 0 : ptend_cp%ls = ptend%ls
908 0 : ptend_cp%lu = ptend%lu
909 0 : ptend_cp%lv = ptend%lv
910 0 : ptend_cp%lq = ptend%lq
911 :
912 0 : call physics_ptend_alloc(ptend_cp, ptend%psetcols)
913 :
914 0 : ptend_cp%top_level = ptend%top_level
915 0 : ptend_cp%bot_level = ptend%bot_level
916 :
917 0 : if (ptend_cp%ls) then
918 0 : ptend_cp%s = ptend%s
919 0 : ptend_cp%hflux_srf = ptend%hflux_srf
920 0 : ptend_cp%hflux_top = ptend%hflux_top
921 : end if
922 :
923 0 : if (ptend_cp%lu) then
924 0 : ptend_cp%u = ptend%u
925 0 : ptend_cp%taux_srf = ptend%taux_srf
926 0 : ptend_cp%taux_top = ptend%taux_top
927 : end if
928 :
929 0 : if (ptend_cp%lv) then
930 0 : ptend_cp%v = ptend%v
931 0 : ptend_cp%tauy_srf = ptend%tauy_srf
932 0 : ptend_cp%tauy_top = ptend%tauy_top
933 : end if
934 :
935 0 : if (any(ptend_cp%lq(:))) then
936 0 : ptend_cp%q = ptend%q
937 0 : ptend_cp%cflx_srf = ptend%cflx_srf
938 0 : ptend_cp%cflx_top = ptend%cflx_top
939 : end if
940 :
941 0 : end subroutine physics_ptend_copy
942 :
943 : !===============================================================================
944 :
945 3083520 : subroutine physics_ptend_reset(ptend)
946 : !-----------------------------------------------------------------------
947 : ! Reset the parameterization tendency structure to "empty"
948 : !-----------------------------------------------------------------------
949 :
950 : !------------------------------Arguments--------------------------------
951 : type(physics_ptend), intent(inout) :: ptend ! Parameterization tendencies
952 : !-----------------------------------------------------------------------
953 :
954 3083520 : if(ptend%ls) then
955 2268426240 : ptend%s = 0._r8
956 32378880 : ptend%hflux_srf = 0._r8
957 32378880 : ptend%hflux_top = 0._r8
958 : endif
959 3083520 : if(ptend%lu) then
960 914688000 : ptend%u = 0._r8
961 13056000 : ptend%taux_srf = 0._r8
962 13056000 : ptend%taux_top = 0._r8
963 : endif
964 3083520 : if(ptend%lv) then
965 827792640 : ptend%v = 0._r8
966 11815680 : ptend%tauy_srf = 0._r8
967 11815680 : ptend%tauy_top = 0._r8
968 : endif
969 152605440 : if(any (ptend%lq(:))) then
970 >61158*10^7 : ptend%q = 0._r8
971 8732044800 : ptend%cflx_srf = 0._r8
972 8732044800 : ptend%cflx_top = 0._r8
973 : end if
974 :
975 3083520 : ptend%top_level = 1
976 3083520 : ptend%bot_level = pver
977 :
978 3083520 : return
979 : end subroutine physics_ptend_reset
980 :
981 : !===============================================================================
982 786535680 : subroutine physics_ptend_init(ptend, psetcols, name, ls, lu, lv, lq)
983 : !-----------------------------------------------------------------------
984 : ! Allocate the fields in the structure which are specified.
985 : ! Initialize the parameterization tendency structure to "empty"
986 : !-----------------------------------------------------------------------
987 :
988 : !------------------------------Arguments--------------------------------
989 : type(physics_ptend), intent(out) :: ptend ! Parameterization tendencies
990 : integer, intent(in) :: psetcols ! maximum number of columns
991 : character(len=*) :: name ! optional name of parameterization which produced tendencies.
992 : logical, optional :: ls ! if true, then fields to support dsdt are allocated
993 : logical, optional :: lu ! if true, then fields to support dudt are allocated
994 : logical, optional :: lv ! if true, then fields to support dvdt are allocated
995 : logical, dimension(pcnst),optional :: lq ! if true, then fields to support dqdt are allocated
996 :
997 : !-----------------------------------------------------------------------
998 :
999 3874560 : if (allocated(ptend%s)) then
1000 0 : call endrun(' physics_ptend_init: ptend should not be allocated before calling this routine')
1001 : end if
1002 :
1003 3874560 : ptend%name = name
1004 3874560 : ptend%psetcols = psetcols
1005 :
1006 : ! If no fields being stored, initialize all values to appropriate nulls and return
1007 3874560 : if (.not. present(ls) .and. .not. present(lu) .and. .not. present(lv) .and. .not. present(lq) ) then
1008 791040 : ptend%ls = .false.
1009 791040 : ptend%lu = .false.
1010 791040 : ptend%lv = .false.
1011 160581120 : ptend%lq(:) = .false.
1012 791040 : ptend%top_level = 1
1013 791040 : ptend%bot_level = pver
1014 791040 : return
1015 : end if
1016 :
1017 3083520 : if (present(ls)) then
1018 1904640 : ptend%ls = ls
1019 : else
1020 1178880 : ptend%ls = .false.
1021 : end if
1022 :
1023 3083520 : if (present(lu)) then
1024 768000 : ptend%lu = lu
1025 : else
1026 2315520 : ptend%lu = .false.
1027 : end if
1028 :
1029 3083520 : if (present(lv)) then
1030 695040 : ptend%lv = lv
1031 : else
1032 2388480 : ptend%lv = .false.
1033 : end if
1034 :
1035 3083520 : if (present(lq)) then
1036 516042240 : ptend%lq(:) = lq(:)
1037 : else
1038 109912320 : ptend%lq(:) = .false.
1039 : end if
1040 :
1041 3083520 : call physics_ptend_alloc(ptend, psetcols)
1042 :
1043 3083520 : call physics_ptend_reset(ptend)
1044 :
1045 3083520 : return
1046 3874560 : end subroutine physics_ptend_init
1047 :
1048 : !===============================================================================
1049 :
1050 7680 : subroutine physics_state_set_grid(lchnk, phys_state)
1051 : !-----------------------------------------------------------------------
1052 : ! Set the grid components of the physics_state object
1053 : !-----------------------------------------------------------------------
1054 :
1055 : integer, intent(in) :: lchnk
1056 : type(physics_state), intent(inout) :: phys_state
1057 :
1058 : ! local variables
1059 : integer :: i, ncol
1060 : real(r8) :: rlon(pcols)
1061 : real(r8) :: rlat(pcols)
1062 :
1063 : !-----------------------------------------------------------------------
1064 : ! get_ncols_p requires a state which does not have sub-columns
1065 7680 : if (phys_state%psetcols .ne. pcols) then
1066 0 : call endrun('physics_state_set_grid: cannot pass in a state which has sub-columns')
1067 : end if
1068 :
1069 7680 : ncol = get_ncols_p(lchnk)
1070 :
1071 7680 : if(ncol<=0) then
1072 0 : write(iulog,*) lchnk, ncol
1073 0 : call endrun('physics_state_set_grid')
1074 : end if
1075 :
1076 7680 : call get_rlon_all_p(lchnk, ncol, rlon)
1077 7680 : call get_rlat_all_p(lchnk, ncol, rlat)
1078 7680 : phys_state%ncol = ncol
1079 7680 : phys_state%lchnk = lchnk
1080 118272 : do i=1,ncol
1081 110592 : phys_state%lat(i) = rlat(i)
1082 118272 : phys_state%lon(i) = rlon(i)
1083 : end do
1084 7680 : call init_geo_unique(phys_state,ncol)
1085 :
1086 7680 : end subroutine physics_state_set_grid
1087 :
1088 :
1089 7680 : subroutine init_geo_unique(phys_state,ncol)
1090 : integer, intent(in) :: ncol
1091 : type(physics_state), intent(inout) :: phys_state
1092 : logical :: match
1093 : integer :: i, j, ulatcnt, uloncnt
1094 :
1095 130560 : phys_state%ulat=-999._r8
1096 130560 : phys_state%ulon=-999._r8
1097 130560 : phys_state%latmapback=0
1098 130560 : phys_state%lonmapback=0
1099 : match=.false.
1100 7680 : ulatcnt=0
1101 7680 : uloncnt=0
1102 7680 : match=.false.
1103 118272 : do i=1,ncol
1104 852480 : do j=1,ulatcnt
1105 852480 : if(phys_state%lat(i) .eq. phys_state%ulat(j)) then
1106 0 : match=.true.
1107 0 : phys_state%latmapback(i)=j
1108 : end if
1109 : end do
1110 110592 : if(.not. match) then
1111 110592 : ulatcnt=ulatcnt+1
1112 110592 : phys_state%ulat(ulatcnt)=phys_state%lat(i)
1113 110592 : phys_state%latmapback(i)=ulatcnt
1114 : end if
1115 :
1116 110592 : match=.false.
1117 671616 : do j=1,uloncnt
1118 671616 : if(phys_state%lon(i) .eq. phys_state%ulon(j)) then
1119 32160 : match=.true.
1120 32160 : phys_state%lonmapback(i)=j
1121 : end if
1122 : end do
1123 110592 : if(.not. match) then
1124 78432 : uloncnt=uloncnt+1
1125 78432 : phys_state%ulon(uloncnt)=phys_state%lon(i)
1126 78432 : phys_state%lonmapback(i)=uloncnt
1127 : end if
1128 118272 : match=.false.
1129 :
1130 : end do
1131 7680 : phys_state%uloncnt=uloncnt
1132 7680 : phys_state%ulatcnt=ulatcnt
1133 :
1134 7680 : call get_gcol_all_p(phys_state%lchnk,pcols,phys_state%cid)
1135 :
1136 :
1137 7680 : end subroutine init_geo_unique
1138 :
1139 : !===============================================================================
1140 0 : subroutine physics_cnst_limit(state)
1141 : type(physics_state), intent(inout) :: state
1142 :
1143 : integer :: i,k, ncol
1144 :
1145 : real(r8) :: mmrSum_O_O2_H ! Sum of mass mixing ratios for O, O2, and H
1146 : real(r8), parameter :: mmrMin=1.e-20_r8 ! lower limit of o2, o, and h mixing ratios
1147 : real(r8), parameter :: N2mmrMin=1.e-6_r8 ! lower limit of N2 mass mixing ratio
1148 : real(r8), parameter :: H2lim=6.e-5_r8 ! H2 limiter: 10x global H2 MMR (Roble, 1995)
1149 : integer :: ixo, ixo2, ixh, ixh2
1150 :
1151 0 : if ( waccmx_is('ionosphere') .or. waccmx_is('neutral') ) then
1152 0 : call cnst_get_ind('O', ixo)
1153 0 : call cnst_get_ind('O2', ixo2)
1154 0 : call cnst_get_ind('H', ixh)
1155 0 : call cnst_get_ind('H2', ixh2)
1156 :
1157 0 : ncol = state%ncol
1158 :
1159 : !------------------------------------------------------------
1160 : ! Ensure N2 = 1-(O2 + O + H) mmr is greater than 0
1161 : ! Check for unusually large H2 values and set to lower value.
1162 : !------------------------------------------------------------
1163 :
1164 0 : do k=1,pver
1165 0 : do i=1,ncol
1166 :
1167 0 : if (state%q(i,k,ixo) < mmrMin) state%q(i,k,ixo) = mmrMin
1168 0 : if (state%q(i,k,ixo2) < mmrMin) state%q(i,k,ixo2) = mmrMin
1169 :
1170 0 : mmrSum_O_O2_H = state%q(i,k,ixo)+state%q(i,k,ixo2)+state%q(i,k,ixh)
1171 :
1172 0 : if ((1._r8-mmrMin-mmrSum_O_O2_H) < 0._r8) then
1173 :
1174 0 : state%q(i,k,ixo) = state%q(i,k,ixo) * (1._r8 - N2mmrMin) / mmrSum_O_O2_H
1175 :
1176 0 : state%q(i,k,ixo2) = state%q(i,k,ixo2) * (1._r8 - N2mmrMin) / mmrSum_O_O2_H
1177 :
1178 0 : state%q(i,k,ixh) = state%q(i,k,ixh) * (1._r8 - N2mmrMin) / mmrSum_O_O2_H
1179 :
1180 : endif
1181 :
1182 0 : if(state%q(i,k,ixh2) > H2lim) then
1183 0 : state%q(i,k,ixh2) = H2lim
1184 : endif
1185 :
1186 : end do
1187 : end do
1188 :
1189 : end if
1190 0 : end subroutine physics_cnst_limit
1191 :
1192 : !===============================================================================
1193 72960 : subroutine physics_dme_adjust(state, tend, qini, liqini, iceini, dt)
1194 : use air_composition, only: dry_air_species_num,thermodynamic_active_species_num
1195 : use air_composition, only: thermodynamic_active_species_idx
1196 : use dycore, only: dycore_is
1197 : !-----------------------------------------------------------------------
1198 : !
1199 : ! Purpose: Adjust the dry mass in each layer back to the value of physics input state
1200 : !
1201 : ! Method: Conserve the integrated mass, momentum and total energy in each layer
1202 : ! by scaling the specific mass of consituents, specific momentum (velocity)
1203 : ! and specific total energy by the relative change in layer mass. Solve for
1204 : ! the new temperature by subtracting the new kinetic energy from total energy
1205 : ! and inverting the hydrostatic equation
1206 : !
1207 : ! The mass in each layer is modified, changing the relationship of the layer
1208 : ! interfaces and midpoints to the surface pressure. The result is no longer in
1209 : ! the original hybrid coordinate.
1210 : !
1211 : ! Author: Byron Boville
1212 :
1213 : ! !REVISION HISTORY:
1214 : ! 03.03.28 Boville Created, partly from code by Lin in p_d_adjust
1215 : !
1216 : !-----------------------------------------------------------------------
1217 :
1218 : implicit none
1219 : !
1220 : ! Arguments
1221 : !
1222 : type(physics_state), intent(inout) :: state
1223 : type(physics_tend ), intent(inout) :: tend
1224 : real(r8), intent(in ) :: qini(pcols,pver) ! initial specific humidity
1225 : real(r8), intent(in ) :: liqini(pcols,pver) ! initial total liquid
1226 : real(r8), intent(in ) :: iceini(pcols,pver) ! initial total ice
1227 : real(r8), intent(in ) :: dt ! model physics timestep
1228 : !
1229 : !---------------------------Local workspace-----------------------------
1230 : !
1231 : integer :: lchnk ! chunk identifier
1232 : integer :: ncol ! number of atmospheric columns
1233 : integer :: k,m ! Longitude, level indices
1234 : real(r8) :: fdq(pcols) ! mass adjustment factor
1235 : real(r8) :: te(pcols) ! total energy in a layer
1236 : real(r8) :: utmp(pcols) ! temp variable for recalculating the initial u values
1237 : real(r8) :: vtmp(pcols) ! temp variable for recalculating the initial v values
1238 :
1239 : real(r8) :: zvirv(pcols,pver) ! Local zvir array pointer
1240 :
1241 : real(r8) :: tot_water (pcols,2) ! total water (initial, present)
1242 : real(r8) :: tot_water_chg(pcols) ! total water change
1243 :
1244 :
1245 : real(r8),allocatable :: cpairv_loc(:,:)
1246 : integer :: m_cnst
1247 : !
1248 : !-----------------------------------------------------------------------
1249 :
1250 72960 : if (state%psetcols .ne. pcols) then
1251 0 : call endrun('physics_dme_adjust: cannot pass in a state which has sub-columns')
1252 : end if
1253 :
1254 72960 : lchnk = state%lchnk
1255 72960 : ncol = state%ncol
1256 :
1257 : ! adjust dry mass in each layer back to input value, while conserving
1258 : ! constituents, momentum, and total energy
1259 1123584 : state%ps(:ncol) = state%pint(:ncol,1)
1260 :
1261 : !
1262 : ! original code for backwards compatability with FV
1263 : !
1264 72960 : if (.not.(dycore_is('MPAS') .or. dycore_is('SE'))) then
1265 5180160 : do k = 1, pver
1266 :
1267 : ! adjusment factor is just change in water vapor
1268 78650880 : fdq(:ncol) = 1._r8 + state%q(:ncol,k,1) - qini(:ncol,k)
1269 :
1270 : ! adjust constituents to conserve mass in each layer
1271 1036761600 : do m = 1, pcnst
1272 15892584960 : state%q(:ncol,k,m) = state%q(:ncol,k,m) / fdq(:ncol)
1273 : end do
1274 : ! compute new total pressure variables
1275 78650880 : state%pdel (:ncol,k ) = state%pdel(:ncol,k ) * fdq(:ncol)
1276 78650880 : state%ps(:ncol) = state%ps(:ncol) + state%pdel(:ncol,k)
1277 78650880 : state%pint (:ncol,k+1) = state%pint(:ncol,k ) + state%pdel(:ncol,k)
1278 78650880 : state%lnpint(:ncol,k+1) = log(state%pint(:ncol,k+1))
1279 78723840 : state%rpdel (:ncol,k ) = 1._r8/ state%pdel(:ncol,k )
1280 : end do
1281 : else
1282 0 : do k = 1, pver
1283 0 : tot_water(:ncol,1) = qini(:ncol,k) +liqini(:ncol,k)+iceini(:ncol,k) !initial total H2O
1284 0 : tot_water(:ncol,2) = 0.0_r8
1285 0 : do m_cnst=dry_air_species_num+1,thermodynamic_active_species_num
1286 0 : m = thermodynamic_active_species_idx(m_cnst)
1287 0 : tot_water(:ncol,2) = tot_water(:ncol,2)+state%q(:ncol,k,m)
1288 : end do
1289 0 : fdq(:ncol) = 1._r8 + tot_water(:ncol,2) - tot_water(:ncol,1)
1290 : ! adjust constituents to conserve mass in each layer
1291 0 : do m = 1, pcnst
1292 0 : state%q(:ncol,k,m) = state%q(:ncol,k,m) / fdq(:ncol)
1293 : end do
1294 : ! compute new total pressure variables
1295 0 : state%pdel (:ncol,k ) = state%pdel(:ncol,k ) * fdq(:ncol)
1296 0 : state%ps(:ncol) = state%ps(:ncol) + state%pdel(:ncol,k)
1297 0 : state%pint (:ncol,k+1) = state%pint(:ncol,k ) + state%pdel(:ncol,k)
1298 0 : state%lnpint(:ncol,k+1) = log(state%pint(:ncol,k+1))
1299 0 : state%rpdel (:ncol,k ) = 1._r8/ state%pdel(:ncol,k )
1300 : !note that mid-level variables (e.g. pmid) are not recomputed
1301 : end do
1302 : endif
1303 72960 : if ( waccmx_is('ionosphere') .or. waccmx_is('neutral') ) then
1304 0 : zvirv(:,:) = shr_const_rwv / rairv(:,:,state%lchnk) - 1._r8
1305 : else
1306 86895360 : zvirv(:,:) = zvir
1307 : endif
1308 :
1309 72960 : end subroutine physics_dme_adjust
1310 :
1311 : !-----------------------------------------------------------------------
1312 :
1313 : !===============================================================================
1314 960000 : subroutine physics_state_copy(state_in, state_out)
1315 :
1316 : use ppgrid, only: pver, pverp
1317 : use constituents, only: pcnst
1318 :
1319 : implicit none
1320 :
1321 : !
1322 : ! Arguments
1323 : !
1324 : type(physics_state), intent(in) :: state_in
1325 : type(physics_state), intent(out) :: state_out
1326 :
1327 : !
1328 : ! Local variables
1329 : !
1330 : integer i, k, m, ncol
1331 :
1332 : ! Allocate state_out with same subcol dimension as state_in
1333 960000 : call physics_state_alloc ( state_out, state_in%lchnk, state_in%psetcols)
1334 :
1335 960000 : ncol = state_in%ncol
1336 :
1337 960000 : state_out%psetcols = state_in%psetcols
1338 960000 : state_out%ngrdcol = state_in%ngrdcol
1339 960000 : state_out%lchnk = state_in%lchnk
1340 960000 : state_out%ncol = state_in%ncol
1341 960000 : state_out%count = state_in%count
1342 :
1343 14784000 : do i = 1, ncol
1344 13824000 : state_out%lat(i) = state_in%lat(i)
1345 13824000 : state_out%lon(i) = state_in%lon(i)
1346 13824000 : state_out%ps(i) = state_in%ps(i)
1347 14784000 : state_out%phis(i) = state_in%phis(i)
1348 : end do
1349 30528000 : state_out%te_ini(:ncol,:) = state_in%te_ini(:ncol,:)
1350 30528000 : state_out%te_cur(:ncol,:) = state_in%te_cur(:ncol,:)
1351 14784000 : state_out%tw_ini(:ncol) = state_in%tw_ini(:ncol)
1352 14784000 : state_out%tw_cur(:ncol) = state_in%tw_cur(:ncol)
1353 :
1354 68160000 : do k = 1, pver
1355 1035840000 : do i = 1, ncol
1356 967680000 : state_out%temp_ini(i,k) = state_in%temp_ini(i,k)
1357 967680000 : state_out%z_ini(i,k) = state_in%z_ini(i,k)
1358 967680000 : state_out%t(i,k) = state_in%t(i,k)
1359 967680000 : state_out%u(i,k) = state_in%u(i,k)
1360 967680000 : state_out%v(i,k) = state_in%v(i,k)
1361 967680000 : state_out%s(i,k) = state_in%s(i,k)
1362 967680000 : state_out%omega(i,k) = state_in%omega(i,k)
1363 967680000 : state_out%pmid(i,k) = state_in%pmid(i,k)
1364 967680000 : state_out%pdel(i,k) = state_in%pdel(i,k)
1365 967680000 : state_out%rpdel(i,k) = state_in%rpdel(i,k)
1366 967680000 : state_out%lnpmid(i,k) = state_in%lnpmid(i,k)
1367 967680000 : state_out%exner(i,k) = state_in%exner(i,k)
1368 1034880000 : state_out%zm(i,k) = state_in%zm(i,k)
1369 : end do
1370 : end do
1371 :
1372 69120000 : do k = 1, pverp
1373 1050624000 : do i = 1, ncol
1374 981504000 : state_out%pint(i,k) = state_in%pint(i,k)
1375 981504000 : state_out%lnpint(i,k) = state_in%lnpint(i,k)
1376 1049664000 : state_out%zi(i,k) = state_in% zi(i,k)
1377 : end do
1378 : end do
1379 :
1380 :
1381 14784000 : do i = 1, ncol
1382 14784000 : state_out%psdry(i) = state_in%psdry(i)
1383 : end do
1384 68160000 : do k = 1, pver
1385 1035840000 : do i = 1, ncol
1386 967680000 : state_out%lnpmiddry(i,k) = state_in%lnpmiddry(i,k)
1387 967680000 : state_out%pmiddry(i,k) = state_in%pmiddry(i,k)
1388 967680000 : state_out%pdeldry(i,k) = state_in%pdeldry(i,k)
1389 1034880000 : state_out%rpdeldry(i,k) = state_in%rpdeldry(i,k)
1390 : end do
1391 : end do
1392 69120000 : do k = 1, pverp
1393 1050624000 : do i = 1, ncol
1394 981504000 : state_out%pintdry(i,k) = state_in%pintdry(i,k)
1395 1049664000 : state_out%lnpintdry(i,k) = state_in%lnpintdry(i,k)
1396 : end do
1397 : end do
1398 :
1399 194880000 : do m = 1, pcnst
1400 13769280000 : do k = 1, pver
1401 >20923*10^7 : do i = 1, ncol
1402 >20904*10^7 : state_out%q(i,k,m) = state_in%q(i,k,m)
1403 : end do
1404 : end do
1405 : end do
1406 :
1407 960000 : end subroutine physics_state_copy
1408 : !===============================================================================
1409 :
1410 0 : subroutine physics_tend_init(tend)
1411 :
1412 : implicit none
1413 :
1414 : !
1415 : ! Arguments
1416 : !
1417 : type(physics_tend), intent(inout) :: tend
1418 :
1419 : !
1420 : ! Local variables
1421 : !
1422 :
1423 0 : if (.not. allocated(tend%dtdt)) then
1424 0 : call endrun('physics_tend_init: tend must be allocated before it can be initialized')
1425 : end if
1426 :
1427 0 : tend%dtdt = 0._r8
1428 0 : tend%dudt = 0._r8
1429 0 : tend%dvdt = 0._r8
1430 0 : tend%flx_net = 0._r8
1431 0 : tend%te_tnd = 0._r8
1432 0 : tend%tw_tnd = 0._r8
1433 :
1434 0 : end subroutine physics_tend_init
1435 :
1436 : !===============================================================================
1437 :
1438 80640 : subroutine set_state_pdry (state,pdeld_calc)
1439 :
1440 : use ppgrid, only: pver
1441 : implicit none
1442 :
1443 : type(physics_state), intent(inout) :: state
1444 : logical, optional, intent(in) :: pdeld_calc ! .true. do calculate pdeld [default]
1445 : ! .false. don't calculate pdeld
1446 : integer ncol
1447 : integer k
1448 : logical do_pdeld_calc
1449 :
1450 80640 : if ( present(pdeld_calc) ) then
1451 0 : do_pdeld_calc = pdeld_calc
1452 : else
1453 : do_pdeld_calc = .true.
1454 : endif
1455 :
1456 80640 : ncol = state%ncol
1457 :
1458 :
1459 1241856 : state%psdry(:ncol) = state%pint(:ncol,1)
1460 1241856 : state%pintdry(:ncol,1) = state%pint(:ncol,1)
1461 :
1462 80640 : if (do_pdeld_calc) then
1463 5725440 : do k = 1, pver
1464 87010560 : state%pdeldry(:ncol,k) = state%pdel(:ncol,k)*(1._r8-state%q(:ncol,k,1))
1465 : end do
1466 : endif
1467 5725440 : do k = 1, pver
1468 86929920 : state%pintdry(:ncol,k+1) = state%pintdry(:ncol,k)+state%pdeldry(:ncol,k)
1469 86929920 : state%pmiddry(:ncol,k) = (state%pintdry(:ncol,k+1)+state%pintdry(:ncol,k))/2._r8
1470 87010560 : state%psdry(:ncol) = state%psdry(:ncol) + state%pdeldry(:ncol,k)
1471 : end do
1472 :
1473 87010560 : state%rpdeldry(:ncol,:) = 1._r8/state%pdeldry(:ncol,:)
1474 87010560 : state%lnpmiddry(:ncol,:) = log(state%pmiddry(:ncol,:))
1475 88252416 : state%lnpintdry(:ncol,:) = log(state%pintdry(:ncol,:))
1476 :
1477 80640 : end subroutine set_state_pdry
1478 :
1479 : !===============================================================================
1480 :
1481 322560 : subroutine set_wet_to_dry(state, convert_cnst_type)
1482 :
1483 : ! Convert mixing ratios from a wet to dry basis for constituents of type
1484 : ! convert_cnst_type. Constituents are given a type when they are added
1485 : ! to the constituent array by a call to cnst_add during the register
1486 : ! phase of initialization. There are two constituent types: 'wet' for
1487 : ! water species and 'dry' for non-water species.
1488 :
1489 : use constituents, only: pcnst, cnst_type
1490 :
1491 : type(physics_state), intent(inout) :: state
1492 : character(len=3), intent(in) :: convert_cnst_type
1493 :
1494 : ! local variables
1495 : integer m, ncol
1496 : character(len=*), parameter :: sub = 'set_wet_to_dry'
1497 : !-----------------------------------------------------------------------------
1498 :
1499 : ! check input
1500 322560 : if (.not.(convert_cnst_type == 'wet' .or. convert_cnst_type == 'dry')) then
1501 0 : write(iulog,*) sub//': FATAL: convert_cnst_type not recognized: '//convert_cnst_type
1502 0 : call endrun(sub//': FATAL: convert_cnst_type not recognized: '//convert_cnst_type)
1503 : end if
1504 :
1505 322560 : ncol = state%ncol
1506 :
1507 65479680 : do m = 1, pcnst
1508 65479680 : if (cnst_type(m) == convert_cnst_type) then
1509 19664386560 : state%q(:ncol,:,m) = state%q(:ncol,:,m)*state%pdel(:ncol,:)/state%pdeldry(:ncol,:)
1510 : end if
1511 : end do
1512 :
1513 322560 : end subroutine set_wet_to_dry
1514 :
1515 : !===============================================================================
1516 :
1517 72960 : subroutine set_dry_to_wet(state, convert_cnst_type)
1518 :
1519 : ! Convert mixing ratios from a dry to wet basis for constituents of type
1520 : ! convert_cnst_type. Constituents are given a type when they are added
1521 : ! to the constituent array by a call to cnst_add during the register
1522 : ! phase of initialization. There are two constituent types: 'wet' for
1523 : ! water species and 'dry' for non-water species.
1524 :
1525 : use constituents, only: pcnst, cnst_type
1526 :
1527 : type(physics_state), intent(inout) :: state
1528 : character(len=3), intent(in) :: convert_cnst_type
1529 :
1530 : ! local variables
1531 : integer m, ncol
1532 : character(len=*), parameter :: sub = 'set_dry_to_wet'
1533 : !-----------------------------------------------------------------------------
1534 :
1535 : ! check input
1536 72960 : if (.not.(convert_cnst_type == 'wet' .or. convert_cnst_type == 'dry')) then
1537 0 : write(iulog,*) sub//': FATAL: convert_cnst_type not recognized: '//convert_cnst_type
1538 0 : call endrun(sub//': FATAL: convert_cnst_type not recognized: '//convert_cnst_type)
1539 : end if
1540 :
1541 72960 : ncol = state%ncol
1542 :
1543 14810880 : do m = 1, pcnst
1544 14810880 : if (cnst_type(m) == convert_cnst_type) then
1545 14957529600 : state%q(:ncol,:,m) = state%q(:ncol,:,m)*state%pdeldry(:ncol,:)/state%pdel(:ncol,:)
1546 : end if
1547 : end do
1548 :
1549 72960 : end subroutine set_dry_to_wet
1550 :
1551 : !===============================================================================
1552 :
1553 967680 : subroutine physics_state_alloc(state,lchnk,psetcols)
1554 :
1555 : use infnan, only: inf, assignment(=)
1556 :
1557 : ! allocate the individual state components
1558 :
1559 : type(physics_state), intent(inout) :: state
1560 : integer,intent(in) :: lchnk
1561 :
1562 : integer, intent(in) :: psetcols
1563 :
1564 : integer :: ierr=0
1565 :
1566 967680 : state%lchnk = lchnk
1567 967680 : state%psetcols = psetcols
1568 967680 : state%ngrdcol = get_ncols_p(lchnk) ! Number of grid columns
1569 :
1570 : !----------------------------------
1571 : ! Following variables will be overwritten by sub-column generator, if sub-columns are being used
1572 :
1573 : ! state%ncol - is initialized in physics_state_set_grid, if not using sub-columns
1574 :
1575 : !----------------------------------
1576 :
1577 2903040 : allocate(state%lat(psetcols), stat=ierr)
1578 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lat')
1579 :
1580 1935360 : allocate(state%lon(psetcols), stat=ierr)
1581 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lon')
1582 :
1583 1935360 : allocate(state%ps(psetcols), stat=ierr)
1584 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%ps')
1585 :
1586 1935360 : allocate(state%psdry(psetcols), stat=ierr)
1587 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%psdry')
1588 :
1589 1935360 : allocate(state%phis(psetcols), stat=ierr)
1590 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%phis')
1591 :
1592 1935360 : allocate(state%ulat(psetcols), stat=ierr)
1593 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%ulat')
1594 :
1595 1935360 : allocate(state%ulon(psetcols), stat=ierr)
1596 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%ulon')
1597 :
1598 2903040 : allocate(state%t(psetcols,pver), stat=ierr)
1599 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%t')
1600 :
1601 1935360 : allocate(state%u(psetcols,pver), stat=ierr)
1602 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%u')
1603 :
1604 1935360 : allocate(state%v(psetcols,pver), stat=ierr)
1605 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%v')
1606 :
1607 1935360 : allocate(state%s(psetcols,pver), stat=ierr)
1608 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%s')
1609 :
1610 1935360 : allocate(state%omega(psetcols,pver), stat=ierr)
1611 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%omega')
1612 :
1613 1935360 : allocate(state%pmid(psetcols,pver), stat=ierr)
1614 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%pmid')
1615 :
1616 1935360 : allocate(state%pmiddry(psetcols,pver), stat=ierr)
1617 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%pmiddry')
1618 :
1619 1935360 : allocate(state%pdel(psetcols,pver), stat=ierr)
1620 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%pdel')
1621 :
1622 1935360 : allocate(state%pdeldry(psetcols,pver), stat=ierr)
1623 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%pdeldry')
1624 :
1625 1935360 : allocate(state%rpdel(psetcols,pver), stat=ierr)
1626 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%rpdel')
1627 :
1628 1935360 : allocate(state%rpdeldry(psetcols,pver), stat=ierr)
1629 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%rpdeldry')
1630 :
1631 1935360 : allocate(state%lnpmid(psetcols,pver), stat=ierr)
1632 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lnpmid')
1633 :
1634 1935360 : allocate(state%lnpmiddry(psetcols,pver), stat=ierr)
1635 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lnpmiddry')
1636 :
1637 1935360 : allocate(state%exner(psetcols,pver), stat=ierr)
1638 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%exner')
1639 :
1640 1935360 : allocate(state%zm(psetcols,pver), stat=ierr)
1641 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%zm')
1642 :
1643 3870720 : allocate(state%q(psetcols,pver,pcnst), stat=ierr)
1644 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%q')
1645 :
1646 2903040 : allocate(state%pint(psetcols,pver+1), stat=ierr)
1647 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%pint')
1648 :
1649 1935360 : allocate(state%pintdry(psetcols,pver+1), stat=ierr)
1650 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%pintdry')
1651 :
1652 1935360 : allocate(state%lnpint(psetcols,pver+1), stat=ierr)
1653 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lnpint')
1654 :
1655 1935360 : allocate(state%lnpintdry(psetcols,pver+1), stat=ierr)
1656 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lnpintdry')
1657 :
1658 1935360 : allocate(state%zi(psetcols,pver+1), stat=ierr)
1659 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%zi')
1660 :
1661 2903040 : allocate(state%te_ini(psetcols,2), stat=ierr)
1662 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%te_ini')
1663 :
1664 1935360 : allocate(state%te_cur(psetcols,2), stat=ierr)
1665 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%te_cur')
1666 :
1667 1935360 : allocate(state%tw_ini(psetcols), stat=ierr)
1668 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%tw_ini')
1669 :
1670 1935360 : allocate(state%tw_cur(psetcols), stat=ierr)
1671 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%tw_cur')
1672 :
1673 1935360 : allocate(state%temp_ini(psetcols,pver), stat=ierr)
1674 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%temp_ini')
1675 :
1676 1935360 : allocate(state%z_ini(psetcols,pver), stat=ierr)
1677 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%z_ini')
1678 :
1679 2903040 : allocate(state%latmapback(psetcols), stat=ierr)
1680 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%latmapback')
1681 :
1682 1935360 : allocate(state%lonmapback(psetcols), stat=ierr)
1683 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%lonmapback')
1684 :
1685 1935360 : allocate(state%cid(psetcols), stat=ierr)
1686 967680 : if ( ierr /= 0 ) call endrun('physics_state_alloc error: allocation error for state%cid')
1687 :
1688 967680 : state%lat(:) = inf
1689 967680 : state%lon(:) = inf
1690 967680 : state%ulat(:) = inf
1691 967680 : state%ulon(:) = inf
1692 967680 : state%ps(:) = inf
1693 967680 : state%psdry(:) = inf
1694 967680 : state%phis(:) = inf
1695 967680 : state%t(:,:) = inf
1696 967680 : state%u(:,:) = inf
1697 967680 : state%v(:,:) = inf
1698 967680 : state%s(:,:) = inf
1699 967680 : state%omega(:,:) = inf
1700 967680 : state%pmid(:,:) = inf
1701 967680 : state%pmiddry(:,:) = inf
1702 967680 : state%pdel(:,:) = inf
1703 967680 : state%pdeldry(:,:) = inf
1704 967680 : state%rpdel(:,:) = inf
1705 967680 : state%rpdeldry(:,:) = inf
1706 967680 : state%lnpmid(:,:) = inf
1707 967680 : state%lnpmiddry(:,:) = inf
1708 967680 : state%exner(:,:) = inf
1709 967680 : state%zm(:,:) = inf
1710 967680 : state%q(:,:,:) = inf
1711 :
1712 967680 : state%pint(:,:) = inf
1713 967680 : state%pintdry(:,:) = inf
1714 967680 : state%lnpint(:,:) = inf
1715 967680 : state%lnpintdry(:,:) = inf
1716 967680 : state%zi(:,:) = inf
1717 :
1718 967680 : state%te_ini(:,:) = inf
1719 967680 : state%te_cur(:,:) = inf
1720 967680 : state%tw_ini(:) = inf
1721 967680 : state%tw_cur(:) = inf
1722 967680 : state%temp_ini(:,:) = inf
1723 967680 : state%z_ini(:,:) = inf
1724 :
1725 967680 : end subroutine physics_state_alloc
1726 :
1727 : !===============================================================================
1728 :
1729 403200 : subroutine physics_state_dealloc(state)
1730 :
1731 : ! deallocate the individual state components
1732 :
1733 : type(physics_state), intent(inout) :: state
1734 : integer :: ierr = 0
1735 :
1736 403200 : deallocate(state%lat, stat=ierr)
1737 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lat')
1738 :
1739 403200 : deallocate(state%lon, stat=ierr)
1740 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lon')
1741 :
1742 403200 : deallocate(state%ps, stat=ierr)
1743 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%ps')
1744 :
1745 403200 : deallocate(state%psdry, stat=ierr)
1746 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%psdry')
1747 :
1748 403200 : deallocate(state%phis, stat=ierr)
1749 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%phis')
1750 :
1751 403200 : deallocate(state%ulat, stat=ierr)
1752 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%ulat')
1753 :
1754 403200 : deallocate(state%ulon, stat=ierr)
1755 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%ulon')
1756 :
1757 403200 : deallocate(state%t, stat=ierr)
1758 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%t')
1759 :
1760 403200 : deallocate(state%u, stat=ierr)
1761 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%u')
1762 :
1763 403200 : deallocate(state%v, stat=ierr)
1764 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%v')
1765 :
1766 403200 : deallocate(state%s, stat=ierr)
1767 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%s')
1768 :
1769 403200 : deallocate(state%omega, stat=ierr)
1770 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%omega')
1771 :
1772 403200 : deallocate(state%pmid, stat=ierr)
1773 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%pmid')
1774 :
1775 403200 : deallocate(state%pmiddry, stat=ierr)
1776 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%pmiddry')
1777 :
1778 403200 : deallocate(state%pdel, stat=ierr)
1779 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%pdel')
1780 :
1781 403200 : deallocate(state%pdeldry, stat=ierr)
1782 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%pdeldry')
1783 :
1784 403200 : deallocate(state%rpdel, stat=ierr)
1785 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%rpdel')
1786 :
1787 403200 : deallocate(state%rpdeldry, stat=ierr)
1788 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%rpdeldry')
1789 :
1790 403200 : deallocate(state%lnpmid, stat=ierr)
1791 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lnpmid')
1792 :
1793 403200 : deallocate(state%lnpmiddry, stat=ierr)
1794 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lnpmiddry')
1795 :
1796 403200 : deallocate(state%exner, stat=ierr)
1797 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%exner')
1798 :
1799 403200 : deallocate(state%zm, stat=ierr)
1800 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%zm')
1801 :
1802 403200 : deallocate(state%q, stat=ierr)
1803 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%q')
1804 :
1805 403200 : deallocate(state%pint, stat=ierr)
1806 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%pint')
1807 :
1808 403200 : deallocate(state%pintdry, stat=ierr)
1809 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%pintdry')
1810 :
1811 403200 : deallocate(state%lnpint, stat=ierr)
1812 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lnpint')
1813 :
1814 403200 : deallocate(state%lnpintdry, stat=ierr)
1815 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lnpintdry')
1816 :
1817 403200 : deallocate(state%zi, stat=ierr)
1818 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%zi')
1819 :
1820 403200 : deallocate(state%te_ini, stat=ierr)
1821 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%te_ini')
1822 :
1823 403200 : deallocate(state%te_cur, stat=ierr)
1824 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%te_cur')
1825 :
1826 403200 : deallocate(state%tw_ini, stat=ierr)
1827 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%tw_ini')
1828 :
1829 403200 : deallocate(state%tw_cur, stat=ierr)
1830 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%tw_cur')
1831 :
1832 403200 : deallocate(state%temp_ini, stat=ierr)
1833 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%temp_ini')
1834 :
1835 403200 : deallocate(state%z_ini, stat=ierr)
1836 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%z_ini')
1837 :
1838 403200 : deallocate(state%latmapback, stat=ierr)
1839 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%latmapback')
1840 :
1841 403200 : deallocate(state%lonmapback, stat=ierr)
1842 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%lonmapback')
1843 :
1844 403200 : deallocate(state%cid, stat=ierr)
1845 403200 : if ( ierr /= 0 ) call endrun('physics_state_dealloc error: deallocation error for state%cid')
1846 :
1847 403200 : end subroutine physics_state_dealloc
1848 :
1849 : !===============================================================================
1850 :
1851 7680 : subroutine physics_tend_alloc(tend,psetcols)
1852 :
1853 : use infnan, only : inf, assignment(=)
1854 : ! allocate the individual tend components
1855 :
1856 : type(physics_tend), intent(inout) :: tend
1857 :
1858 : integer, intent(in) :: psetcols
1859 :
1860 : integer :: ierr = 0
1861 :
1862 7680 : tend%psetcols = psetcols
1863 :
1864 23040 : allocate(tend%dtdt(psetcols,pver), stat=ierr)
1865 7680 : if ( ierr /= 0 ) call endrun('physics_tend_alloc error: allocation error for tend%dtdt')
1866 :
1867 15360 : allocate(tend%dudt(psetcols,pver), stat=ierr)
1868 7680 : if ( ierr /= 0 ) call endrun('physics_tend_alloc error: allocation error for tend%dudt')
1869 :
1870 15360 : allocate(tend%dvdt(psetcols,pver), stat=ierr)
1871 7680 : if ( ierr /= 0 ) call endrun('physics_tend_alloc error: allocation error for tend%dvdt')
1872 :
1873 23040 : allocate(tend%flx_net(psetcols), stat=ierr)
1874 7680 : if ( ierr /= 0 ) call endrun('physics_tend_alloc error: allocation error for tend%flx_net')
1875 :
1876 15360 : allocate(tend%te_tnd(psetcols), stat=ierr)
1877 7680 : if ( ierr /= 0 ) call endrun('physics_tend_alloc error: allocation error for tend%te_tnd')
1878 :
1879 15360 : allocate(tend%tw_tnd(psetcols), stat=ierr)
1880 7680 : if ( ierr /= 0 ) call endrun('physics_tend_alloc error: allocation error for tend%tw_tnd')
1881 :
1882 7680 : tend%dtdt(:,:) = inf
1883 7680 : tend%dudt(:,:) = inf
1884 7680 : tend%dvdt(:,:) = inf
1885 7680 : tend%flx_net(:) = inf
1886 7680 : tend%te_tnd(:) = inf
1887 7680 : tend%tw_tnd(:) = inf
1888 :
1889 7680 : end subroutine physics_tend_alloc
1890 :
1891 : !===============================================================================
1892 :
1893 0 : subroutine physics_tend_dealloc(tend)
1894 :
1895 : ! deallocate the individual tend components
1896 :
1897 : type(physics_tend), intent(inout) :: tend
1898 : integer :: ierr = 0
1899 :
1900 0 : deallocate(tend%dtdt, stat=ierr)
1901 0 : if ( ierr /= 0 ) call endrun('physics_tend_dealloc error: deallocation error for tend%dtdt')
1902 :
1903 0 : deallocate(tend%dudt, stat=ierr)
1904 0 : if ( ierr /= 0 ) call endrun('physics_tend_dealloc error: deallocation error for tend%dudt')
1905 :
1906 0 : deallocate(tend%dvdt, stat=ierr)
1907 0 : if ( ierr /= 0 ) call endrun('physics_tend_dealloc error: deallocation error for tend%dvdt')
1908 :
1909 0 : deallocate(tend%flx_net, stat=ierr)
1910 0 : if ( ierr /= 0 ) call endrun('physics_tend_dealloc error: deallocation error for tend%flx_net')
1911 :
1912 0 : deallocate(tend%te_tnd, stat=ierr)
1913 0 : if ( ierr /= 0 ) call endrun('physics_tend_dealloc error: deallocation error for tend%te_tnd')
1914 :
1915 0 : deallocate(tend%tw_tnd, stat=ierr)
1916 0 : if ( ierr /= 0 ) call endrun('physics_tend_dealloc error: deallocation error for tend%tw_tnd')
1917 0 : end subroutine physics_tend_dealloc
1918 :
1919 : !===============================================================================
1920 :
1921 3083520 : subroutine physics_ptend_alloc(ptend,psetcols)
1922 :
1923 : ! allocate the individual ptend components
1924 :
1925 : type(physics_ptend), intent(inout) :: ptend
1926 :
1927 : integer, intent(in) :: psetcols
1928 :
1929 : integer :: ierr = 0
1930 :
1931 3083520 : ptend%psetcols = psetcols
1932 :
1933 3083520 : if (ptend%ls) then
1934 5713920 : allocate(ptend%s(psetcols,pver), stat=ierr)
1935 1904640 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%s')
1936 :
1937 5713920 : allocate(ptend%hflux_srf(psetcols), stat=ierr)
1938 1904640 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%hflux_srf')
1939 :
1940 3809280 : allocate(ptend%hflux_top(psetcols), stat=ierr)
1941 1904640 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%hflux_top')
1942 : end if
1943 :
1944 3083520 : if (ptend%lu) then
1945 2304000 : allocate(ptend%u(psetcols,pver), stat=ierr)
1946 768000 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%u')
1947 :
1948 2304000 : allocate(ptend%taux_srf(psetcols), stat=ierr)
1949 768000 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%taux_srf')
1950 :
1951 1536000 : allocate(ptend%taux_top(psetcols), stat=ierr)
1952 768000 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%taux_top')
1953 : end if
1954 :
1955 3083520 : if (ptend%lv) then
1956 2085120 : allocate(ptend%v(psetcols,pver), stat=ierr)
1957 695040 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%v')
1958 :
1959 2085120 : allocate(ptend%tauy_srf(psetcols), stat=ierr)
1960 695040 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%tauy_srf')
1961 :
1962 1390080 : allocate(ptend%tauy_top(psetcols), stat=ierr)
1963 695040 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%tauy_top')
1964 : end if
1965 :
1966 152605440 : if (any(ptend%lq)) then
1967 10168320 : allocate(ptend%q(psetcols,pver,pcnst), stat=ierr)
1968 2542080 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%q')
1969 :
1970 7626240 : allocate(ptend%cflx_srf(psetcols,pcnst), stat=ierr)
1971 2542080 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%cflx_srf')
1972 :
1973 5084160 : allocate(ptend%cflx_top(psetcols,pcnst), stat=ierr)
1974 2542080 : if ( ierr /= 0 ) call endrun('physics_ptend_alloc error: allocation error for ptend%cflx_top')
1975 : end if
1976 :
1977 3083520 : end subroutine physics_ptend_alloc
1978 :
1979 : !===============================================================================
1980 :
1981 3728640 : subroutine physics_ptend_dealloc(ptend)
1982 :
1983 : ! deallocate the individual ptend components
1984 :
1985 : type(physics_ptend), intent(inout) :: ptend
1986 : integer :: ierr = 0
1987 :
1988 3728640 : ptend%psetcols = 0
1989 :
1990 3728640 : if (allocated(ptend%s)) deallocate(ptend%s, stat=ierr)
1991 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%s')
1992 :
1993 3728640 : if (allocated(ptend%hflux_srf)) deallocate(ptend%hflux_srf, stat=ierr)
1994 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%hflux_srf')
1995 :
1996 3728640 : if (allocated(ptend%hflux_top)) deallocate(ptend%hflux_top, stat=ierr)
1997 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%hflux_top')
1998 :
1999 3728640 : if (allocated(ptend%u)) deallocate(ptend%u, stat=ierr)
2000 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%u')
2001 :
2002 3728640 : if (allocated(ptend%taux_srf)) deallocate(ptend%taux_srf, stat=ierr)
2003 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%taux_srf')
2004 :
2005 3728640 : if (allocated(ptend%taux_top)) deallocate(ptend%taux_top, stat=ierr)
2006 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%taux_top')
2007 :
2008 3728640 : if (allocated(ptend%v)) deallocate(ptend%v, stat=ierr)
2009 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%v')
2010 :
2011 3728640 : if (allocated(ptend%tauy_srf)) deallocate(ptend%tauy_srf, stat=ierr)
2012 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%tauy_srf')
2013 :
2014 3728640 : if (allocated(ptend%tauy_top)) deallocate(ptend%tauy_top, stat=ierr)
2015 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%tauy_top')
2016 :
2017 3728640 : if (allocated(ptend%q)) deallocate(ptend%q, stat=ierr)
2018 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%q')
2019 :
2020 3728640 : if (allocated(ptend%cflx_srf)) deallocate(ptend%cflx_srf, stat=ierr)
2021 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%cflx_srf')
2022 :
2023 3728640 : if(allocated(ptend%cflx_top)) deallocate(ptend%cflx_top, stat=ierr)
2024 3728640 : if ( ierr /= 0 ) call endrun('physics_ptend_dealloc error: deallocation error for ptend%cflx_top')
2025 :
2026 3728640 : end subroutine physics_ptend_dealloc
2027 :
2028 0 : end module physics_types
|