CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_debug.F90
1
2
3
4!!
5!! This is the place to switch between different debug outputs.
6!! - The default behavior for Intel (or any compiler other than GNU)
7!! is to print mininmum, maximum and 32-bit Adler checksum for arrays.
8!! - The default behavior for GNU is to mininmum, maximum and
9!! mean value of arrays, because calculating the checksum leads
10!! to segmentation faults with gfortran (bug in malloc?).
11!! - If none of the #define preprocessor statements is used,
12!! arrays are printed in full (this is often unpractical).
13!! - All output to stdout/stderr from these routines are prefixed
14!! with 'XXX: ' so that they can be easily removed from the log files
15!! using "grep -ve 'XXX: ' ..." if needed.
16!! - Only one #define statement can be active at any time
17!!
18!! Available options for debug output:
19!!
20!! #define PRINT_SUM: print mininmum, maximum and mean value of arrays
21!!
22!! #define PRINT_CHKSUM: mininmum, maximum and 32-bit Adler checksum for arrays
23!!
24
25#ifdef __GFORTRAN__
26#define PRINT_SUM
27#else
28#define PRINT_CHKSUM
29#endif
30
31!!
32!!
33!!
34
35 module print_var_chksum
36
37 use machine, only: kind_phys
38
39 implicit none
40
41 private
42
43 public chksum_int, chksum_real, print_var
44
45 interface print_var
46 module procedure print_logic_0d
47 module procedure print_logic_1d
48 module procedure print_int_0d
49 module procedure print_int_1d
50 module procedure print_int_2d
51 module procedure print_real_0d
52 module procedure print_real_1d
53 module procedure print_real_2d
54 module procedure print_real_3d
55 module procedure print_real_4d
56 end interface
57
58 contains
59
60 subroutine print_logic_0d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
61
62 integer, intent(in) :: mpirank, omprank, blkno
63 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
64 character(len=*), intent(in) :: name
65 logical, intent(in) :: var
66
67 write(0,'(2a,3i6,1x,l)') 'XXX: ', trim(name), mpirank, omprank, blkno, var
68
69 end subroutine print_logic_0d
70
71 subroutine print_logic_1d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
72
73 integer, intent(in) :: mpirank, omprank, blkno
74 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
75 character(len=*), intent(in) :: name
76 logical, intent(in) :: var(:)
77
78 integer :: i
79
80#ifdef PRINT_SUM
81 write(0,'(2a,3i6,2i8)') 'XXX: ', trim(name), mpirank, omprank, blkno, size(var), count(var)
82#elif defined(PRINT_CHKSUM)
83 write(0,'(2a,3i6,2i8)') 'XXX: ', trim(name), mpirank, omprank, blkno, size(var), count(var)
84#else
85 do i=lbound(var,1),ubound(var,1)
86 write(0,'(2a,3i6,i6,2e16.7,1x,l)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, lat_d(i), lon_d(i), var(i)
87 end do
88#endif
89
90 end subroutine print_logic_1d
91
92 subroutine print_int_0d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
93
94 integer, intent(in) :: mpirank, omprank, blkno
95 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
96 character(len=*), intent(in) :: name
97 integer, intent(in) :: var
98
99 write(0,'(2a,3i6,i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, var
100
101 end subroutine print_int_0d
102
103 subroutine print_int_1d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
104
105 integer, intent(in) :: mpirank, omprank, blkno
106 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
107 character(len=*), intent(in) :: name
108 integer, intent(in) :: var(:)
109
110 integer :: i
111
112#ifdef PRINT_SUM
113 write(0,'(2a,3i6,3i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, sum(var), minval(var), maxval(var)
114#elif defined(PRINT_CHKSUM)
115 write(0,'(2a,3i6,i20,2i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, chksum_int(size(var),var), minval(var), maxval(var)
116#else
117 do i=lbound(var,1),ubound(var,1)
118 write(0,'(2a,3i6,i6,2e16.7,i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, lat_d(i), lon_d(i), var(i)
119 end do
120#endif
121
122 end subroutine print_int_1d
123
124 subroutine print_int_2d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
125
126 integer, intent(in) :: mpirank, omprank, blkno
127 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
128 character(len=*), intent(in) :: name
129 integer, intent(in) :: var(:,:)
130
131 integer :: i, k
132
133#ifdef PRINT_SUM
134 write(0,'(2a,3i6,3i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, sum(var), minval(var), maxval(var)
135#elif defined(PRINT_CHKSUM)
136 write(0,'(2a,3i6,i20,2i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, chksum_int(size(var),var), minval(var), maxval(var)
137#else
138 do i=lbound(var,1),ubound(var,1)
139 do k=lbound(var,2),ubound(var,2)
140 write(0,'(2a,3i6,2i6,2e16.7,i15)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, k, lat_d(i), lon_d(i), var(i,k)
141 end do
142 end do
143#endif
144
145 end subroutine print_int_2d
146
147 subroutine print_real_0d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
148
149 integer, intent(in) :: mpirank, omprank, blkno
150 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
151 character(len=*), intent(in) :: name
152 real(kind_phys), intent(in) :: var
153
154 write(0,'(2a,3i6,e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, var
155
156 end subroutine print_real_0d
157
158 subroutine print_real_1d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
159
160 integer, intent(in) :: mpirank, omprank, blkno
161 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
162 character(len=*), intent(in) :: name
163 real(kind_phys), intent(in) :: var(:)
164
165 integer :: i
166
167#ifdef PRINT_SUM
168 write(0,'(2a,3i6,3e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, sum(var), minval(var), maxval(var)
169#elif defined(PRINT_CHKSUM)
170 write(0,'(2a,3i6,i20,2e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, chksum_real(size(var),var), minval(var), maxval(var)
171#else
172 do i=lbound(var,1),ubound(var,1)
173 write(0,'(2a,3i6,i6,2e16.7,e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, lat_d(i), lon_d(i), var(i)
174 end do
175#endif
176
177 end subroutine print_real_1d
178
179 subroutine print_real_2d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
180
181 integer, intent(in) :: mpirank, omprank, blkno
182 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
183 character(len=*), intent(in) :: name
184 real(kind_phys), intent(in) :: var(:,:)
185
186 integer :: k, i
187
188#ifdef PRINT_SUM
189 write(0,'(2a,3i6,3e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, sum(var), minval(var), maxval(var)
190#elif defined(PRINT_CHKSUM)
191 write(0,'(2a,3i6,i20,2e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, chksum_real(size(var),reshape(var,(/size(var)/))), minval(var), maxval(var)
192#else
193 do i=lbound(var,1),ubound(var,1)
194 do k=lbound(var,2),ubound(var,2)
195 write(0,'(2a,3i6,2i6,2e16.7,e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, k, lat_d(i), lon_d(i), var(i,k)
196 end do
197 end do
198#endif
199
200 end subroutine print_real_2d
201
202 subroutine print_real_3d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
203
204 integer, intent(in) :: mpirank, omprank, blkno
205 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
206 character(len=*), intent(in) :: name
207 real(kind_phys), intent(in) :: var(:,:,:)
208
209 integer :: k, i, l
210
211#ifdef PRINT_SUM
212 write(0,'(2a,3i6,3e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, sum(var), minval(var), maxval(var)
213#elif defined(PRINT_CHKSUM)
214 write(0,'(2a,3i6,i20,2e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, chksum_real(size(var),reshape(var,(/size(var)/))), minval(var), maxval(var)
215#else
216 do i=lbound(var,1),ubound(var,1)
217 do k=lbound(var,2),ubound(var,2)
218 do l=lbound(var,3),ubound(var,3)
219 write(0,'(2a,3i6,3i6,2e16.7,e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, k, l, lat_d(i), lon_d(i), var(i,k,l)
220 end do
221 end do
222 end do
223#endif
224
225 end subroutine print_real_3d
226
227 subroutine print_real_4d(mpirank, omprank, blkno, lat_d, lon_d, name, var)
228
229 integer, intent(in) :: mpirank, omprank, blkno
230 real(kind_phys), intent(in) :: lat_d(:), lon_d(:)
231 character(len=*), intent(in) :: name
232 real(kind_phys), intent(in) :: var(:,:,:,:)
233
234 integer :: k, i, l, m
235
236#ifdef PRINT_SUM
237 write(0,'(2a,3i6,3e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, sum(var), minval(var), maxval(var)
238#elif defined(PRINT_CHKSUM)
239 write(0,'(2a,3i6,i20,2e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, chksum_real(size(var),reshape(var,(/size(var)/))), minval(var), maxval(var)
240#else
241 do i=lbound(var,1),ubound(var,1)
242 do k=lbound(var,2),ubound(var,2)
243 do l=lbound(var,3),ubound(var,3)
244 do m=lbound(var,4),ubound(var,4)
245 write(0,'(2a,3i6,4i6,2e16.7,e35.25)') 'XXX: ', trim(name), mpirank, omprank, blkno, i, k, l, m, lat_d(i), lon_d(i), var(i,k,l,m)
246 end do
247 end do
248 end do
249 end do
250#endif
251
252 end subroutine print_real_4d
253
254 function chksum_int(N, var) result(hash)
255
256 integer, intent(in) :: N
257 integer, dimension(1:N), intent(in) :: var
258 integer*8, dimension(1:N) :: int_var
259 integer*8 :: a, b, i, hash
260 integer*8, parameter :: mod_adler=65521
261
262 a=1
263 b=0
264 i=1
265 hash = 0
266 int_var = transfer(var, a, n)
267
268 do i= 1, n
269 a = mod(a + int_var(i), mod_adler)
270 b = mod(b+a, mod_adler)
271 end do
272
273 hash = ior(b * 65536, a)
274
275 end function chksum_int
276
277 function chksum_real(N, var) result(hash)
278
279 integer, intent(in) :: N
280 real(kind_phys), dimension(1:N), intent(in) :: var
281 integer*8, dimension(1:N) :: int_var
282 integer*8 :: a, b, i, hash
283 integer*8, parameter :: mod_adler=65521
284
285 a=1
286 b=0
287 i=1
288 hash = 0
289 int_var = transfer(var, a, n)
290
291 do i= 1, n
292 a = mod(a + int_var(i), mod_adler)
293 b = mod(b+a, mod_adler)
294 end do
295
296 hash = ior(b * 65536, a)
297
298 end function chksum_real
299
300 end module print_var_chksum
301
303
304 use print_var_chksum, only: print_var
305
306 implicit none
307
308 private
309
310 public gfs_diagtoscreen_init, gfs_diagtoscreen_timestep_init, gfs_diagtoscreen_run
311
312 contains
313
317 subroutine gfs_diagtoscreen_init (Model, Data, Interstitial, errmsg, errflg)
318
319 use gfs_typedefs, only: gfs_control_type, gfs_data_type
320 use ccpp_typedefs, only: gfs_interstitial_type
321
322 implicit none
323
324 !--- interface variables
325 type(gfs_control_type), intent(in) :: model
326 type(gfs_data_type), intent(in) :: data(:)
327 type(gfs_interstitial_type), intent(in) :: interstitial(:)
328 character(len=*), intent(out) :: errmsg
329 integer, intent(out) :: errflg
330
331 !--- local variables
332 integer :: i
333
334 ! Initialize CCPP error handling variables
335 errmsg = ''
336 errflg = 0
337
338 do i=1,size(data)
339 call gfs_diagtoscreen_run (model, Data(i)%Statein, Data(i)%Stateout, Data(i)%Sfcprop, &
340 Data(i)%Coupling, Data(i)%Grid, Data(i)%Tbd, Data(i)%Cldprop, &
341 Data(i)%Radtend, Data(i)%Intdiag, interstitial(1), &
342 size(interstitial), i, errmsg, errflg)
343 end do
344
345 end subroutine gfs_diagtoscreen_init
346
350 subroutine gfs_diagtoscreen_timestep_init (Model, Data, Interstitial, errmsg, errflg)
351
352 use gfs_typedefs, only: gfs_control_type, gfs_data_type
353 use ccpp_typedefs, only: gfs_interstitial_type
354
355 implicit none
356
357 !--- interface variables
358 type(gfs_control_type), intent(in) :: model
359 type(gfs_data_type), intent(in) :: data(:)
360 type(gfs_interstitial_type), intent(in) :: interstitial(:)
361 character(len=*), intent(out) :: errmsg
362 integer, intent(out) :: errflg
363
364 !--- local variables
365 integer :: i
366
367 ! Initialize CCPP error handling variables
368 errmsg = ''
369 errflg = 0
370
371 do i=1,size(data)
372 call gfs_diagtoscreen_run (model, Data(i)%Statein, Data(i)%Stateout, Data(i)%Sfcprop, &
373 Data(i)%Coupling, Data(i)%Grid, Data(i)%Tbd, Data(i)%Cldprop, &
374 Data(i)%Radtend, Data(i)%Intdiag, interstitial(1), &
375 size(interstitial), i, errmsg, errflg)
376 end do
377
378 end subroutine gfs_diagtoscreen_timestep_init
379
383 subroutine gfs_diagtoscreen_run (Model, Statein, Stateout, Sfcprop, Coupling, &
384 Grid, Tbd, Cldprop, Radtend, Diag, Interstitial, &
385 nthreads, blkno, errmsg, errflg)
386
387#ifdef MPI
388 use mpi_f08
389#endif
390#ifdef _OPENMP
391 use omp_lib
392#endif
393 use gfs_typedefs, only: gfs_control_type, gfs_statein_type, &
394 gfs_stateout_type, gfs_sfcprop_type, &
395 gfs_coupling_type, gfs_grid_type, &
396 gfs_tbd_type, gfs_cldprop_type, &
397 gfs_radtend_type, gfs_diag_type
398 use ccpp_typedefs, only: gfs_interstitial_type
399
400 implicit none
401
402 !--- interface variables
403 type(gfs_control_type), intent(in ) :: model
404 type(gfs_statein_type), intent(in ) :: statein
405 type(gfs_stateout_type), intent(in ) :: stateout
406 type(gfs_sfcprop_type), intent(in ) :: sfcprop
407 type(gfs_coupling_type), intent(in ) :: coupling
408 type(gfs_grid_type), intent(in ) :: grid
409 type(gfs_tbd_type), intent(in ) :: tbd
410 type(gfs_cldprop_type), intent(in ) :: cldprop
411 type(gfs_radtend_type), intent(in ) :: radtend
412 type(gfs_diag_type), intent(in ) :: diag
413 type(gfs_interstitial_type), intent(in) :: interstitial
414 integer, intent(in ) :: nthreads
415 integer, intent(in ) :: blkno
416 character(len=*), intent(out) :: errmsg
417 integer, intent(out) :: errflg
418
419 !--- local variables
420 integer :: impi, iomp, ierr, n, idtend, iprocess, itracer
421 integer :: mpirank, mpisize, mpicomm
422 integer :: omprank, ompsize
423
424 ! Initialize CCPP error handling variables
425 errmsg = ''
426 errflg = 0
427
428#ifdef MPI
429 mpicomm = model%communicator
430 mpirank = model%me
431 mpisize = model%ntasks
432#else
433 mpirank = 0
434 mpisize = 1
435 mpicomm = 0
436#endif
437#ifdef _OPENMP
438 omprank = omp_get_thread_num()
439 ompsize = nthreads
440#else
441 omprank = 0
442 ompsize = 1
443#endif
444
445#ifdef _OPENMP
446!$OMP BARRIER
447#endif
448#ifdef MPI
449! call MPI_BARRIER(mpicomm,ierr)
450#endif
451
452 do impi=0,mpisize-1
453 do iomp=0,ompsize-1
454 if (mpirank==impi .and. omprank==iomp) then
455 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Model%kdt' , model%kdt)
456 ! Sfcprop
457 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%slmsk' , sfcprop%slmsk)
458 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%oceanfrac' , sfcprop%oceanfrac)
459 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%landfrac' , sfcprop%landfrac)
460 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%lakefrac' , sfcprop%lakefrac)
461 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tsfc' , sfcprop%tsfc)
462 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tsfco' , sfcprop%tsfco)
463 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tsfcl' , sfcprop%tsfcl)
464 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tisfc' , sfcprop%tisfc)
465 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%snowd' , sfcprop%snowd)
466 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%zorl' , sfcprop%zorl)
467 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%zorlw' , sfcprop%zorlw)
468 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%zorll' , sfcprop%zorll)
469 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%zorli' , sfcprop%zorli)
470 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%zorlwav' , sfcprop%zorlwav)
471 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%fice' , sfcprop%fice)
472 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%hprime' , sfcprop%hprime)
473 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%sncovr' , sfcprop%sncovr)
474 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%snoalb' , sfcprop%snoalb)
475 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%alvsf' , sfcprop%alvsf)
476 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%alnsf' , sfcprop%alnsf)
477 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%alvwf' , sfcprop%alvwf)
478 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%alnwf' , sfcprop%alnwf)
479 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%facsf' , sfcprop%facsf)
480 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%facwf' , sfcprop%facwf)
481 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%slope' , sfcprop%slope)
482 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%slope_save', sfcprop%slope_save)
483 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%shdmin' , sfcprop%shdmin)
484 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%shdmax' , sfcprop%shdmax)
485 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tg3' , sfcprop%tg3)
486 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%vfrac' , sfcprop%vfrac)
487 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%vtype' , sfcprop%vtype)
488 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%vtype_save', sfcprop%vtype_save)
489 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%stype' , sfcprop%stype)
490 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%stype_save', sfcprop%stype_save)
491
492 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%scolor' , sfcprop%scolor)
493 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%scolore_save', sfcprop%scolor_save)
494
495 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%uustar' , sfcprop%uustar)
496 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%oro' , sfcprop%oro)
497 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%oro_uf' , sfcprop%oro_uf)
498 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%hice' , sfcprop%hice)
499 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%weasd' , sfcprop%weasd)
500 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%weasdl' , sfcprop%weasdl)
501 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%%weasdi' , sfcprop%weasdi)
502 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%canopy' , sfcprop%canopy)
503 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%ffmm' , sfcprop%ffmm)
504 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%ffhh' , sfcprop%ffhh)
505 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%f10m' , sfcprop%f10m)
506 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tprcp' , sfcprop%tprcp)
507 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%srflag' , sfcprop%srflag)
508 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%slc' , sfcprop%slc)
509 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%smc' , sfcprop%smc)
510 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%stc' , sfcprop%stc)
511 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%t2m' , sfcprop%t2m)
512 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%q2m' , sfcprop%q2m)
513 if (model%nstf_name(1)>0) then
514 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tref ', sfcprop%tref)
515 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%z_c ', sfcprop%z_c)
516 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%c_0 ', sfcprop%c_0)
517 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%c_d ', sfcprop%c_d)
518 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%w_0 ', sfcprop%w_0)
519 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%w_d ', sfcprop%w_d)
520 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xt ', sfcprop%xt)
521 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xs ', sfcprop%xs)
522 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xu ', sfcprop%xu)
523 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xv ', sfcprop%xv)
524 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xz ', sfcprop%xz)
525 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%zm ', sfcprop%zm)
526 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xtts ', sfcprop%xtts)
527 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%xzts ', sfcprop%xzts)
528 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%d_conv ', sfcprop%d_conv)
529 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%ifd ', sfcprop%ifd)
530 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%dt_cool ', sfcprop%dt_cool)
531 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%qrain ', sfcprop%qrain)
532 end if
533 ! CCPP/RUC only
534 if (model%lsm == model%lsm_ruc) then
535 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%sh2o', sfcprop%sh2o)
536 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%smois', sfcprop%smois)
537 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tslb', sfcprop%tslb)
538 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%clw_surf_land', sfcprop%clw_surf_land)
539 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%clw_surf_ice', sfcprop%clw_surf_ice)
540 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%qwv_surf_land', sfcprop%qwv_surf_land)
541 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%qwv_surf_ice', sfcprop%qwv_surf_ice)
542 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%flag_frsoil', sfcprop%flag_frsoil)
543 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%rhofr', sfcprop%rhofr)
544 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tsnow_land', sfcprop%tsnow_land)
545 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%tsnow_ice', sfcprop%tsnow_ice)
546 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%snowfallac_land', sfcprop%snowfallac_land)
547 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%snowfallac_ice', sfcprop%snowfallac_ice)
548 end if
549 ! Revised surface albedo and emissivity calculation
550 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%emis_lnd', sfcprop%emis_lnd)
551 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%emis_ice', sfcprop%emis_ice)
552 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%emis_wat', sfcprop%emis_wat)
553 ! NoahMP and RUC
554 if (model%lsm == model%lsm_ruc .or. model%lsm == model%lsm_noahmp) then
555 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdirvis_lnd', sfcprop%albdirvis_lnd)
556 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdirnir_lnd', sfcprop%albdirnir_lnd)
557 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdifvis_lnd', sfcprop%albdifvis_lnd)
558 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdifnir_lnd', sfcprop%albdifnir_lnd)
559 end if
560 ! RUC only
561 if (model%lsm == model%lsm_ruc) then
562 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdirvis_ice', sfcprop%albdirvis_ice)
563 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdifvis_ice', sfcprop%albdifvis_ice)
564 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdirnir_ice', sfcprop%albdirnir_ice)
565 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%albdifnir_ice', sfcprop%albdifnir_ice)
566 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%sfalb_lnd', sfcprop%sfalb_lnd)
567 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%sfalb_ice', sfcprop%sfalb_ice)
568 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%sfalb_lnd_bck', sfcprop%sfalb_lnd_bck)
569 end if
570 ! Radtend
571 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcfsw%upfxc', radtend%sfcfsw(:)%upfxc)
572 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcfsw%dnfxc', radtend%sfcfsw(:)%dnfxc)
573 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcfsw%upfx0', radtend%sfcfsw(:)%upfx0)
574 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcfsw%dnfx0', radtend%sfcfsw(:)%dnfx0)
575 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcflw%upfxc', radtend%sfcflw(:)%upfxc)
576 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcflw%upfx0', radtend%sfcflw(:)%upfx0)
577 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcflw%dnfxc', radtend%sfcflw(:)%dnfxc)
578 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfcflw%dnfx0', radtend%sfcflw(:)%dnfx0)
579 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%htrsw', radtend%htrsw)
580 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%htrlw', radtend%htrlw)
581 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%sfalb', radtend%sfalb)
582 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%coszen', radtend%coszen)
583 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%tsflw', radtend%tsflw)
584 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%semis', radtend%semis)
585 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%coszdg', radtend%coszdg)
586 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%swhc', radtend%swhc)
587 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%lwhc', radtend%lwhc)
588 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Radtend%lwhd', radtend%lwhd)
589 ! Tbd
590 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%icsdsw' , tbd%icsdsw)
591 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%icsdlw' , tbd%icsdlw)
592 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%ozpl' , tbd%ozpl)
593 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%h2opl' , tbd%h2opl)
594 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%rann' , tbd%rann)
595 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%acv' , tbd%acv)
596 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%acvb' , tbd%acvb)
597 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%acvt' , tbd%acvt)
598 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%hpbl' , tbd%hpbl)
599 if(model%imfdeepcnv>0 .or. model%imfshalcnv>0) then
600 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%ud_mf' , tbd%ud_mf)
601 endif
602 if (model%do_sppt) then
603 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%dtdtnp' , tbd%dtdtnp)
604 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%dtotprcp' , tbd%dtotprcp)
605 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%dcnvprcp' , tbd%dcnvprcp)
606 end if
607 if (model%cplflx .or. model%cplchm) then
608 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%drain_cpl' , tbd%drain_cpl)
609 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%dsnow_cpl' , tbd%dsnow_cpl)
610 end if
611 if (model%nctp > 0 .and. model%cscnv) then
612 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%phy_fctd' , tbd%phy_fctd)
613 end if
614 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%phy_f2d' , tbd%phy_f2d)
615 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%phy_f3d' , tbd%phy_f3d)
616 do n=1,size(tbd%phy_f3d(1,1,:))
617 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%phy_f3d_n' , tbd%phy_f3d(:,:,n))
618 end do
619 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%in_nm' , tbd%in_nm)
620 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%ccn_nm' , tbd%ccn_nm)
621 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%aer_nm' , tbd%aer_nm)
622 if (model%imfdeepcnv == model%imfdeepcnv_gf .or. model%imfdeepcnv == model%imfdeepcnv_unified) then
623 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%cactiv' , tbd%cactiv)
624 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%cactiv_m' , tbd%cactiv_m)
625 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Tbd%aod_gf' , tbd%aod_gf)
626 end if
627 ! Diag
628 !call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Diag%fluxr ', Diag%fluxr)
629 !do n=1,size(Diag%fluxr(1,:))
630 ! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Diag%fluxr_n ', Diag%fluxr(:,n))
631 !end do
632 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%srunoff ', diag%srunoff)
633 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%evbs ', diag%evbs)
634 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%evcw ', diag%evcw)
635 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%sbsno ', diag%sbsno)
636 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%evbsa ', diag%evbsa)
637 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%evcwa ', diag%evcwa)
638 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%snohfa ', diag%snohfa)
639 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%transa ', diag%transa)
640 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%sbsnoa ', diag%sbsnoa)
641 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%snowca ', diag%snowca)
642 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%soilm ', diag%soilm)
643 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%tmpmin ', diag%tmpmin)
644 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%tmpmax ', diag%tmpmax)
645 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dusfc ', diag%dusfc)
646 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvsfc ', diag%dvsfc)
647 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dtsfc ', diag%dtsfc)
648 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dqsfc ', diag%dqsfc)
649 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totprcp ', diag%totprcp)
650 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totice ', diag%totice)
651 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totsnw ', diag%totsnw)
652 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totgrp ', diag%totgrp)
653 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totprcpb ', diag%totprcpb)
654 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%toticeb ', diag%toticeb)
655 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totsnwb ', diag%totsnwb)
656 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%totgrpb ', diag%totgrpb)
657 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%suntim ', diag%suntim)
658 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%runoff ', diag%runoff)
659 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%ep ', diag%ep)
660 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%cldwrk ', diag%cldwrk)
661 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dugwd ', diag%dugwd)
662 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvgwd ', diag%dvgwd)
663 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%psmean ', diag%psmean)
664 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%cnvprcp ', diag%cnvprcp)
665 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%cnvprcpb ', diag%cnvprcpb)
666 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%spfhmin ', diag%spfhmin)
667 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%spfhmax ', diag%spfhmax)
668 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%u10mmax ', diag%u10mmax)
669 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%v10mmax ', diag%v10mmax)
670 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%wind10mmax ', diag%wind10mmax)
671 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%rain ', diag%rain)
672 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%rainc ', diag%rainc)
673 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%ice ', diag%ice)
674 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%snow ', diag%snow)
675 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%graupel ', diag%graupel)
676 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%u10m ', diag%u10m)
677 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%v10m ', diag%v10m)
678 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dpt2m ', diag%dpt2m)
679 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%zlvl ', diag%zlvl)
680 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%psurf ', diag%psurf)
681 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%pwat ', diag%pwat)
682 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%t1 ', diag%t1)
683 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%q1 ', diag%q1)
684 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%u1 ', diag%u1)
685 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%v1 ', diag%v1)
686 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%chh ', diag%chh)
687 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%cmm ', diag%cmm)
688 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dlwsfci ', diag%dlwsfci)
689 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%ulwsfci ', diag%ulwsfci)
690 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dswsfci ', diag%dswsfci)
691 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%nswsfci ', diag%nswsfci)
692 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%uswsfci ', diag%uswsfci)
693 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dusfci ', diag%dusfci)
694 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvsfci ', diag%dvsfci)
695 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dtsfci ', diag%dtsfci)
696 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dqsfci ', diag%dqsfci)
697 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%gfluxi ', diag%gfluxi)
698 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%epi ', diag%epi)
699 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%smcwlt2 ', diag%smcwlt2)
700 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%smcref2 ', diag%smcref2)
701 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%sr ', diag%sr)
702 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%tdomr ', diag%tdomr)
703 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%tdomzr ', diag%tdomzr)
704 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%tdomip ', diag%tdomip)
705 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%tdoms ', diag%tdoms)
706 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%trans ', diag%trans)
707 ! CCPP/RUC only
708 if (model%lsm == model%lsm_ruc) then
709 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Sfcprop%wetness ', sfcprop%wetness)
710 else
711 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%wet1 ', diag%wet1)
712 end if
713 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%zmtnblck ', diag%zmtnblck)
714 if (model%ldiag3d) then
715 !do itracer=2,Model%ntracp100
716 ! do iprocess=1,Model%nprocess
717 ! idtend = Model%dtidx(itracer,iprocess)
718 ! if(idtend>=1) then
719 ! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, &
720 ! 'dtend_'//Model%dtend_tracer_labels(itracer)//'_' &
721 ! //Model%dtend_cause_labels(iprocess), Diag%dtend(1,1,idtend))
722 ! endif
723 ! enddo
724 !enddo
725 if (model%qdiag3d) then
726 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%upd_mf ', diag%upd_mf)
727 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dwn_mf ', diag%dwn_mf)
728 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%det_mf ', diag%det_mf)
729 end if
730 end if
731 if(model%lradar) then
732 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%refl_10cm ', diag%refl_10cm)
733 end if
734 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dkt ', diag%dkt)
735 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dku ', diag%dku)
736 ! CCPP/MYNNEDMF only
737 if (model%do_mynnedmf) then
738 if (model%bl_mynn_output .ne. 0) then
739 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%edmf_a ', diag%edmf_a)
740 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%edmf_w ', diag%edmf_w)
741 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%edmf_qt ', diag%edmf_qt)
742 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%edmf_thl ', diag%edmf_thl)
743 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%edmf_ent ', diag%edmf_ent)
744 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%edmf_qc ', diag%edmf_qc)
745 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%sub_thl ', diag%sub_thl)
746 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%sub_sqv ', diag%sub_sqv)
747 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%det_thl ', diag%det_thl)
748 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%det_sqv ', diag%det_sqv)
749 end if
750 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%maxwidth ', diag%maxwidth)
751 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%maxMF ', diag%maxMF)
752 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%ztop_plume ', diag%ztop_plume)
753 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%exch_h ', diag%exch_h)
754 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%exch_m ', diag%exch_m)
755 end if
756 ! UGWP - incomplete list
757 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dudt_gw ', diag%dudt_gw)
758 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvdt_gw ', diag%dvdt_gw)
759 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dtdt_gw ', diag%dtdt_gw)
760 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%kdis_gw ', diag%kdis_gw)
761 if (model%do_ugwp_v1 .or. model%gwd_opt==33 .or. model%gwd_opt==22) then
762 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dudt_ogw ', diag%dudt_ogw )
763 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvdt_ogw ', diag%dvdt_ogw )
764 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dudt_obl ', diag%dudt_obl )
765 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvdt_obl ', diag%dvdt_obl )
766 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dudt_oss ', diag%dudt_oss )
767 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvdt_oss ', diag%dvdt_oss )
768 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dudt_ofd ', diag%dudt_ofd )
769 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dvdt_ofd ', diag%dvdt_ofd )
770 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%du_ogwcol ', diag%du_ogwcol)
771 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dv_ogwcol ', diag%dv_ogwcol)
772 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%du_oblcol ', diag%du_oblcol)
773 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dv_oblcol ', diag%dv_oblcol)
774 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%du_osscol ', diag%du_osscol)
775 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dv_osscol ', diag%dv_osscol)
776 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%du_ofdcol ', diag%du_ofdcol)
777 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dv_ofdcol ', diag%dv_ofdcol)
778 else
779 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Diag%dudt_ogw ', diag%dudt_ogw)
780 end if
781 ! Statein
782 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%phii' , statein%phii)
783 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%prsi' , statein%prsi)
784 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%prsik' , statein%prsik)
785 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%phil' , statein%phil)
786 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%prsl' , statein%prsl)
787 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%prslk' , statein%prslk)
788 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%pgr' , statein%pgr)
789 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%ugrs' , statein%ugrs)
790 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%vgrs' , statein%vgrs)
791 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%vvl' , statein%vvl)
792 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%tgrs' , statein%tgrs)
793 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%qgrs' , statein%qgrs)
794 do n=1,size(statein%qgrs(1,1,:))
795 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%qgrs_n', statein%qgrs(:,:,n))
796 end do
797 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%diss_est', statein%diss_est)
798 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%smc' , statein%smc)
799 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%stc' , statein%stc)
800 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Statein%slc' , statein%slc)
801 ! Stateout
802 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Stateout%gu0', stateout%gu0)
803 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Stateout%gv0', stateout%gv0)
804 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Stateout%gt0', stateout%gt0)
805 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Stateout%gq0', stateout%gq0)
806 do n=1,size(stateout%gq0(1,1,:))
807 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Stateout%gq0_n', stateout%gq0(:,:,n))
808 end do
809 ! Coupling
810 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nirbmdi', coupling%nirbmdi)
811 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nirdfdi', coupling%nirdfdi)
812 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%visbmdi', coupling%visbmdi)
813 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%visdfdi', coupling%visdfdi)
814 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nirbmui', coupling%nirbmui)
815 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nirdfui', coupling%nirdfui)
816 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%visbmui', coupling%visbmui)
817 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%visdfui', coupling%visdfui)
818 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%sfcdsw ', coupling%sfcdsw )
819 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%sfcnsw ', coupling%sfcnsw )
820 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%sfcdlw ', coupling%sfcdlw )
821 if (model%cplflx .or. model%do_sppt .or. model%cplchm) then
822 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%rain_cpl', coupling%rain_cpl)
823 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%snow_cpl', coupling%snow_cpl)
824 end if
825! if (Model%cplwav2atm) then
826! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Coupling%zorlwav_cpl' , Coupling%zorlwav_cpl )
827! end if
828 if (model%cplflx) then
829 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%oro_cpl' , coupling%oro_cpl )
830 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%slmsk_cpl' , coupling%slmsk_cpl )
831 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%slimskin_cpl', coupling%slimskin_cpl )
832 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dusfcin_cpl ', coupling%dusfcin_cpl )
833 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvsfcin_cpl ', coupling%dvsfcin_cpl )
834 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dtsfcin_cpl ', coupling%dtsfcin_cpl )
835 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dqsfcin_cpl ', coupling%dqsfcin_cpl )
836 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ulwsfcin_cpl', coupling%ulwsfcin_cpl )
837! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Coupling%tseain_cpl ', Coupling%tseain_cpl )
838! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Coupling%tisfcin_cpl ', Coupling%tisfcin_cpl )
839! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Coupling%ficein_cpl ', Coupling%ficein_cpl )
840! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Coupling%hicein_cpl ', Coupling%hicein_cpl )
841 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%hsnoin_cpl ', coupling%hsnoin_cpl )
842 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dusfc_cpl ', coupling%dusfc_cpl )
843 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvsfc_cpl ', coupling%dvsfc_cpl )
844 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dtsfc_cpl ', coupling%dtsfc_cpl )
845 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dqsfc_cpl ', coupling%dqsfc_cpl )
846 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dlwsfc_cpl ', coupling%dlwsfc_cpl )
847 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dswsfc_cpl ', coupling%dswsfc_cpl )
848 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dnirbm_cpl ', coupling%dnirbm_cpl )
849 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dnirdf_cpl ', coupling%dnirdf_cpl )
850 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvisbm_cpl ', coupling%dvisbm_cpl )
851 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvisdf_cpl ', coupling%dvisdf_cpl )
852 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nlwsfc_cpl ', coupling%nlwsfc_cpl )
853 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nswsfc_cpl ', coupling%nswsfc_cpl )
854 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nnirbm_cpl ', coupling%nnirbm_cpl )
855 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nnirdf_cpl ', coupling%nnirdf_cpl )
856 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nvisbm_cpl ', coupling%nvisbm_cpl )
857 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nvisdf_cpl ', coupling%nvisdf_cpl )
858 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dusfci_cpl ', coupling%dusfci_cpl )
859 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvsfci_cpl ', coupling%dvsfci_cpl )
860 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dtsfci_cpl ', coupling%dtsfci_cpl )
861 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dqsfci_cpl ', coupling%dqsfci_cpl )
862 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dlwsfci_cpl ', coupling%dlwsfci_cpl )
863 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dswsfci_cpl ', coupling%dswsfci_cpl )
864 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dnirbmi_cpl ', coupling%dnirbmi_cpl )
865 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dnirdfi_cpl ', coupling%dnirdfi_cpl )
866 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvisbmi_cpl ', coupling%dvisbmi_cpl )
867 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvisdfi_cpl ', coupling%dvisdfi_cpl )
868 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nlwsfci_cpl ', coupling%nlwsfci_cpl )
869 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nswsfci_cpl ', coupling%nswsfci_cpl )
870 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nnirbmi_cpl ', coupling%nnirbmi_cpl )
871 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nnirdfi_cpl ', coupling%nnirdfi_cpl )
872 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nvisbmi_cpl ', coupling%nvisbmi_cpl )
873 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nvisdfi_cpl ', coupling%nvisdfi_cpl )
874 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%t2mi_cpl ', coupling%t2mi_cpl )
875 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%q2mi_cpl ', coupling%q2mi_cpl )
876 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%u10mi_cpl ', coupling%u10mi_cpl )
877 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%v10mi_cpl ', coupling%v10mi_cpl )
878 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%tsfci_cpl ', coupling%tsfci_cpl )
879 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%psurfi_cpl ', coupling%psurfi_cpl )
880 if (model%use_med_flux) then
881 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dusfcino_cpl ', coupling%dusfcino_cpl )
882 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dvsfcino_cpl ', coupling%dvsfcino_cpl )
883 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dtsfcino_cpl ', coupling%dtsfcino_cpl )
884 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%dqsfcino_cpl ', coupling%dqsfcino_cpl )
885 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ulwsfcino_cpl', coupling%ulwsfcino_cpl )
886 end if
887 end if
888 if (model%cplchm) then
889 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%rainc_cpl', coupling%rainc_cpl)
890 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ushfsfci ', coupling%ushfsfci )
891 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%pfi_lsan', coupling%pfi_lsan )
892 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%pfl_lsan', coupling%pfl_lsan )
893 end if
894 if (model%do_sppt) then
895 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%sppt_wts', coupling%sppt_wts)
896 end if
897 if (model%do_shum) then
898 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%shum_wts', coupling%shum_wts)
899 end if
900 if (model%do_skeb) then
901 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%skebu_wts', coupling%skebu_wts )
902 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%skebv_wts', coupling%skebv_wts )
903 end if
904 if (model%lndp_type /= 0) then
905 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%sfc_wts' , coupling%sfc_wts )
906 end if
907 if (model%do_ca) then
908 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ca1 ', coupling%ca1 )
909 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ca_deep ', coupling%ca_deep )
910 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ca_turb ', coupling%ca_turb )
911 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ca_shal ', coupling%ca_shal )
912 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ca_rad ', coupling%ca_rad )
913 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%ca_micro ', coupling%ca_micro )
914 end if
915 if(model%imp_physics == model%imp_physics_thompson .and. model%ltaerosol) then
916 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nwfa2d', coupling%nwfa2d)
917 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%nifa2d', coupling%nifa2d)
918 end if
919 if (model%do_RRTMGP) then
920 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%fluxlwUP_jac', coupling%fluxlwUP_jac)
921 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Coupling%htrlw', coupling%htrlw)
922 end if
923 !
924 ! Grid
925 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%xlon ', grid%xlon )
926 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%xlat ', grid%xlat )
927 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%xlat_d', grid%xlat_d)
928 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%sinlat', grid%sinlat)
929 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%coslat', grid%coslat)
930 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%area ', grid%area )
931 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%dx ', grid%dx )
932 if (model%kdt>0 .and. model%ntoz>0) then
933 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%ddy_o3 ', grid%ddy_o3 )
934 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%jindx1_o3', grid%jindx1_o3)
935 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%jindx2_o3', grid%jindx2_o3)
936 endif
937 if (model%kdt>0 .and. model%h2o_phys) then
938 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%ddy_h ', grid%ddy_h )
939 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%jindx1_h', grid%jindx1_h)
940 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%jindx2_h', grid%jindx2_h)
941 endif
942 if (model%kdt>0 .and. model%do_ugwp_v1) then
943 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%ddy_j1tau ', grid%ddy_j1tau )
944 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%ddy_j2tau ', grid%ddy_j2tau )
945 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%jindx1_tau', grid%jindx1_tau )
946 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Grid%jindx2_tau', grid%jindx2_tau )
947 endif
948 end if
949#ifdef _OPENMP
950!$OMP BARRIER
951#endif
952 end do
953#ifdef MPI
954! call MPI_BARRIER(mpicomm,ierr)
955#endif
956 end do
957
958#ifdef _OPENMP
959!$OMP BARRIER
960#endif
961#ifdef MPI
962! call MPI_BARRIER(mpicomm,ierr)
963#endif
964
965 end subroutine gfs_diagtoscreen_run
966
967 end module gfs_diagtoscreen
968
969
971
972 use print_var_chksum, only: print_var
973
974 implicit none
975
976 private
977
978 public gfs_interstitialtoscreen_init, gfs_interstitialtoscreen_timestep_init, gfs_interstitialtoscreen_run
979
980 contains
981
985 subroutine gfs_interstitialtoscreen_init (Model, Data, Interstitial, errmsg, errflg)
986
987 use gfs_typedefs, only: gfs_control_type, gfs_data_type
988 use ccpp_typedefs, only: gfs_interstitial_type
989
990 implicit none
991
992 !--- interface variables
993 type(gfs_control_type), intent(in) :: model
994 type(gfs_data_type), intent(in) :: data(:)
995 type(gfs_interstitial_type), intent(in) :: interstitial(:)
996 character(len=*), intent(out) :: errmsg
997 integer, intent(out) :: errflg
998
999 !--- local variables
1000 integer :: i
1001
1002 ! Initialize CCPP error handling variables
1003 errmsg = ''
1004 errflg = 0
1005
1006
1007 do i=1,size(interstitial)
1008 call gfs_interstitialtoscreen_run (model, Data(1)%Statein, Data(1)%Stateout, Data(1)%Sfcprop, &
1009 Data(1)%Coupling, Data(1)%Grid, Data(1)%Tbd, Data(1)%Cldprop, &
1010 Data(1)%Radtend, Data(1)%Intdiag, interstitial(i), &
1011 size(interstitial), -999, errmsg, errflg)
1012 end do
1013
1014 end subroutine gfs_interstitialtoscreen_init
1015
1019 subroutine gfs_interstitialtoscreen_timestep_init (Model, Data, Interstitial, errmsg, errflg)
1020
1021 use gfs_typedefs, only: gfs_control_type, gfs_data_type
1022 use ccpp_typedefs, only: gfs_interstitial_type
1023
1024 implicit none
1025
1026 !--- interface variables
1027 type(gfs_control_type), intent(in) :: model
1028 type(gfs_data_type), intent(in) :: data(:)
1029 type(gfs_interstitial_type), intent(in) :: interstitial(:)
1030 character(len=*), intent(out) :: errmsg
1031 integer, intent(out) :: errflg
1032
1033 !--- local variables
1034 integer :: i
1035
1036 ! Initialize CCPP error handling variables
1037 errmsg = ''
1038 errflg = 0
1039
1040
1041 do i=1,size(interstitial)
1042 call gfs_interstitialtoscreen_run (model, Data(1)%Statein, Data(1)%Stateout, Data(1)%Sfcprop, &
1043 Data(1)%Coupling, Data(1)%Grid, Data(1)%Tbd, Data(1)%Cldprop, &
1044 Data(1)%Radtend, Data(1)%Intdiag, interstitial(i), &
1045 size(interstitial), -999, errmsg, errflg)
1046 end do
1047
1048 end subroutine gfs_interstitialtoscreen_timestep_init
1049
1053 subroutine gfs_interstitialtoscreen_run (Model, Statein, Stateout, Sfcprop, Coupling, &
1054 Grid, Tbd, Cldprop, Radtend, Diag, Interstitial, &
1055 nthreads, blkno, errmsg, errflg)
1056
1057#ifdef MPI
1058 use mpi_f08
1059#endif
1060#ifdef _OPENMP
1061 use omp_lib
1062#endif
1063 use machine, only: kind_phys
1064 use gfs_typedefs, only: gfs_control_type, gfs_statein_type, &
1065 gfs_stateout_type, gfs_sfcprop_type, &
1066 gfs_coupling_type, gfs_grid_type, &
1067 gfs_tbd_type, gfs_cldprop_type, &
1068 gfs_radtend_type, gfs_diag_type
1069 use ccpp_typedefs, only: gfs_interstitial_type
1070
1071 implicit none
1072
1073 !--- interface variables
1074 type(gfs_control_type), intent(in ) :: model
1075 type(gfs_statein_type), intent(in ) :: statein
1076 type(gfs_stateout_type), intent(in ) :: stateout
1077 type(gfs_sfcprop_type), intent(in ) :: sfcprop
1078 type(gfs_coupling_type), intent(in ) :: coupling
1079 type(gfs_grid_type), intent(in ) :: grid
1080 type(gfs_tbd_type), intent(in ) :: tbd
1081 type(gfs_cldprop_type), intent(in ) :: cldprop
1082 type(gfs_radtend_type), intent(in ) :: radtend
1083 type(gfs_diag_type), intent(in ) :: diag
1084 type(gfs_interstitial_type), intent(in) :: interstitial
1085 integer, intent(in ) :: nthreads
1086 integer, intent(in ) :: blkno
1087 character(len=*), intent( out) :: errmsg
1088 integer, intent( out) :: errflg
1089
1090 !--- local variables
1091 integer :: impi, iomp, ierr
1092 integer :: mpirank, mpisize, mpicomm
1093 integer :: omprank, ompsize
1094 integer :: istart, iend, kstart, kend
1095
1096 ! Initialize CCPP error handling variables
1097 errmsg = ''
1098 errflg = 0
1099
1100#ifdef MPI
1101 mpicomm = model%communicator
1102 mpirank = model%me
1103 call mpi_comm_size(mpicomm, mpisize, ierr)
1104#else
1105 mpirank = 0
1106 mpisize = 1
1107 mpicomm = 0
1108#endif
1109#ifdef _OPENMP
1110 omprank = omp_get_thread_num()
1111 ompsize = nthreads
1112#else
1113 omprank = 0
1114 ompsize = 1
1115#endif
1116
1117#ifdef _OPENMP
1118!$OMP BARRIER
1119#endif
1120#ifdef MPI
1121! call MPI_BARRIER(mpicomm,ierr)
1122#endif
1123
1124 do impi=0,mpisize-1
1125 do iomp=0,ompsize-1
1126 if (mpirank==impi .and. omprank==iomp) then
1127 ! Print static variables
1128 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ipr ', interstitial%ipr )
1129 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%itc ', interstitial%itc )
1130 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%latidxprnt ', interstitial%latidxprnt )
1131 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%levi ', interstitial%levi )
1132 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%lmk ', interstitial%lmk )
1133 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%lmp ', interstitial%lmp )
1134 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nbdlw ', interstitial%nbdlw )
1135 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nbdsw ', interstitial%nbdsw )
1136 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nf_aelw ', interstitial%nf_aelw )
1137 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nf_aesw ', interstitial%nf_aesw )
1138 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nsamftrac ', interstitial%nsamftrac )
1139 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nscav ', interstitial%nscav )
1140 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nspc1 ', interstitial%nspc1 )
1141 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ntiwx ', interstitial%ntiwx )
1142 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nvdiff ', interstitial%nvdiff )
1143 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%phys_hydrostatic ', interstitial%phys_hydrostatic )
1144 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%skip_macro ', interstitial%skip_macro )
1145 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%trans_aero ', interstitial%trans_aero )
1146 ! Print all other variables
1147 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjsfculw_land ', interstitial%adjsfculw_land )
1148 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjsfculw_ice ', interstitial%adjsfculw_ice )
1149 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjsfculw_water ', interstitial%adjsfculw_water )
1150 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjnirbmd ', interstitial%adjnirbmd )
1151 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjnirbmu ', interstitial%adjnirbmu )
1152 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjnirdfd ', interstitial%adjnirdfd )
1153 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjnirdfu ', interstitial%adjnirdfu )
1154 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjvisbmd ', interstitial%adjvisbmd )
1155 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjvisbmu ', interstitial%adjvisbmu )
1156 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjvisdfu ', interstitial%adjvisdfu )
1157 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%adjvisdfd ', interstitial%adjvisdfd )
1158 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%aerodp ', interstitial%aerodp )
1159 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%alb1d ', interstitial%alb1d )
1160 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%bexp1d ', interstitial%bexp1d )
1161 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cd ', interstitial%cd )
1162 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cd_ice ', interstitial%cd_ice )
1163 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cd_land ', interstitial%cd_land )
1164 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cd_water ', interstitial%cd_water )
1165 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cdq ', interstitial%cdq )
1166 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cdq_ice ', interstitial%cdq_ice )
1167 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cdq_land ', interstitial%cdq_land )
1168 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cdq_water ', interstitial%cdq_water )
1169 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%chh_ice ', interstitial%chh_ice )
1170 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%chh_land ', interstitial%chh_land )
1171 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%chh_water ', interstitial%chh_water )
1172 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cldf ', interstitial%cldf )
1173 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cldsa ', interstitial%cldsa )
1174 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cldtaulw ', interstitial%cldtaulw )
1175 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cldtausw ', interstitial%cldtausw )
1176 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld1d ', interstitial%cld1d )
1177 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%clw ', interstitial%clw )
1178 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%clx ', interstitial%clx )
1179 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%clouds ', interstitial%clouds )
1180 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cmm_ice ', interstitial%cmm_ice )
1181 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cmm_land ', interstitial%cmm_land )
1182 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cmm_water ', interstitial%cmm_water )
1183 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnvc ', interstitial%cnvc )
1184 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnvw ', interstitial%cnvw )
1185 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ctei_r ', interstitial%ctei_r )
1186 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ctei_rml ', interstitial%ctei_rml )
1187 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cumabs ', interstitial%cumabs )
1188 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dd_mf ', interstitial%dd_mf )
1189 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%de_lgth ', interstitial%de_lgth )
1190 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%del ', interstitial%del )
1191 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%del_gz ', interstitial%del_gz )
1192 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%delr ', interstitial%delr )
1193 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dlength ', interstitial%dlength )
1194 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dqdt ', interstitial%dqdt )
1195 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dqsfc1 ', interstitial%dqsfc1 )
1196 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%drain ', interstitial%drain )
1197 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dtdt ', interstitial%dtdt )
1198 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dtsfc1 ', interstitial%dtsfc1 )
1199 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dtzm ', interstitial%dtzm )
1200 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dt_mf ', interstitial%dt_mf )
1201 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dudt ', interstitial%dudt )
1202 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dusfcg ', interstitial%dusfcg )
1203 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dusfc1 ', interstitial%dusfc1 )
1204 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dvdftra ', interstitial%dvdftra )
1205 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dvdt ', interstitial%dvdt )
1206 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dvsfcg ', interstitial%dvsfcg )
1207 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dvsfc1 ', interstitial%dvsfc1 )
1208 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dzlyr ', interstitial%dzlyr )
1209 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%elvmax ', interstitial%elvmax )
1210 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ep1d ', interstitial%ep1d )
1211 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ep1d_ice ', interstitial%ep1d_ice )
1212 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ep1d_land ', interstitial%ep1d_land )
1213 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ep1d_water ', interstitial%ep1d_water )
1214 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%evap_ice ', interstitial%evap_ice )
1215 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%evap_land ', interstitial%evap_land )
1216 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%evap_water ', interstitial%evap_water )
1217 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ext_diag_thompson_reset', interstitial%ext_diag_thompson_reset)
1218 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%faerlw ', interstitial%faerlw )
1219 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%faersw ', interstitial%faersw )
1220 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ffhh_ice ', interstitial%ffhh_ice )
1221 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ffhh_land ', interstitial%ffhh_land )
1222 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ffhh_water ', interstitial%ffhh_water )
1223 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fh2 ', interstitial%fh2 )
1224 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fh2_ice ', interstitial%fh2_ice )
1225 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fh2_land ', interstitial%fh2_land )
1226 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fh2_water ', interstitial%fh2_water )
1227 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%flag_cice ', interstitial%flag_cice )
1228 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%flag_guess ', interstitial%flag_guess )
1229 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%flag_iter ', interstitial%flag_iter )
1230 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ffmm_ice ', interstitial%ffmm_ice )
1231 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ffmm_land ', interstitial%ffmm_land )
1232 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ffmm_water ', interstitial%ffmm_water )
1233 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fm10 ', interstitial%fm10 )
1234 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fm10_ice ', interstitial%fm10_ice )
1235 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fm10_land ', interstitial%fm10_land )
1236 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fm10_water ', interstitial%fm10_water )
1237 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%frain ', interstitial%frain )
1238 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%frland ', interstitial%frland )
1239 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fscav ', interstitial%fscav )
1240 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fswtr ', interstitial%fswtr )
1241 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gabsbdlw ', interstitial%gabsbdlw )
1242 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gabsbdlw_ice ', interstitial%gabsbdlw_ice )
1243 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gabsbdlw_land ', interstitial%gabsbdlw_land )
1244 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gabsbdlw_water ', interstitial%gabsbdlw_water )
1245 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gamma ', interstitial%gamma )
1246 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gamq ', interstitial%gamq )
1247 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gamt ', interstitial%gamt )
1248 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gasvmr ', interstitial%gasvmr )
1249 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gflx ', interstitial%gflx )
1250 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gflx_ice ', interstitial%gflx_ice )
1251 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gflx_land ', interstitial%gflx_land )
1252 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gflx_water ', interstitial%gflx_water )
1253 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gwdcu ', interstitial%gwdcu )
1254 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%gwdcv ', interstitial%gwdcv )
1255 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%zvfun ', interstitial%zvfun )
1256 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%hffac ', interstitial%hffac )
1257 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%hflxq ', interstitial%hflxq )
1258 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%hflx_ice ', interstitial%hflx_ice )
1259 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%hflx_land ', interstitial%hflx_land )
1260 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%hflx_water ', interstitial%hflx_water )
1261 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%htlwc ', interstitial%htlwc )
1262 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%htlw0 ', interstitial%htlw0 )
1263 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%htswc ', interstitial%htswc )
1264 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%htsw0 ', interstitial%htsw0 )
1265 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dry ', interstitial%dry )
1266 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%idxday ', interstitial%idxday )
1267 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%icy ', interstitial%icy )
1268 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%lake ', interstitial%lake )
1269 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ocean ', interstitial%ocean )
1270 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%islmsk ', interstitial%islmsk )
1271 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%islmsk_cice ', interstitial%islmsk_cice )
1272 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%wet ', interstitial%wet )
1273 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kb ', interstitial%kb )
1274 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kbot ', interstitial%kbot )
1275 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kcnv ', interstitial%kcnv )
1276 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kd ', interstitial%kd )
1277 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kinver ', interstitial%kinver )
1278 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kpbl ', interstitial%kpbl )
1279 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kt ', interstitial%kt )
1280 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ktop ', interstitial%ktop )
1281 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%max_hourly_reset ', interstitial%max_hourly_reset )
1282 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%mbota ', interstitial%mbota )
1283 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%mtopa ', interstitial%mtopa )
1284 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%nday ', interstitial%nday )
1285 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%oa4 ', interstitial%oa4 )
1286 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%oc ', interstitial%oc )
1287 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%olyr ', interstitial%olyr )
1288 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%plvl ', interstitial%plvl )
1289 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%plyr ', interstitial%plyr )
1290 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%prcpmp ', interstitial%prcpmp )
1291 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%prnum ', interstitial%prnum )
1292 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qlyr ', interstitial%qlyr )
1293 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qss_ice ', interstitial%qss_ice )
1294 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qss_land ', interstitial%qss_land )
1295 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qss_water ', interstitial%qss_water )
1296 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%radar_reset ', interstitial%radar_reset )
1297 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%raddt ', interstitial%raddt )
1298 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%raincd ', interstitial%raincd )
1299 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%raincs ', interstitial%raincs )
1300 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rainmcadj ', interstitial%rainmcadj )
1301 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rainp ', interstitial%rainp )
1302 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rb ', interstitial%rb )
1303 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rb_ice ', interstitial%rb_ice )
1304 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rb_land ', interstitial%rb_land )
1305 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rb_water ', interstitial%rb_water )
1306 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rhc ', interstitial%rhc )
1307 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%runoff ', interstitial%runoff )
1308 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%save_q ', interstitial%save_q )
1309 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%save_t ', interstitial%save_t )
1310 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%save_tcp ', interstitial%save_tcp )
1311 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%save_u ', interstitial%save_u )
1312 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%save_v ', interstitial%save_v )
1313 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%scmpsw%uvbfc ', interstitial%scmpsw%uvbfc )
1314 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%scmpsw%uvbf0 ', interstitial%scmpsw%uvbf0 )
1315 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%scmpsw%nirbm ', interstitial%scmpsw%nirbm )
1316 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%scmpsw%nirdf ', interstitial%scmpsw%nirdf )
1317 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%scmpsw%visbm ', interstitial%scmpsw%visbm )
1318 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%scmpsw%visdf ', interstitial%scmpsw%visdf )
1319 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%sfcalb ', interstitial%sfcalb )
1320 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%sigma ', interstitial%sigma )
1321 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%sigmaf ', interstitial%sigmaf )
1322 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%sigmafrac ', interstitial%sigmafrac )
1323 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%sigmatot ', interstitial%sigmatot )
1324 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%snowc ', interstitial%snowc )
1325 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%snowd_ice ', interstitial%snowd_ice )
1326 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%snowd_land ', interstitial%snowd_land )
1327 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%snohf ', interstitial%snohf )
1328 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%snowmt ', interstitial%snowmt )
1329 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%stress ', interstitial%stress )
1330 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%stress_ice ', interstitial%stress_ice )
1331 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%stress_land ', interstitial%stress_land )
1332 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%stress_water ', interstitial%stress_water )
1333 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%theta ', interstitial%theta )
1334 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tlvl ', interstitial%tlvl )
1335 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tlyr ', interstitial%tlyr )
1336 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tprcp_ice ', interstitial%tprcp_ice )
1337 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tprcp_land ', interstitial%tprcp_land )
1338 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tprcp_water ', interstitial%tprcp_water )
1339 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tseal ', interstitial%tseal )
1340 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tsfa ', interstitial%tsfa )
1341 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tsfc_water ', interstitial%tsfc_water )
1342 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tsfg ', interstitial%tsfg )
1343 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tsurf_ice ', interstitial%tsurf_ice )
1344 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tsurf_land ', interstitial%tsurf_land )
1345 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tsurf_water ', interstitial%tsurf_water )
1346 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%uustar_ice ', interstitial%uustar_ice )
1347 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%uustar_land ', interstitial%uustar_land )
1348 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%uustar_water ', interstitial%uustar_water )
1349 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%vdftra ', interstitial%vdftra )
1350 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%vegf1d ', interstitial%vegf1d )
1351 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%wcbmax ', interstitial%wcbmax )
1352! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%weasd_ice ', Interstitial%weasd_ice )
1353! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%weasd_land ', Interstitial%weasd_land )
1354! call print_var(mpirank, omprank, blkno, Grid%xlat_d, Grid%xlon_d, 'Interstitial%weasd_water ', Interstitial%weasd_water )
1355 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%wind ', interstitial%wind )
1356 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%work1 ', interstitial%work1 )
1357 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%work2 ', interstitial%work2 )
1358 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%work3 ', interstitial%work3 )
1359 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%xcosz ', interstitial%xcosz )
1360 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%xlai1d ', interstitial%xlai1d )
1361 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%xmu ', interstitial%xmu )
1362 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%z01d ', interstitial%z01d )
1363 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%zt1d ', interstitial%zt1d )
1364 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ztmax_ice ', interstitial%ztmax_ice )
1365 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ztmax_land ', interstitial%ztmax_land )
1366 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ztmax_water ', interstitial%ztmax_water )
1367 ! UGWP
1368 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tau_mtb ', interstitial%tau_mtb )
1369 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tau_ogw ', interstitial%tau_ogw )
1370 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tau_tofd ', interstitial%tau_tofd )
1371 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tau_ngw ', interstitial%tau_ngw )
1372 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tau_oss ', interstitial%tau_oss )
1373 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dudt_mtb ', interstitial%dudt_mtb )
1374 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dudt_tms ', interstitial%dudt_tms )
1375 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%zmtb ', interstitial%zmtb )
1376 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%zlwb ', interstitial%zlwb )
1377 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%zogw ', interstitial%zogw )
1378 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%zngw ', interstitial%zngw )
1379 ! UGWP v1
1380 if (model%do_ugwp_v1) then
1381 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dudt_ngw ', interstitial%dudt_ngw )
1382 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dvdt_ngw ', interstitial%dvdt_ngw )
1383 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%dtdt_ngw ', interstitial%dtdt_ngw )
1384 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%kdis_ngw ', interstitial%kdis_ngw )
1385 end if
1386 !-- GSD drag suite
1387 if (model%gwd_opt==3 .or. model%gwd_opt==33 .or. &
1388 model%gwd_opt==2 .or. model%gwd_opt==22) then
1389 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%varss ', interstitial%varss )
1390 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ocss ', interstitial%ocss )
1391 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%oa4ss ', interstitial%oa4ss )
1392 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%clxss ', interstitial%clxss )
1393 end if
1394 ! GFDL and Thompson MP
1395 if (model%imp_physics == model%imp_physics_gfdl .or. model%imp_physics == model%imp_physics_thompson .or. model%imp_physics == model%imp_physics_nssl) then
1396 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%graupelmp ', interstitial%graupelmp )
1397 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%icemp ', interstitial%icemp )
1398 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%rainmp ', interstitial%rainmp )
1399 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%snowmp ', interstitial%snowmp )
1400 ! Morrison-Gettelman
1401 else if (model%imp_physics == model%imp_physics_mg) then
1402 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ncgl ', interstitial%ncgl )
1403 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ncpr ', interstitial%ncpr )
1404 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ncps ', interstitial%ncps )
1405 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qgl ', interstitial%qgl )
1406 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qrn ', interstitial%qrn )
1407 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qsnw ', interstitial%qsnw )
1408 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qlcn ', interstitial%qlcn )
1409 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qicn ', interstitial%qicn )
1410 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%w_upi ', interstitial%w_upi )
1411 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cf_upi ', interstitial%cf_upi )
1412 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnv_mfd ', interstitial%cnv_mfd )
1413 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnv_dqldt ', interstitial%cnv_dqldt )
1414 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%clcn ', interstitial%clcn )
1415 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnv_fice ', interstitial%cnv_fice )
1416 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnv_ndrop ', interstitial%cnv_ndrop )
1417 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cnv_nice ', interstitial%cnv_nice )
1418 end if
1419 ! SHOC
1420 if (model%do_shoc) then
1421 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ncgl ', interstitial%ncgl )
1422 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qrn ', interstitial%qrn )
1423 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qsnw ', interstitial%qsnw )
1424 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qgl ', interstitial%qgl )
1425 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ncpi ', interstitial%ncpi )
1426 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%ncpl ', interstitial%ncpl )
1427 end if
1428 ! Noah MP
1429 if (model%lsm == model%lsm_noahmp) then
1430 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%t2mmp ', interstitial%t2mmp )
1431 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%q2mp ', interstitial%q2mp )
1432 end if
1433 ! RRTMGP
1434 if (model%do_RRTMGP) then
1435 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%aerosolslw ', interstitial%aerosolslw )
1436 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%aerosolssw ', interstitial%aerosolssw )
1437 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_frac ', interstitial%cld_frac )
1438 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_lwp ', interstitial%cld_lwp )
1439 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_reliq ', interstitial%cld_reliq )
1440 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_iwp ', interstitial%cld_iwp )
1441 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_reice ', interstitial%cld_reice )
1442 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_swp ', interstitial%cld_swp )
1443 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_resnow ', interstitial%cld_resnow )
1444 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_rwp ', interstitial%cld_rwp )
1445 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cld_rerain ', interstitial%cld_rerain )
1446 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%precip_frac ', interstitial%precip_frac )
1447 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxlwUP_allsky ', interstitial%fluxlwUP_allsky )
1448 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxlwDOWN_allsky ', interstitial%fluxlwDOWN_allsky )
1449 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxlwUP_clrsky ', interstitial%fluxlwUP_clrsky )
1450 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxlwDOWN_clrsky ', interstitial%fluxlwDOWN_clrsky )
1451 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxswUP_allsky ', interstitial%fluxswUP_allsky )
1452 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxswDOWN_allsky ', interstitial%fluxswDOWN_allsky )
1453 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxswUP_clrsky ', interstitial%fluxswUP_clrsky )
1454 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%fluxswDOWN_clrsky ', interstitial%fluxswDOWN_clrsky )
1455 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%relhum ', interstitial%relhum )
1456 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%q_lay ', interstitial%q_lay )
1457 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%qs_lay ', interstitial%qs_lay )
1458 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%deltaZ ', interstitial%deltaZ )
1459 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%p_lay ', interstitial%p_lay )
1460 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%p_lev ', interstitial%p_lev )
1461 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%t_lay ', interstitial%t_lay )
1462 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%t_lev ', interstitial%t_lev )
1463 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%tv_lay ', interstitial%tv_lay )
1464 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%cloud_overlap_param ', interstitial%cloud_overlap_param )
1465 call print_var(mpirank, omprank, blkno, grid%xlat_d, grid%xlon_d, 'Interstitial%precip_overlap_param', interstitial%precip_overlap_param )
1466 end if
1467 end if
1468#ifdef _OPENMP
1469!$OMP BARRIER
1470#endif
1471 end do
1472#ifdef MPI
1473! call MPI_BARRIER(mpicomm,ierr)
1474#endif
1475 end do
1476
1477#ifdef _OPENMP
1478!$OMP BARRIER
1479#endif
1480#ifdef MPI
1481! call MPI_BARRIER(mpicomm,ierr)
1482#endif
1483
1484 end subroutine gfs_interstitialtoscreen_run
1485
1486 end module gfs_interstitialtoscreen
1487
1489
1490 private
1491
1492 public gfs_abort_run
1493
1494 contains
1495
1499 subroutine gfs_abort_run (Model, blkno, errmsg, errflg)
1500
1501 use machine, only: kind_phys
1502 use gfs_typedefs, only: gfs_control_type
1503
1504 implicit none
1505
1506 !--- interface variables
1507 type(gfs_control_type), intent(in ) :: model
1508 integer, intent(in ) :: blkno
1509 character(len=*), intent( out) :: errmsg
1510 integer, intent( out) :: errflg
1511
1512 ! Initialize CCPP error handling variables
1513 errmsg = ''
1514 errflg = 0
1515
1516 if (model%kdt==1 .and. blkno==size(model%blksz)) then
1517 if (model%me==model%master) write(0,*) "GFS_abort_run: ABORTING MODEL"
1518 call sleep(10)
1519 stop
1520 end if
1521
1522 end subroutine gfs_abort_run
1523
1524 end module gfs_abort
1525
1527
1528 private
1529
1530 public gfs_checkland_run
1531
1532 contains
1533
1537 subroutine gfs_checkland_run (me, master, blkno, im, kdt, iter, flag_iter, flag_guess, &
1538 flag_init, flag_restart, frac_grid, isot, ivegsrc, stype,scolor, vtype, slope, &
1539 dry, icy, wet, lake, ocean, oceanfrac, landfrac, lakefrac, slmsk, islmsk, &
1540 zorl, zorlw, zorll, zorli, fice, errmsg, errflg )
1541
1542 use machine, only: kind_phys
1543
1544 implicit none
1545
1546 ! Interface variables
1547 integer, intent(in ) :: me
1548 integer, intent(in ) :: master
1549 integer, intent(in ) :: blkno
1550 integer, intent(in ) :: im
1551 integer, intent(in ) :: kdt
1552 integer, intent(in ) :: iter
1553 logical, intent(in ) :: flag_iter(:)
1554 logical, intent(in ) :: flag_guess(:)
1555 logical, intent(in ) :: flag_init
1556 logical, intent(in ) :: flag_restart
1557 logical, intent(in ) :: frac_grid
1558 integer, intent(in ) :: isot
1559 integer, intent(in ) :: ivegsrc
1560 integer, intent(in ) :: stype(:)
1561 integer, intent(in ) :: scolor(:)
1562
1563 integer, intent(in ) :: vtype(:)
1564 integer, intent(in ) :: slope(:)
1565 logical, intent(in ) :: dry(:)
1566 logical, intent(in ) :: icy(:)
1567 logical, intent(in ) :: wet(:)
1568 logical, intent(in ) :: lake(:)
1569 logical, intent(in ) :: ocean(:)
1570 real(kind_phys), intent(in ) :: oceanfrac(:)
1571 real(kind_phys), intent(in ) :: landfrac(:)
1572 real(kind_phys), intent(in ) :: lakefrac(:)
1573 real(kind_phys), intent(in ) :: slmsk(:)
1574 integer, intent(in ) :: islmsk(:)
1575 real(kind_phys), intent(in ) :: zorl(:)
1576 real(kind_phys), intent(in ) :: zorlw(:)
1577 real(kind_phys), intent(in ) :: zorll(:)
1578 real(kind_phys), intent(in ) :: zorli(:)
1579 real(kind_phys), intent(in ) :: fice(:)
1580 character(len=*), intent( out) :: errmsg
1581 integer, intent( out) :: errflg
1582
1583 ! Local variables
1584 integer :: i
1585
1586 errflg = 0
1587 errmsg = ''
1588
1589 write(0,'(a,i5)') 'YYY: me :', me
1590 write(0,'(a,i5)') 'YYY: master :', master
1591 write(0,'(a,i5)') 'YYY: blkno :', blkno
1592 write(0,'(a,i5)') 'YYY: im :', im
1593 write(0,'(a,i5)') 'YYY: kdt :', kdt
1594 write(0,'(a,i5)') 'YYY: iter :', iter
1595 write(0,'(a,1x,l)') 'YYY: flag_init :', flag_init
1596 write(0,'(a,1x,l)') 'YYY: flag_restart :', flag_restart
1597 write(0,'(a,1x,l)') 'YYY: frac_grid :', frac_grid
1598 write(0,'(a,i5)') 'YYY: isot :', isot
1599 write(0,'(a,i5)') 'YYY: ivegsrc :', ivegsrc
1600
1601 do i=1,im
1602 !if (fice(i)>0.999) then
1603 !if (vegtype(i)==15) then
1604 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, flag_iter(i) :', i, blkno, flag_iter(i)
1605 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, flag_guess(i) :', i, blkno, flag_guess(i)
1606 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, stype(i) :', i, blkno, stype(i)
1607
1608 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, scolor(i) :', i, blkno, scolor(i)
1609 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, vtype(i) :', i, blkno, vtype(i)
1610 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, slope(i) :', i, blkno, slope(i)
1611 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, dry(i) :', i, blkno, dry(i)
1612 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, icy(i) :', i, blkno, icy(i)
1613 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, wet(i) :', i, blkno, wet(i)
1614 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, lake(i) :', i, blkno, lake(i)
1615 write(0,'(a,2i5,1x,1x,l)') 'YYY: i, blk, ocean(i) :', i, blkno, ocean(i)
1616 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, oceanfrac(i) :', i, blkno, oceanfrac(i)
1617 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, landfrac(i) :', i, blkno, landfrac(i)
1618 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, lakefrac(i) :', i, blkno, lakefrac(i)
1619 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, fice(i) :', i, blkno, fice(i)
1620 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, slmsk(i) :', i, blkno, slmsk(i)
1621 write(0,'(a,2i5,1x,i5)') 'YYY: i, blk, islmsk(i) :', i, blkno, islmsk(i)
1622 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, zorl(i) :', i, blkno, zorl(i)
1623 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, zorlw(i) :', i, blkno, zorlw(i)
1624 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, zorli(i) :', i, blkno, zorli(i)
1625 write(0,'(a,2i5,1x,e16.7)')'YYY: i, blk, zorll(i) :', i, blkno, zorll(i)
1626 !end if
1627 end do
1628
1629 end subroutine gfs_checkland_run
1630 end module gfs_checkland
1631
1633
1634 private
1635
1636 public gfs_checktracers_init, gfs_checktracers_timestep_init, gfs_checktracers_run
1637
1638 contains
1639
1643 subroutine gfs_checktracers_init (me, master, im, levs, ntracer, kdt, qgrs, gq0, errmsg, errflg)
1644
1645 use machine, only: kind_phys
1646
1647 implicit none
1648
1649 ! Interface variables
1650 integer, intent(in ) :: me
1651 integer, intent(in ) :: master
1652 integer, intent(in ) :: im
1653 integer, intent(in ) :: levs
1654 integer, intent(in ) :: ntracer
1655 integer, intent(in ) :: kdt
1656 real(kind_phys), intent(in ) :: qgrs(:,:,:)
1657 real(kind_phys), intent(in ) :: gq0(:,:,:)
1658 character(len=*), intent( out) :: errmsg
1659 integer, intent( out) :: errflg
1660
1661 call gfs_checktracers_timestep_init (me, master, im, levs, ntracer, kdt, qgrs, gq0, errmsg, errflg)
1662
1663 end subroutine gfs_checktracers_init
1664
1668 subroutine gfs_checktracers_timestep_init (me, master, im, levs, ntracer, kdt, qgrs, gq0, errmsg, errflg)
1669
1670 use machine, only: kind_phys
1671
1672 implicit none
1673
1674 ! Interface variables
1675 integer, intent(in ) :: me
1676 integer, intent(in ) :: master
1677 integer, intent(in ) :: im
1678 integer, intent(in ) :: levs
1679 integer, intent(in ) :: ntracer
1680 integer, intent(in ) :: kdt
1681 real(kind_phys), intent(in ) :: qgrs(:,:,:)
1682 real(kind_phys), intent(in ) :: gq0(:,:,:)
1683 character(len=*), intent( out) :: errmsg
1684 integer, intent( out) :: errflg
1685
1686 ! Local variables
1687 integer :: i, k, n
1688
1689 errflg = 0
1690 errmsg = ''
1691
1692 write(0,'(a,i5)') 'YYY: me :', me
1693 write(0,'(a,i5)') 'YYY: master :', master
1694 write(0,'(a,i5)') 'YYY: im :', im
1695 write(0,'(a,i5)') 'YYY: levs :', levs
1696 write(0,'(a,i5)') 'YYY: ntracer :', ntracer
1697 write(0,'(a,i5)') 'YYY: kdt :', kdt
1698
1699 do n=1,ntracer
1700 do i=1,im
1701 do k=1,levs
1702 if (qgrs(i,k,n)<0 .or. gq0(i,k,n)<0) then
1703 write(0,'(a,4i5,1x,2e16.7)') 'YYY: blk, n, i, k, qgrs, gq0 :', -999, n, i, k, qgrs(i,k,n), gq0(i,k,n)
1704 end if
1705 end do
1706 end do
1707 end do
1708
1709 end subroutine gfs_checktracers_timestep_init
1710
1714 subroutine gfs_checktracers_run (me, master, blkno, im, levs, ntracer, kdt, qgrs, gq0, errmsg, errflg)
1715
1716 use machine, only: kind_phys
1717
1718 implicit none
1719
1720 ! Interface variables
1721 integer, intent(in ) :: me
1722 integer, intent(in ) :: master
1723 integer, intent(in ) :: blkno
1724 integer, intent(in ) :: im
1725 integer, intent(in ) :: levs
1726 integer, intent(in ) :: ntracer
1727 integer, intent(in ) :: kdt
1728 real(kind_phys), intent(in ) :: qgrs(:,:,:)
1729 real(kind_phys), intent(in ) :: gq0(:,:,:)
1730 character(len=*), intent( out) :: errmsg
1731 integer, intent( out) :: errflg
1732
1733 ! Local variables
1734 integer :: i, k, n
1735
1736 errflg = 0
1737 errmsg = ''
1738
1739 write(0,'(a,i5)') 'YYY: me :', me
1740 write(0,'(a,i5)') 'YYY: master :', master
1741 write(0,'(a,i5)') 'YYY: blkno :', blkno
1742 write(0,'(a,i5)') 'YYY: im :', im
1743 write(0,'(a,i5)') 'YYY: levs :', levs
1744 write(0,'(a,i5)') 'YYY: ntracer :', ntracer
1745 write(0,'(a,i5)') 'YYY: kdt :', kdt
1746
1747 do n=1,ntracer
1748 do i=1,im
1749 do k=1,levs
1750 if (qgrs(i,k,n)<0 .or. gq0(i,k,n)<0) then
1751 write(0,'(a,4i5,1x,2e16.7)') 'YYY: blk, n, i, k, qgrs, gq0 :', blkno, n, i, k, qgrs(i,k,n), gq0(i,k,n)
1752 end if
1753 end do
1754 end do
1755 end do
1756
1757 end subroutine gfs_checktracers_run
1758
1759 end module gfs_checktracers
subroutine count(slimsk, sno, ijmax)
This subroutine counts number of points for the four surface conditions.
Definition sfcsub.F:2631