Line data Source code
1 : module restart_dynamics
2 :
3 : !-----------------------------------------------------------------------
4 : !
5 : ! Read and write dynamics fields to restart file
6 : !
7 : ! !HISTORY:
8 : ! 2006.04.13 Sawyer Removed dependency on prognostics
9 : ! 2006.07.01 Sawyer Transitioned q3 tracers to T_TRACERS
10 : ! 2008.10.02 Edwards Added pio support
11 : !
12 : !-----------------------------------------------------------------------
13 :
14 : use shr_kind_mod, only: r8=>shr_kind_r8, r4=>shr_kind_r4
15 : use pmgrid, only: plon, plat, plev
16 : use constituents, only: pcnst, cnst_name
17 : use time_manager, only: get_curr_time
18 :
19 : use hycoef, only: init_restart_hycoef, write_restart_hycoef
20 :
21 : use dyn_comp, only: dyn_import_t, dyn_export_t
22 : use dyn_grid, only: get_horiz_grid_dim_d, get_dyn_grid_parm
23 : use dynamics_vars, only: t_fvdycore_state, t_fvdycore_grid
24 : use dyn_internal_state, only: get_dyn_state
25 :
26 : use cam_grid_support, only: cam_grid_write_attr, cam_grid_id
27 : use cam_grid_support, only: cam_grid_header_info_t
28 : use cam_pio_utils, only: pio_subsystem
29 : use pio, only: file_desc_t, io_desc_t, var_desc_t, &
30 : pio_double, pio_unlimited, pio_offset_kind, &
31 : pio_setdebuglevel, pio_setframe, &
32 : pio_def_var, pio_def_dim, &
33 : pio_inq_varid, &
34 : pio_put_var, pio_get_var, &
35 : pio_write_darray, pio_read_darray, &
36 : pio_initdecomp, pio_freedecomp
37 :
38 : use cam_logfile, only: iulog
39 : use cam_abortutils, only: endrun
40 :
41 : #if ( defined OFFLINE_DYN )
42 : use metdata, only: write_met_restart, read_met_restart
43 : #endif
44 :
45 : implicit none
46 : private
47 : save
48 :
49 : public :: &
50 : init_restart_dynamics, &
51 : write_restart_dynamics, &
52 : read_restart_dynamics
53 :
54 : integer, parameter :: namlen=16
55 :
56 : type restart_var_t
57 : real(r8), pointer :: v2d(:,:)
58 : real(r8), pointer :: v3d(:, :, :)
59 : real(r8), pointer :: v4d(:, :, :, :)
60 : real(r8), pointer :: v5d(:, :, :, :, :)
61 : real(r4), pointer :: v2dr4(:,:)
62 : real(r4), pointer :: v3dr4(:, :, :)
63 : real(r4), pointer :: v4dr4(:, :, :, :)
64 : real(r4), pointer :: v5dr4(:, :, :, :, :)
65 :
66 : type(var_desc_t), pointer :: vdesc
67 : integer :: ndims
68 : integer :: timelevels
69 : character(len=namlen) :: name
70 : end type restart_var_t
71 :
72 : integer, parameter :: restartvarcnt = 6+pcnst
73 :
74 : type(var_desc_t) :: timedesc
75 :
76 : type(restart_var_t) :: restartvars(restartvarcnt)
77 : logical :: restart_varlistw_initialized = .false.
78 : logical :: restart_varlistr_initialized = .false.
79 :
80 : !=========================================================================================
81 : contains
82 : !=========================================================================================
83 :
84 3072 : subroutine init_restart_dynamics(File, dyn_out)
85 :
86 : ! Input arguments
87 : type(File_desc_t), intent(inout) :: File
88 : type(Dyn_export_t), intent(in) :: dyn_out
89 :
90 : integer :: vdimids(2)
91 : integer :: ierr
92 : integer :: hdim1, hdim2
93 : integer :: hdimids(2)
94 : character(len=namlen) :: name
95 :
96 : integer :: alldims(3), alldims2d(2)
97 : integer :: i, timelevels
98 : integer :: ndims
99 : type(var_desc_t), pointer :: vdesc
100 : type(cam_grid_header_info_t) :: info
101 : !----------------------------------------------------------------------------
102 :
103 1536 : call init_restart_hycoef(File, vdimids)
104 :
105 1536 : call get_horiz_grid_dim_d(hdim1, hdim2)
106 1536 : ierr = PIO_Def_Dim(File, 'lon', hdim1, hdimids(1))
107 1536 : ierr = PIO_Def_Dim(File, 'lat', hdim2, hdimids(2))
108 :
109 1536 : ierr = PIO_Def_Var(File, 'time', pio_double, timedesc)
110 :
111 4608 : alldims(1:2) = hdimids(1:2)
112 1536 : alldims(3) = vdimids(1)
113 :
114 1536 : alldims2d(1:2) = hdimids(1:2)
115 :
116 1536 : call init_restart_varlistw(dyn_out)
117 :
118 379392 : do i = 1, restartvarcnt
119 :
120 377856 : call get_restart_var(File, i, name, timelevels, ndims, vdesc)
121 :
122 379392 : if (timelevels > 1) then
123 0 : call endrun('not expecting timelevels>1 in fv dycore')
124 : else
125 377856 : if (ndims==1) then
126 : ! broken i think
127 0 : ierr = PIO_Def_Var(File, name, pio_double, hdimids(2:2), vdesc)
128 377856 : else if (ndims==2) then
129 3072 : ierr = PIO_Def_Var(File, name, pio_double, alldims2d(1:2), vdesc)
130 374784 : else if (ndims==3) then
131 374784 : ierr = PIO_Def_Var(File, name, pio_double, alldims(1:3), vdesc)
132 : end if
133 377856 : call pio_setframe(File, vdesc, int(-1,kind=pio_offset_kind))
134 : end if
135 : end do
136 :
137 : #if ( defined OFFLINE_DYN )
138 : ! This write happens here rather than from write_restart_dynamics because
139 : ! only attributes are written (two filenames) and that should happen during
140 : ! the file definition phase.
141 : call write_met_restart( File )
142 : #endif
143 :
144 1536 : end subroutine init_restart_dynamics
145 :
146 : !=========================================================================================
147 :
148 4608 : subroutine write_restart_dynamics(File, dyn_out)
149 :
150 : ! arguments
151 : type(File_desc_t), intent(inout) :: File
152 : type(Dyn_export_t), intent(in) :: dyn_out
153 :
154 : ! Local workspace
155 : logical :: use_transfer
156 : integer :: ndcur, nscur
157 : integer :: hdim1, hdim2
158 : type(io_desc_t) :: iodesc2d, iodesc3d
159 : integer :: ierr
160 : real(r8) :: time, mold(1), null(0)
161 : integer :: i, s3d(1), s2d(1), ct
162 : integer(kind=pio_offset_kind) :: t
163 : integer :: ndims, isize(1), timelevels
164 : type(var_desc_t), pointer :: vdesc
165 1536 : integer, pointer :: ldof(:)
166 : character(len=namlen) :: name
167 : type (T_FVDYCORE_STATE), pointer :: dyn_state
168 : !----------------------------------------------------------------------------
169 :
170 1536 : call write_restart_hycoef(File)
171 :
172 : ! transfer is the fastest method to flatten the multidimensional arrays into the 1d needed by pio
173 : ! but it doesn't work correctly if the array is zero length...
174 :
175 1536 : use_transfer = .true.
176 : #if ( defined SPMD )
177 1536 : dyn_state => get_dyn_state()
178 1536 : if (dyn_state%grid%iam >= dyn_state%grid%npes_xy) then
179 0 : use_transfer = .false.
180 : end if
181 : #endif
182 :
183 1536 : call get_curr_time(ndcur, nscur)
184 1536 : call get_horiz_grid_dim_d(hdim1, hdim2)
185 :
186 1536 : ldof => get_restart_decomp(hdim1, hdim2, plev)
187 6144 : call pio_initdecomp(pio_subsystem, pio_double, (/hdim1, hdim2, plev/), ldof, iodesc3d)
188 1536 : deallocate(ldof)
189 :
190 1536 : ldof => get_restart_decomp(hdim1, hdim2, 1)
191 4608 : call pio_initdecomp(pio_subsystem, pio_double, (/hdim1, hdim2/), ldof, iodesc2d)
192 1536 : deallocate(ldof)
193 :
194 : !++bee it does not appear that this time variable is read by read_restart_dynamics. is it read somewhere else?
195 1536 : time = ndcur+(real(nscur,kind=r8))/86400._r8
196 1536 : ierr = pio_put_var(File, timedesc%varid, time)
197 :
198 379392 : do i=1,restartvarcnt
199 377856 : call get_restart_var(File, i, name, timelevels, ndims, vdesc)
200 379392 : if(ndims==2) then
201 3072 : if(use_transfer) then
202 3072 : call pio_write_darray(File, vdesc, iodesc2d, transfer(restartvars(i)%v2d(:,:), mold), ierr)
203 : else
204 0 : call pio_write_darray(File, vdesc, iodesc2d, null, ierr)
205 : end if
206 374784 : else if(ndims==3) then
207 374784 : if(use_transfer) then
208 374784 : call pio_write_darray(File, vdesc, iodesc3d, transfer(restartvars(i)%v3d(:,:,:), mold), ierr)
209 : else
210 0 : call pio_write_darray(File, vdesc, iodesc3d, null, ierr)
211 : end if
212 : end if
213 : end do
214 :
215 1536 : call pio_freedecomp(File, iodesc2d)
216 1536 : call pio_freedecomp(File, iodesc3d)
217 :
218 1536 : end subroutine write_restart_dynamics
219 :
220 : !=========================================================================================
221 :
222 768 : subroutine read_restart_dynamics(File, dyn_in, dyn_out)
223 :
224 : use dyn_comp, only: dyn_init
225 :
226 : ! arguments
227 : type(file_desc_t) :: File
228 : type(dyn_export_t) :: dyn_out
229 : type(dyn_import_t) :: dyn_in
230 :
231 : ! local variables
232 : type(t_fvdycore_state), pointer :: dyn_state
233 : type(t_fvdycore_grid), pointer :: grid
234 :
235 : integer :: beglonxy, endlonxy, beglatxy, endlatxy
236 : integer :: i, s2d, s3d
237 :
238 768 : real(r8), allocatable :: tmp(:)
239 : integer :: dims3d(3)
240 : integer :: ndims
241 : integer :: ierr
242 : character(len=namlen) :: name
243 768 : integer, pointer :: ldof(:)
244 :
245 : integer :: timelevels ! not used in fv
246 : type(var_desc_t), pointer :: vdesc
247 : type(io_desc_t) :: iodesc2d, iodesc3d
248 : !----------------------------------------------------------------------------
249 :
250 : #if ( defined OFFLINE_DYN )
251 : call read_met_restart(File)
252 : #endif
253 :
254 1536 : dyn_state => get_dyn_state()
255 768 : grid => dyn_state%grid
256 :
257 768 : call dyn_init(dyn_in, dyn_out)
258 :
259 768 : beglonxy = grid%ifirstxy
260 768 : endlonxy = grid%ilastxy
261 768 : beglatxy = grid%jfirstxy
262 768 : endlatxy = grid%jlastxy
263 :
264 768 : dims3d(1) = endlonxy - beglonxy + 1
265 768 : dims3d(2) = endlatxy - beglatxy + 1
266 768 : dims3d(3) = plev
267 :
268 768 : s2d = dims3d(1)*dims3d(2)
269 768 : s3d = s2d*dims3d(3)
270 2304 : allocate(tmp(s3d))
271 :
272 768 : call init_restart_varlistr(dyn_in)
273 :
274 768 : ldof => get_restart_decomp(plon, plat, plev)
275 768 : call pio_initdecomp(pio_subsystem, pio_double, (/plon, plat, plev/), ldof, iodesc3d)
276 768 : deallocate(ldof)
277 :
278 768 : ldof => get_restart_decomp(plon, plat, 1)
279 768 : call pio_initdecomp(pio_subsystem, pio_double, (/plon, plat/), ldof, iodesc2d)
280 768 : deallocate(ldof)
281 :
282 189696 : do i = 1, restartvarcnt
283 188928 : call get_restart_var(File, i, name, timelevels, ndims, vdesc)
284 188928 : ierr = PIO_Inq_varid(File, name, vdesc)
285 188928 : call pio_setframe(File, vdesc, int(-1,kind=pio_offset_kind))
286 378624 : if(ndims==2) then
287 1536 : call pio_read_darray(File, vdesc, iodesc2d, tmp(1:s2d), ierr)
288 119808 : restartvars(i)%v2d(:,:) = reshape(tmp(1:s2d), dims3d(1:2))
289 187392 : else if(ndims==3) then
290 187392 : call pio_read_darray(File, restartvars(i)%vdesc, iodesc3d, tmp(1:s3d), ierr)
291 456486912 : restartvars(i)%v3d(:,:,:) = reshape(tmp(1:s3d), dims3d)
292 : end if
293 : end do
294 :
295 768 : deallocate(tmp)
296 768 : call pio_freedecomp(File, iodesc2d)
297 768 : call pio_freedecomp(File, iodesc3d)
298 :
299 1536 : end subroutine read_restart_dynamics
300 :
301 : !=========================================================================================
302 : ! Private
303 : !=========================================================================================
304 :
305 566784 : subroutine set_r_var(name, timelevels, index, v2, v3, v4, v5, v2r4, v3r4, v4r4, v5r4)
306 :
307 : ! arguments
308 : character(len=*), intent(in) :: name
309 : integer, intent(in) :: timelevels, index
310 : real(r8), target, optional :: v2(:,:), v3(:,:,:), v4(:,:,:,:), v5(:,:,:,:,:)
311 : real(r4), target, optional :: v2r4(:,:), v3r4(:,:,:), v4r4(:,:,:,:), v5r4(:,:,:,:,:)
312 : !----------------------------------------------------------------------------
313 :
314 566784 : restartvars(index)%name = name
315 566784 : restartvars(index)%timelevels = timelevels
316 :
317 566784 : if (present(v2)) then
318 4608 : restartvars(index)%ndims = 2
319 4608 : restartvars(index)%v2d => v2
320 562176 : else if (present(v3)) then
321 562176 : restartvars(index)%ndims = 3
322 562176 : restartvars(index)%v3d => v3
323 0 : else if (present(v4)) then
324 0 : restartvars(index)%ndims = 4
325 0 : restartvars(index)%v4d => v4
326 0 : else if (present(v5)) then
327 0 : restartvars(index)%ndims = 5
328 0 : restartvars(index)%v5d => v5
329 0 : else if (present(v2r4)) then
330 0 : restartvars(index)%ndims = 2
331 0 : restartvars(index)%v2dr4 => v2r4
332 0 : else if (present(v3r4)) then
333 0 : restartvars(index)%ndims = 3
334 0 : restartvars(index)%v3dr4 => v3r4
335 0 : else if (present(v4r4)) then
336 0 : restartvars(index)%ndims = 4
337 0 : restartvars(index)%v4dr4 => v4r4
338 0 : else if (present(v5r4)) then
339 0 : restartvars(index)%ndims = 5
340 0 : restartvars(index)%v5dr4 => v5r4
341 : else
342 0 : call endrun('set_r_var: ERROR: bad ndims')
343 : end if
344 :
345 566784 : allocate(restartvars(index)%vdesc)
346 :
347 567552 : end subroutine set_r_var
348 :
349 : !=========================================================================================
350 :
351 1536 : subroutine init_restart_varlistw(dyn_out)
352 :
353 : type(dyn_export_t) :: dyn_out
354 :
355 : integer :: vcnt=1
356 : integer :: i, m
357 : !----------------------------------------------------------------------------
358 :
359 : ! Should only be called once
360 1536 : if (restart_varlistw_initialized) return
361 :
362 1536 : restart_varlistw_initialized = .true.
363 :
364 1536 : call set_r_var('PHIS', 1, vcnt, v2=dyn_out%phis)
365 :
366 1536 : vcnt = vcnt + 1
367 1536 : call set_r_var('U', 1, vcnt, v3=dyn_out%u3s)
368 :
369 1536 : vcnt = vcnt + 1
370 1536 : call set_r_var('V', 1, vcnt, v3=dyn_out%v3s)
371 :
372 1536 : vcnt = vcnt + 1
373 1536 : call set_r_var('DELP', 1, vcnt, v3=dyn_out%delp)
374 :
375 1536 : vcnt = vcnt + 1
376 1536 : call set_r_var('PT', 1, vcnt, v3=dyn_out%pt)
377 :
378 370176 : do m = 1, pcnst
379 368640 : vcnt = vcnt + 1
380 370176 : call set_r_var(cnst_name(m), 1, vcnt, v3=dyn_out%tracer(:,:,:,m))
381 : end do
382 :
383 1536 : vcnt = vcnt + 1
384 1536 : call set_r_var('PS', 1, vcnt, v2=dyn_out%ps )
385 :
386 :
387 1536 : if (vcnt /= restartvarcnt) then
388 0 : write(iulog,*) 'init_restart_varlistw: ERROR: vcnt= ', vcnt, &
389 0 : ' restartvarcnt=', restartvarcnt
390 0 : call endrun('init_restart_varlistw: ERROR: bad restartvarcnt')
391 : end if
392 :
393 : end subroutine init_restart_varlistw
394 :
395 : !=========================================================================================
396 :
397 768 : subroutine init_restart_varlistr(dyn_in)
398 :
399 : type(dyn_import_t) :: dyn_in
400 :
401 : integer :: vcnt=1
402 : integer :: i, m
403 : !----------------------------------------------------------------------------
404 :
405 : ! Should only be called once
406 768 : if (restart_varlistr_initialized) return
407 :
408 768 : restart_varlistr_initialized = .true.
409 :
410 768 : call set_r_var('PHIS', 1, vcnt, v2=dyn_in%phis)
411 :
412 768 : vcnt = vcnt + 1
413 768 : call set_r_var('U', 1, vcnt, v3=dyn_in%u3s)
414 :
415 768 : vcnt = vcnt + 1
416 768 : call set_r_var('V', 1, vcnt, v3=dyn_in%v3s)
417 :
418 768 : vcnt = vcnt + 1
419 768 : call set_r_var('DELP', 1, vcnt, v3=dyn_in%delp)
420 :
421 768 : vcnt = vcnt + 1
422 768 : call set_r_var('PT', 1, vcnt, v3=dyn_in%pt)
423 :
424 185088 : do m = 1, pcnst
425 184320 : vcnt = vcnt + 1
426 185088 : call set_r_var(cnst_name(m), 1, vcnt, v3=dyn_in%tracer(:,:,:,m))
427 : end do
428 :
429 768 : vcnt = vcnt + 1
430 768 : call set_r_var('PS', 1, vcnt, v2=dyn_in%ps)
431 :
432 768 : if (vcnt /= restartvarcnt) then
433 0 : write(iulog,*) 'init_restart_varlistr: ERROR: vcnt= ', vcnt, &
434 0 : ' restartvarcnt=',restartvarcnt
435 0 : call endrun('init_restart_varlistr: ERROR: bad restartvarcnt')
436 : end if
437 :
438 : end subroutine init_restart_varlistr
439 :
440 : !=========================================================================================
441 :
442 4608 : function get_restart_decomp(hdim1, hdim2, nlev) result(ldof)
443 :
444 : ! Get the integer mapping of a variable in the dynamics decomp in memory.
445 : ! The canonical ordering is as on the file. A 0 value indicates that the
446 : ! variable is not on the file (eg halo or boundary values)
447 :
448 : ! arguments
449 : integer, intent(in) :: hdim1, hdim2, nlev
450 : integer, pointer :: ldof(:)
451 :
452 : ! local variables
453 : integer :: i, k, j
454 : integer :: lcnt
455 : integer :: beglatxy, beglonxy, endlatxy, endlonxy
456 : !----------------------------------------------------------------------------
457 :
458 4608 : beglonxy = get_dyn_grid_parm('beglonxy')
459 4608 : endlonxy = get_dyn_grid_parm('endlonxy')
460 4608 : beglatxy = get_dyn_grid_parm('beglatxy')
461 4608 : endlatxy = get_dyn_grid_parm('endlatxy')
462 :
463 4608 : lcnt = (endlatxy-beglatxy+1)*nlev*(endlonxy-beglonxy+1)
464 13824 : allocate(ldof(lcnt))
465 5478912 : ldof(:) = 0
466 :
467 : lcnt = 0
468 80640 : do k = 1, nlev
469 308736 : do j = beglatxy, endlatxy
470 5778432 : do i = beglonxy, endlonxy
471 5474304 : lcnt = lcnt + 1
472 5702400 : ldof(lcnt) = i + (j-(plat-hdim2+1))*hdim1+(k-1)*hdim1*hdim2
473 : end do
474 : end do
475 : end do
476 :
477 4608 : end function get_restart_decomp
478 :
479 : !=========================================================================================
480 :
481 944640 : subroutine get_restart_var(File, i, name, timelevels, ndims, vdesc)
482 :
483 : type(file_desc_t) :: File
484 : integer, intent(in) :: i
485 : character(len=namlen), intent(out) :: name
486 : integer, intent(out) :: timelevels
487 : integer, intent(out) :: ndims
488 : type(var_desc_t), pointer :: vdesc
489 :
490 944640 : name = restartvars(i)%name
491 944640 : timelevels = restartvars(i)%timelevels
492 944640 : ndims = restartvars(i)%ndims
493 944640 : if (.not.associated(restartvars(i)%vdesc)) then
494 0 : allocate(restartvars(i)%vdesc)
495 : end if
496 944640 : vdesc => restartvars(i)%vdesc
497 :
498 944640 : end subroutine get_restart_var
499 :
500 : !=========================================================================================
501 :
502 0 : end module restart_dynamics
|
