!#define DEBUG 1 module physics_buffer !----------------------------------------------------------------------- ! ! Purpose: ! Buffer managment for persistant variables ! ! Author: J. Edwards ! ! This file is used with genf90.pl to generate buffer.F90 ! !----------------------------------------------------------------------- use shr_kind_mod, only: r8 => shr_kind_r8, r4=> shr_kind_r4, i4=> shr_kind_i4 use ppgrid, only: pcols, begchunk, endchunk, psubcols use cam_logfile, only: iulog use pio, only: var_desc_t use dyn_grid, only: ptimelevels use cam_abortutils, only: endrun use buffer, only: buffer_field_allocate, buffer_field_deallocate, buffer_get_field_ptr, buffer_set_field, & dtype_i4, dtype_r4, dtype_r8, buffer_field_default_type, buffer_field_is_alloc implicit none private ! ngrid_types is a private parameter denoting how many types of a field ! (e.g., grid, subcol) a pbuf can hold integer, parameter :: ngrid_types = 2 ! --col_type parameters -- see pbuf_add_field and pbuf_get_field ! -- These indices should start at 1 and increment by 1 so that they ! -- may be used in a loop from 1 to ngrid_types integer, parameter, public :: col_type_grid = 1 integer, parameter, public :: col_type_subcol = 2 ! ! PBUF field suffix strings for different grid types for restart files ! NB: There should be ngrid_types entries ! character(len=5), parameter :: field_grid_suff(ngrid_types) = (/ " ", "_SCOL" /) integer :: old_time_idx = 1 ! The API has strings 'GLOBAL' and 'PHYSPKG' which correspond to these ! integer constants if global_allocate_all is false fields allocated with ! PHYSPKG persistence are deallocated at the end of each physics time step ! and reallocated at the beginning of the next. ! If global_allocate_all is true then all fields are allocated at model ! start and persist until model completion. integer, parameter :: persistence_global = 1 integer, parameter :: persistence_physpkg = 2 logical :: global_allocate_all = .true. ! ! physics_buffer_hdr carries the description info for each buffer, only one header ! is allocated for each field and each chunk of the field points to it. ! It is carried as a linkedlist for initialization only. ! type physics_buffer_hdr character(len=16) :: name = '' integer :: dtype = -1 integer :: persistence = -1 integer :: dimsizes(6,ngrid_types) = 0 type(var_desc_t) :: vardesc(ngrid_types) ! var id for restart files logical :: is_copy(ngrid_types) = .false. logical :: added = .false. type(physics_buffer_hdr), pointer :: nexthdr => null() end type physics_buffer_hdr ! ! The default type for a buffer field is buffer_field_double (r8) since that is the ! type most often used in the model. The F90 transfer function is used to move ! data of other types into and out of the pbuf2d ! type physics_buffer_desc private integer :: lchnk type(physics_buffer_hdr), pointer :: hdr => null() type(buffer_field_default_type) :: bfg(ngrid_types) end type physics_buffer_desc interface pbuf_get_field ! TYPE int,double,real ! DIMS 1,2,3,4,5 module procedure get_pbuf1d_field_by_index_{DIMS}d_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 module procedure get_pbuf2d_field_by_index_{DIMS}d_{TYPE} end interface interface pbuf_get_field_restart ! TYPE int,double,real module procedure get_pbuf2d_field_restart_{TYPE} end interface interface pbuf_set_field ! TYPE int,double,real ! DIMS 1,2,3,4,5 module procedure set_pbuf2d_field_by_index_{DIMS}d_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 module procedure set_pbuf1d_field_by_index_{DIMS}d_{TYPE} ! TYPE int,double,real module procedure set_pbuf1d_field_const_by_index_{TYPE} ! TYPE int,double,real module procedure set_pbuf2d_field_const_by_index_{TYPE} end interface interface pbuf_add_field ! TYPE int,double,real module procedure pbuf_add_field_{TYPE} end interface public :: pbuf_initialize, & pbuf_readnl, &! read namelist options pbuf_init_time, &! Initialize dyn_time_lvls pbuf_old_tim_idx, &! return the index for the oldest time pbuf_update_tim_idx, &! update the index for the oldest time pbuf_col_type_index, & pbuf_get_field_name, & pbuf_get_field, & pbuf_add_field, & pbuf_register_subcol, & physics_buffer_desc, & pbuf_get_index, & pbuf_get_chunk, & pbuf_allocate, & pbuf_deallocate, & pbuf_set_field, & pbuf_init_restart, & pbuf_write_restart, & pbuf_read_restart, & dtype_r8, dtype_r4, dtype_i4 ! For help debugging code public :: pbuf_dump_pbuf integer, public :: dyn_time_lvls ! number of time levels in physics buffer (dycore dependent) ! ! Currentpbufflds is incremented in calls to pbuf_add_field and determines the size ! of the allocated pbuf2d ! integer :: currentpbufflds=0 type(physics_buffer_hdr), pointer :: hdrbuffertop => null() ! ! Insures that we do not attempt to allocate physics_buffer more than once ! logical :: buffer_initialized =.false. ! ! private pio descriptor for time ! type(var_desc_t) :: timeidx_desc !=============================================================================== CONTAINS !=============================================================================== subroutine pbuf_readnl(nlfile) use namelist_utils, only: find_group_name use units, only: getunit, freeunit use spmd_utils, only: masterproc, mpicom, mstrid=>masterprocid, mpi_logical character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input ! Local variables integer :: unitn, ierr character(len=*), parameter :: sub = 'pbuf_readnl' logical :: pbuf_global_allocate namelist /pbuf_nl/ pbuf_global_allocate !----------------------------------------------------------------------------- pbuf_global_allocate = global_allocate_all ! Read namelist if (masterproc) then unitn = getunit() open( unitn, file=trim(nlfile), status='old' ) call find_group_name(unitn, 'pbuf_nl', status=ierr) if (ierr == 0) then read(unitn, pbuf_nl, iostat=ierr) if (ierr /= 0) then call endrun(sub//': FATAL: reading namelist') end if end if close(unitn) call freeunit(unitn) end if call mpi_bcast(pbuf_global_allocate, 1, mpi_logical, mstrid, mpicom, ierr) if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: pbuf_global_allocate") global_allocate_all = pbuf_global_allocate if (masterproc) then write(iulog,*) 'physics buffer options:' write(iulog,*) ' pbuf_global_allocate =', global_allocate_all end if end subroutine pbuf_readnl !=============================================================================== ! ! Initialize dyn_time_lvls ! subroutine pbuf_init_time() dyn_time_lvls = ptimelevels - 1 end subroutine pbuf_init_time ! ! Return index of oldest time sample in the physics buffer. ! function pbuf_old_tim_idx() integer :: pbuf_old_tim_idx pbuf_old_tim_idx = old_time_idx end function pbuf_old_tim_idx ! ! Update index of old time sample in the physics buffer. ! subroutine pbuf_update_tim_idx() old_time_idx = mod(old_time_idx, dyn_time_lvls) + 1 end subroutine pbuf_update_tim_idx ! ! pbuf_col_type_index returns an index for use with pbuf calls ! ! * col_type: is set to col_type_grid (if use_subcol=.false.) or ! col_type_subcol (if (use_subcol=.true.). ! subroutine pbuf_col_type_index(use_subcol, col_type) logical, intent(in) :: use_subcol integer, intent(out) :: col_type if (use_subcol) then col_type = col_type_subcol else col_type = col_type_grid end if end subroutine pbuf_col_type_index ! ! Return a pointer to the current chunks physics_buffer. ! function pbuf_get_field_name(index) integer, intent(in) :: index character(len=16) :: pbuf_get_field_name integer :: i type(physics_buffer_hdr), pointer :: hdrbuffer hdrbuffer => hdrbuffertop do i=2,index hdrbuffer=>hdrbuffer%nexthdr end do pbuf_get_field_name = hdrbuffer%name end function pbuf_get_field_name function pbuf_get_chunk(pbuf2d, lchnk) type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: lchnk type(physics_buffer_desc), pointer :: pbuf_get_chunk(:) pbuf_get_chunk => pbuf2d(:,lchnk) end function pbuf_get_chunk ! ! Return .true. iff pbuf has an allocated grid field ! logical function pbuf_field_has_gridcols(pbuf, index) result(rval) type(physics_buffer_desc), pointer :: pbuf(:) integer, intent(in) :: index ! If the field is a copy, return .false. even if it is allocated if (pbuf(index)%hdr%is_copy(col_type_grid)) then rval = .false. else rval = buffer_field_is_alloc(pbuf(index)%bfg(col_type_grid)) end if end function pbuf_field_has_gridcols ! ! Return .true. iff pbuf has an allocated subcolumn field ! logical function pbuf_field_has_subcols(pbuf, index) result(rval) type(physics_buffer_desc), pointer :: pbuf(:) integer, intent(in) :: index ! If the field is a copy, return .false. even if it is allocated if (pbuf(index)%hdr%is_copy(col_type_subcol)) then rval = .false. else rval = buffer_field_is_alloc(pbuf(index)%bfg(col_type_subcol)) end if end function pbuf_field_has_subcols ! ! Return .true. iff pbuf has an allocated col_type column field ! logical function pbuf_field_has_col_type(pbuf, index, col_type) result(rval) type(physics_buffer_desc), pointer :: pbuf(:) integer, intent(in) :: index integer, intent(in) :: col_type if (col_type == col_type_grid) then rval = pbuf_field_has_gridcols(pbuf, index) else if(col_type == col_type_subcol) then rval = pbuf_field_has_subcols(pbuf, index) else call endrun('pbuf_field_has_col_type: Invalid col_type') end if end function pbuf_field_has_col_type ! ! Initialize the buffer, should be called after all pbuf_add_field calls ! have been completed and should only be called once in a run ! subroutine pbuf_initialize(pbuf2d) type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer :: i, c type(physics_buffer_hdr), pointer :: hdrbuffer ! ! Allocate memory ! if(buffer_initialized) return ! Allocate at least 1 to avoid unallocated error in ideal physics allocate(pbuf2d(max(1,currentpbufflds),begchunk:endchunk)) if(currentpbufflds<1) return do c = begchunk, endchunk hdrbuffer => hdrbuffertop do i = 1, currentpbufflds if(.not. hdrbuffer%added) then call endrun('pbuf_initialize: pbuf, '//trim(hdrbuffer%name)//', not added') end if pbuf2d(i,c)%lchnk = c pbuf2d(i,c)%hdr => hdrbuffer hdrbuffer => hdrbuffer%nexthdr end do end do buffer_initialized=.true. call pbuf_allocate(pbuf2d, 'global') #ifdef DEBUG call pbuf2d_print(pbuf2d) #endif return end subroutine pbuf_initialize subroutine pbuf_allocate(pbuf2d, persistence) type(physics_buffer_desc), pointer :: pbuf2d(:,:) character(len=*), intent(in) :: persistence integer :: i logical :: allocate_all ! allocate_all is used to force allocation of all fields at same time as allocation ! for global scope. allocate_all = .false. if ( global_allocate_all ) then if ( persistence == 'global' ) then allocate_all = .true. else return end if end if if(allocate_all) then do i=1,currentpbufflds select case(pbuf2d(i,begchunk)%hdr%dtype) case(TYPEDOUBLE) call pbuf_allocate_field_double(pbuf2d, i, dtype_r8) case(TYPEREAL) call pbuf_allocate_field_real(pbuf2d, i, dtype_r4) ! case(i8) ! call pbuf_allocate_field_long(pbuf2d, i) case(TYPEINT) call pbuf_allocate_field_int(pbuf2d, i, dtype_i4) end select end do else if(persistence .eq. 'physpkg') then do i=1,currentpbufflds if(pbuf2d(i,begchunk)%hdr%persistence==persistence_physpkg) then select case(pbuf2d(i,begchunk)%hdr%dtype) case(TYPEDOUBLE) call pbuf_allocate_field_double(pbuf2d, i, dtype_r8) case(TYPEREAL) call pbuf_allocate_field_real(pbuf2d, i, dtype_r4) ! case(i8) ! call pbuf_allocate_field_long(pbuf2d, i) case(TYPEINT) call pbuf_allocate_field_int(pbuf2d, i, dtype_i4) end select end if end do else if(persistence .eq. 'global') then do i=1,currentpbufflds if(pbuf2d(i,begchunk)%hdr%persistence==persistence_global) then select case(pbuf2d(i,begchunk)%hdr%dtype) case(TYPEDOUBLE) call pbuf_allocate_field_double(pbuf2d, i, dtype_r8) case(TYPEREAL) call pbuf_allocate_field_real(pbuf2d, i, dtype_r4) ! case(i8) ! call pbuf_allocate_field_long(pbuf2d, i) case(TYPEINT) call pbuf_allocate_field_int(pbuf2d, i, dtype_i4) end select end if end do end if end subroutine pbuf_allocate ! TYPE int,double,real subroutine pbuf_allocate_field_{TYPE}(pbuf2d, index, dtype) {VTYPE}, intent(in) :: dtype type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: index integer, pointer :: dimsizes(:) integer :: c, i logical :: is_copy do i = 1, ngrid_types ! Note - dimsizes(:)=0 is special case to indicate "do not allocate" and is not fatal ! Since this is called by a single thread, setting the dimsizes pointer to first chunk is okay dimsizes => pbuf2d(index,begchunk)%hdr%dimsizes(:,i) ! Note - We may have dimsizes /= 0 but still don't allocate if copy is_copy = pbuf2d(index,begchunk)%hdr%is_copy(i) if (any(dimsizes(:) < 0)) & call endrun('pbuf_allocate_field: dimsizes must be greater than 0 for pbuf field '& //trim(pbuf2d(index,begchunk)%hdr%name)) if (all(dimsizes(:) /= 0) .and. (.not. is_copy)) then do c = begchunk, endchunk call buffer_field_allocate(pbuf2d(index,c)%bfg(i), dimsizes, dtype) end do end if end do end subroutine pbuf_allocate_field_{TYPE} pure logical function pbuf_do_deallocate(hdr, persistence, col_type) result(rval) type(physics_buffer_hdr), pointer :: hdr character(len=*), intent(in) :: persistence integer, intent(in) :: col_type if (persistence .eq. 'physpkg') then if (hdr%is_copy(col_type)) then ! If this is a copied field, it is always deallocated rval = .true. else if (global_allocate_all) then ! This is an all-global run so no deallocation rval = .false. else if(hdr%persistence == persistence_physpkg) then ! This pbuf is a physpkg type, deallocate rval = .true. else ! This pbuf is global, do not deallocate rval = .false. end if else ! We are doing a global deallocate, everything must go! rval = .true. end if end function pbuf_do_deallocate subroutine pbuf_deallocate(pbuf2d, persistence) type(physics_buffer_desc), pointer :: pbuf2d(:,:) character(len=*) :: persistence integer :: i, j, c do i = 1, currentpbufflds do j = 1, ngrid_types if (pbuf_do_deallocate(pbuf2d(i,begchunk)%hdr, persistence, j)) then if (buffer_field_is_alloc(pbuf2d(i,begchunk)%bfg(j))) then select case(pbuf2d(i,begchunk)%hdr%dtype) case(TYPEDOUBLE) do c = begchunk,endchunk call buffer_field_deallocate(pbuf2d(i,c)%bfg(j), dtype_r8) end do case(TYPEREAL) do c = begchunk,endchunk call buffer_field_deallocate(pbuf2d(i,c)%bfg(j), dtype_r4) end do case(TYPEINT) do c = begchunk,endchunk call buffer_field_deallocate(pbuf2d(i,c)%bfg(j), dtype_i4) end do end select end if end if end do ! ngrid_types end do ! currentpbufflds end subroutine pbuf_deallocate ! Find a pbuf header pointer based on the name input ! Automatically tacks on an extra header if is not found. subroutine find_pbuf_header(name, index, bufptr) character(len=*), intent(in) :: name type(physics_buffer_hdr), pointer :: bufptr integer, intent(out) :: index ! Local Variables logical :: buf_found if(.not. associated(hdrbuffertop)) then ! This is the very first pbuf, allocate and set allocate(hdrbuffertop) currentpbufflds = 1 hdrbuffertop%name = name end if bufptr=>hdrbuffertop buf_found = .false. index = 1 do while(.not. buf_found) if(trim(bufptr%name) == trim(name)) then buf_found = .true. else if (associated(bufptr%nexthdr)) then bufptr=>bufptr%nexthdr index = index + 1 else ! We ran off the end of the buffers, make a new one for ! Sanity check, we should have checked exactly currentpbufflds if (index /= currentpbufflds) then call endrun("find_pbuf_header: currentpbufflds indexing off") end if currentpbufflds = currentpbufflds + 1 index = currentpbufflds allocate(bufptr%nexthdr) bufptr=>bufptr%nexthdr bufptr%name = trim(name) buf_found = .true. end if end if end do end subroutine find_pbuf_header ! ! Register a field in the pbuf ! This should be called from physics register routines. ! persistence must be 'global' or 'physpkg' ! dtype can be any of r8, r4, i4 as defined in shr_kinds_mod.F90 ! col_type is either col_type_grid or col_type_subcol. ! If no col_type, then grid field is defined (i.e., dimsizes set) subroutine pbuf_register_field_int(name, pname, index, persistence, & dtype, dimsizes, col_type, pbuf_add) ! Dummy Arguments character(len=*), intent(in) :: name character(len=*), intent(in) :: pname ! name of calling parameterization integer, intent(out) :: index character(len=*), optional, intent(in) :: persistence integer, optional, intent(in) :: dtype ! used to differentiate specific calls integer, optional, intent(in) :: dimsizes(:) ! dimension sizes of grid field integer, optional, intent(in) :: col_type logical, optional, intent(in) :: pbuf_add ! Local Variables type(physics_buffer_hdr), pointer :: bufptr integer :: dimcnt, col_type_use character(len=128) :: errmsg if(buffer_initialized) then call endrun('Attempt to register pbuf field after buffer initialized') end if call find_pbuf_header(name, index, bufptr) if (present(persistence)) then if(persistence .eq. "global") then bufptr%persistence = persistence_global else bufptr%persistence = persistence_physpkg end if end if if (present(dtype)) then bufptr%dtype = dtype end if if (present(dimsizes)) then ! Normally, we only set buffer dimsizes if dimsizes is passed dimcnt = size(dimsizes) if (present(col_type)) then col_type_use = col_type else col_type_use = col_type_grid end if ! Only allow up to 5 dimensions. #6 is reserved for subcolumn umpacking ! and #7 is for the physics chunk index. if (dimcnt > 5) then call endrun('pbuf_register_field: Attempt to exceed maximum of 5 dimensions for '//trim(name)) end if ! Assign dimensions dependent on col_types input ! Note that dimensions are initialized to zero and set if being used if (col_type_use == col_type_grid) then ! grid is requested bufptr%dimsizes(:,col_type_grid) = 1 bufptr%dimsizes(1:dimcnt,col_type_grid) = dimsizes ! If someone previously registered the subcol, reset those dims if (bufptr%dimsizes(1,col_type_subcol) == pcols*psubcols) then bufptr%dimsizes(2:dimcnt,col_type_subcol) = dimsizes(2:dimcnt) end if else if (col_type_use == col_type_subcol) then bufptr%dimsizes(:,col_type_subcol) = 1 bufptr%dimsizes(1,col_type_subcol) = pcols*psubcols ! This case should only be used for a pbuf_add_field call adding a subcolumn-only field if (dimcnt > 1) then bufptr%dimsizes(2:dimcnt,col_type_subcol) = dimsizes(2:dimcnt) else if (bufptr%dimsizes(1,col_type_grid) > 0) then ! If grid field previously registered or added, update dims bufptr%dimsizes(2:,col_type_subcol) = bufptr%dimsizes(2:,col_type_grid) end if end if end if if (present(pbuf_add)) then if (pbuf_add) then ! An add request has been made but we have to make sure we have all info if (all(bufptr%dimsizes(:,:) == 0)) then write(errmsg, *) 'pbuf_add_field: trying to add field with no dimensions' call endrun(errmsg) end if if (bufptr%dtype < 0) then write(errmsg, *) 'pbuf_add_field: trying to add field with no type' call endrun(errmsg) end if if (bufptr%persistence < 0) then write(errmsg, *) 'pbuf_add_field: trying to add field with bad persistence' call endrun(errmsg) end if bufptr%added = .true. end if end if end subroutine pbuf_register_field_int ! ! Add a field to the pbuf, this should be called from physics register routines ! name is required to be unique, ! persistence must be 'global' or 'physpkg' ! dtype can be any of r8, r4, i4 as defined in shr_kinds_mod.F90 ! If present, col_type must be either col_type_grid or col_type_subcol ! TYPE int,double,real subroutine pbuf_add_field_{TYPE}(name, persistence, dtype, dimsizes, index, col_type) character(len=*), intent(in) :: name, persistence {VTYPE}, intent(in) :: dtype ! used only to differentiate specific calls integer, intent(in) :: dimsizes(:) ! dimension sizes of grid field integer, intent(in), optional :: col_type integer, intent(out) :: index ! Local Variables integer :: col_type_use if(buffer_initialized) then call endrun('Attempt to add pbuf field after buffer initialized') end if if (present(col_type)) then col_type_use = col_type else col_type_use = col_type_grid end if call pbuf_register_field_int(trim(name), '', index, & persistence=persistence, dtype={ITYPE}, & dimsizes=dimsizes, col_type=col_type_use, pbuf_add=.true.) end subroutine pbuf_add_field_{TYPE} subroutine pbuf_register_subcol(name, pname, index) use subcol_utils, only: is_subcol_on ! Dummy Arguments character(len=*), intent(in) :: name character(len=*), intent(in) :: pname ! name of calling parameterization integer, intent(out) :: index ! Local variables integer :: dimsizes(1) ! You really should not call this routine if subcolumns are not on if (.not. is_subcol_on()) then call endrun('pbuf_register_subcol: subcolumns are not active') end if ! Create and pass dimsizes so that subcolums are registered dimsizes(1) = pcols * psubcols call pbuf_register_field_int(trim(name), trim(pname), index, & dimsizes=dimsizes, col_type=col_type_subcol) end subroutine pbuf_register_subcol subroutine pbuf2d_print(pbuf2d) type(physics_buffer_desc), pointer :: pbuf2d(:,:) call pbuf1d_print(pbuf_get_chunk(pbuf2d,begchunk)) end subroutine pbuf2d_print subroutine pbuf1d_print(pbuf) type(physics_buffer_desc), pointer :: pbuf(:) integer :: i type(physics_buffer_desc), pointer :: pbufPtr print *,__FILE__,__LINE__,currentpbufflds,size(pbuf) do i=1,currentpbufflds pbufPtr => pbuf(i) print *,__FILE__,__LINE__,i,trim(pbufPtr%hdr%name),pbufPtr%hdr%dtype,pbufPtr%hdr%persistence,pbufPtr%hdr%dimsizes end do end subroutine pbuf1d_print ! ! Given a pbuf field name return an integer index to the field. ! This index can be used to retrieve the field and is much faster ! than using the name in most cases ! function pbuf_get_index(name, errcode) result(index) character(len=*), intent(in) :: name integer, intent(inout), optional :: errcode integer :: index integer :: i type(physics_buffer_hdr), pointer :: bufptr bufptr=>hdrbuffertop index = -1 do i=1,currentpbufflds if(bufptr%name .eq. name) then index=i exit end if bufptr=>bufptr%nexthdr end do if (present(errcode)) then errcode = index else if(index<0) then call endrun('Attempt to find undefined name in pbuf '//trim(name)) end if end function pbuf_get_index !========================================================================================= ! ! Given a pbuf2d chunk and an index return a pointer to a field chunk ! ! ! TYPE int,double,real ! DIMS 1,2,3,4,5 subroutine get_pbuf1d_field_by_index_{DIMS}d_{TYPE}(pbuf, index, field, start,kount, col_type, copy_if_needed, errcode) ! Get the data based on the col_type which is specified. If no col_type, then grid field is returned use subcol_utils, only: subcol_field_copy type(physics_buffer_desc), pointer:: pbuf(:) integer, intent(in) :: index {VTYPE}, pointer :: field{DIMSTR} integer, intent(in), optional :: start(:),kount(:) integer, intent(in), optional :: col_type logical, intent(in), optional :: copy_if_needed integer, intent(out), optional :: errcode integer, pointer :: dimsizes(:) {VTYPE}, pointer :: gfield{DIMSTR} {VTYPE}, pointer :: sfield{DIMSTR} integer :: col_type_use integer, allocatable :: kount_grid(:) ! For use in copy_if_needed logical :: subset ! Copy the generic type to one compatable with the request if(index<1 .or. index>size(pbuf)) then print *,__FILE__,__LINE__,index call endrun('index out of range') end if ! Default col_type is grid columns if (present(col_type)) then col_type_use = col_type else col_type_use = col_type_grid end if ! If there is an errcode, start with an OK status (zero) if (present(errcode)) then errcode = 0 end if ! Check whether subset of data requested (default is false) subset = .false. if (present(start) .and. present(kount)) subset = .true. ! Check for ill-formed request if ( (present(start) .and. .not. present(kount)) .or. & (.not. present(start) .and. present(kount)) ) then call endrun('pbuf_get_field: Both start and kount must be present for '//trim(pbuf(index)%hdr%name)) end if ! See if we need to copy the grid field to subcolumns if (present(copy_if_needed)) then if (copy_if_needed) then if (col_type_use == col_type_subcol) then ! If a subcolumn field buffer does not exist, allocate one. ! Even if start and kount are being passed, allocate and copy the ! full grid field buffer so that a future access will succeed. if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_use))) then dimsizes => pbuf(index)%hdr%dimsizes(:,col_type_use) dimsizes(2:) = pbuf(index)%hdr%dimsizes(2:,col_type_grid) dimsizes(1) = pcols * psubcols select case(pbuf(index)%hdr%dtype) case(TYPEDOUBLE) call buffer_field_allocate(pbuf(index)%bfg(col_type_subcol), & dimsizes, dtype_r8) case(TYPEREAL) call buffer_field_allocate(pbuf(index)%bfg(col_type_subcol), & dimsizes, dtype_r4) case(TYPEINT) call buffer_field_allocate(pbuf(index)%bfg(col_type_subcol), & dimsizes, dtype_i4) end select ! Set copy only if we did the allocation after init time pbuf(index)%hdr%is_copy(col_type_subcol) = .true. end if else call endrun('pbuf_get_field: copy_if_needed only supported for subcolumns') end if if (pbuf(index)%hdr%is_copy(col_type_subcol)) then ! Only do the copy if we did the alloc (i.e., set the is_copy flag) ! Only copy the portion we are going to hand back (i.e., start, kount) ! Chances are that kount(1) = pcols*psubcols because we are looking ! for a subcolumn field (or we wouldn't be here). Now, ! subcol_field_copy requires kount(1) = pcols for the input and ! therefore, kount(1) = pcols*psubcols for the output. Check and ! make it work if (subset) then ! Create kount array for grid field ! Use input subcol kount array, replacing the first dimension with pcols if (size(kount) > size(pbuf(index)%hdr%dimsizes(:,col_type_subcol))) then call endrun('pbuf_get_field: kount input has too many dimensions') end if if (kount(1) /= pcols * psubcols) then call endrun('pbuf_get_field: kount(1) must be pcols*psubcols when using copy_if_needed=.true.') endif allocate(kount_grid(size(kount))) kount_grid(2:) = kount(2:) kount_grid(1) = pcols ! Don't need to create start array for grid field as start array for subcolumn field is identical if (size(start) > size(pbuf(index)%hdr%dimsizes(:,col_type_subcol))) then call endrun('pbuf_get_field: start input has too many dimensions') end if if (start(1) /= 1) then call endrun('pbuf_get_field: start(1) must be 1 when using copy_if_needed=.true.') end if ! Get the grid field call buffer_get_field_ptr(pbuf(index)%bfg(col_type_grid), & gfield, start, kount_grid) deallocate(kount_grid) else ! Get the grid field call buffer_get_field_ptr(pbuf(index)%bfg(col_type_grid), & gfield) end if ! Get the subcol field pointer (note optional start/kount retain their status in this call) call buffer_get_field_ptr(pbuf(index)%bfg(col_type_subcol), & sfield, start, kount) call subcol_field_copy(gfield, pbuf(index)%lchnk, sfield) end if end if end if ! Copy or not, retrieve the requested field pointer if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_use))) then if (present(errcode)) then errcode = -1 else if (col_type_use == col_type_grid) then call endrun('pbuf_get_field: probably missing a pbuf_add_field call. field not allocated for '& //trim(pbuf(index)%hdr%name)) else if (col_type_use == col_type_subcol) then call endrun('pbuf_get_field: probably missing a pbuf_register_subcol. field not allocated for '& //trim(pbuf(index)%hdr%name)) else call endrun('pbuf_get_field: field not allocated for '//trim(pbuf(index)%hdr%name)) end if end if else ! Get the field pointer (note optional start/kount retain their status in this call) call buffer_get_field_ptr(pbuf(index)%bfg(col_type_use),field,start,kount ) end if end subroutine get_pbuf1d_field_by_index_{DIMS}d_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 subroutine get_pbuf2d_field_by_index_{DIMS}d_{TYPE}(pbuf2d, lchnk, index, field, start, kount, col_type, errcode) ! Get the data based on the col_type which is specified. If no col_type, then grid field is returned type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: lchnk integer, intent(in) :: index integer, intent(in), optional :: start(:),kount(:) integer, intent(in), optional :: col_type integer, intent(out), optional :: errcode {VTYPE}, pointer :: field{DIMSTR} ! Check for ill-formed request if ( (present(start) .and. .not. present(kount)) .or.& (.not. present(start) .and. present(kount)) ) then call endrun('pbuf_get_field: Both start and kount must be present for '//trim(pbuf2d(index,begchunk)%hdr%name)) end if call pbuf_get_field(pbuf_get_chunk(pbuf2d,lchnk), index, field, & start=start, kount=kount, col_type=col_type, errcode=errcode) end subroutine get_pbuf2d_field_by_index_{DIMS}d_{TYPE} ! TYPE int,double,real subroutine get_pbuf2d_field_restart_{TYPE}(pbuf2d, lchnk, index, field, mdimsize, col_type) use subcol_utils, only: subcol_unpack ! NB: This routine is really only useful for write_restart_field ! Get the data based on the col_type which is specified. ! If no col_type, then grid field is assumed ! For grid field, reference buffer and copy into chunk (field) ! If col_type is a subcol, unpack subcolumn data ! Field is then reshaped into (/pcols, psubcols*mdims/). ! Dummy variables type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: lchnk integer, intent(in) :: index integer, intent(in) :: mdimsize integer, intent(in) :: col_type {VTYPE}, intent(inout) :: field(:,:) ! Local variables {VTYPE}, allocatable :: exp_fld(:,:,:,:,:,:) {VTYPE}, pointer :: buf(:,:,:,:,:) {VTYPE} :: fillvalue character(len=128) :: errmsg integer, pointer :: dimsizes(:) #if ({ITYPE}==TYPEREAL) fillvalue = 0._r4 #elif ({ITYPE}==TYPEDOUBLE) fillvalue = 0._r8 #else fillvalue = 0 #endif if (col_type == col_type_grid) then call pbuf_get_field(pbuf2d, lchnk, index, buf, col_type=col_type) field(:,:) = reshape(buf, (/pcols, mdimsize/)) else if (col_type == col_type_subcol) then ! Don't initialize field, unpack will fill in unused slots call pbuf_get_field(pbuf2d, lchnk, index, buf, col_type=col_type) dimsizes => pbuf2d(index, lchnk)%hdr%dimsizes(:,col_type_subcol) allocate(exp_fld(pcols, psubcols, dimsizes(2), dimsizes(3), dimsizes(4), dimsizes(5))) ! unpack the subcolumns into their own dimension call subcol_unpack(lchnk, buf, exp_fld, fillvalue) ! Reshape back to pcols for outputting field(:,:) = reshape(exp_fld, (/pcols, mdimsize/)) deallocate(exp_fld) else write(errmsg, *) "get_pbuf2d_field_restart_{TYPE}: Bad col_type:",col_type call endrun(errmsg) end if end subroutine get_pbuf2d_field_restart_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 subroutine set_pbuf2d_field_const_by_index_{TYPE}(pbuf2d,index,const, col_type) ! Set the field specified by the col_type type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: index {VTYPE},intent(in) :: const integer,intent(in) ,optional :: col_type integer :: c do c=begchunk,endchunk if(present(col_type)) then call set_pbuf1d_field_const_by_index_{TYPE}(pbuf_get_chunk(pbuf2d,c),index,const, col_type=col_type) else call set_pbuf1d_field_const_by_index_{TYPE}(pbuf_get_chunk(pbuf2d,c),index,const) end if end do end subroutine set_pbuf2d_field_const_by_index_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 subroutine set_pbuf1d_field_const_by_index_{TYPE}(pbuf,index,const,start,kount, col_type) ! Set the field(s) specified by the col_type type(physics_buffer_desc), pointer :: pbuf(:) integer, intent(in) :: index {VTYPE},intent(in) :: const integer, intent(in), optional :: start(:),kount(:) integer, intent(in), optional :: col_type integer :: col_type_use logical :: subset ! Default col_type is grid if (present(col_type)) then col_type_use = col_type else col_type_use = col_type_grid end if ! Check whether subset of data requested (default is false) subset = .false. if (present(start) .and. present(kount)) subset = .true. ! Check for ill-formed request if ( (present(start) .and. .not. present(kount)) .or.& (.not. present(start) .and. present(kount)) ) then call endrun('pbuf_set_field: Both start and kount must be present for '//trim(pbuf(index)%hdr%name)) end if ! Set the appropriate grid or sub-column field. Check that the fields have been allocated. if(subset) then if (col_type_use == col_type_subcol) then ! Set sub-column field if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_subcol))) & call endrun('pbuf_set_field: sub-column field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_subcol),const,start,kount) else if (col_type_use == col_type_grid) then ! Set grid field if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_grid))) & call endrun('pbuf_set_field: grid field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_grid),const,start,kount) else call endrun('pbuf_set_field: Trying to set '//trim(pbuf(index)%hdr%name)//& ' but col_type is neither col_type_grid nor col_type_subcol') end if else if (col_type_use == col_type_subcol) then ! Set sub-column field if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_subcol))) & call endrun('pbuf_set_field: sub-column field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_subcol),const) else if (col_type_use == col_type_grid) then ! Set grid field if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_grid))) & call endrun('pbuf_set_field: grid field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_grid),const) else call endrun('pbuf_set_field: Trying to set '//trim(pbuf(index)%hdr%name)//& ' but col_type is neither col_type_grid nor col_type_subcol') end if end if end subroutine set_pbuf1d_field_const_by_index_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 subroutine set_pbuf2d_field_by_index_{DIMS}d_{TYPE}(pbuf2d,index,field, start, kount, col_type) ! Set the field(s) specified by the col_type type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: index integer,intent(in),optional :: start(:), kount(:) integer,intent(in),optional :: col_type logical :: subset integer :: c {VTYPE}, pointer :: fld{DIMSTR} #if ({DIMS}==1) {VTYPE},pointer :: field(:,:) #elif ({DIMS}==2) {VTYPE},pointer :: field(:,:,:) #elif ({DIMS}==3) {VTYPE},pointer :: field(:,:,:,:) #elif ({DIMS}==4) {VTYPE},pointer :: field(:,:,:,:,:) #elif ({DIMS}==5) {VTYPE},pointer :: field(:,:,:,:,:,:) #endif ! Check whether subset of data requested (default is false) subset = .false. if (present(start) .and. present(kount)) subset = .true. ! Check for ill-formed request if ( (present(start) .and. .not. present(kount)) .or.& (.not. present(start) .and. present(kount)) ) then call endrun('pbuf_set_field: Both start and kount must be present for '//trim(pbuf2d(index,begchunk)%hdr%name)) end if do c=begchunk,endchunk fld => get_field_chunk_{DIMS}d_{TYPE}(field,c) if(subset .and. present(col_type)) then call pbuf_set_field(pbuf_get_chunk(pbuf2d,c),index,fld,start,kount, col_type) else if(subset) then call pbuf_set_field(pbuf_get_chunk(pbuf2d,c),index,fld,start,kount) else if(present(col_type)) then call pbuf_set_field(pbuf_get_chunk(pbuf2d,c),index,fld,col_type=col_type) else call pbuf_set_field(pbuf_get_chunk(pbuf2d,c),index,fld) end if end do end subroutine set_pbuf2d_field_by_index_{DIMS}d_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 function get_field_chunk_{DIMS}d_{TYPE}(field, c) result(fld) ! module private helper function {VTYPE}, pointer :: fld{DIMSTR} integer, intent(in) :: c #if ({DIMS}==1) {VTYPE},pointer :: field(:,:) fld => field(:,c) #elif ({DIMS}==2) {VTYPE},pointer :: field(:,:,:) fld => field(:,:,c) #elif ({DIMS}==3) {VTYPE},pointer :: field(:,:,:,:) fld => field(:,:,:,c) #elif ({DIMS}==4) {VTYPE},pointer :: field(:,:,:,:,:) fld => field(:,:,:,:,c) #elif ({DIMS}==5) {VTYPE},pointer :: field(:,:,:,:,:,:) fld => field(:,:,:,:,:,c) #endif end function get_field_chunk_{DIMS}d_{TYPE} ! TYPE int,double,real ! DIMS 1,2,3,4,5 subroutine set_pbuf1d_field_by_index_{DIMS}d_{TYPE}(pbuf,index,field, start, kount, col_type) type(physics_buffer_desc), pointer :: pbuf(:) integer, intent(in) :: index {VTYPE}, intent(in) :: field{DIMSTR} integer,intent(in),optional :: start(:), kount(:) integer,intent(in),optional :: col_type integer :: col_type_use logical :: subset ! Default col_type is grid only if (present(col_type)) then col_type_use = col_type else col_type_use = col_type_grid end if ! Check whether subset of data requested (default is false) subset = .false. if (present(start) .and. present(kount)) subset = .true. ! Check for ill-formed request if ( (present(start) .and. .not. present(kount)) .or.& (.not. present(start) .and. present(kount)) ) then call endrun('pbuf_set_field: Both start and kount must be present for '//trim(pbuf(index)%hdr%name)) end if if(subset) then ! Set sub-column field if (col_type_use == col_type_subcol) then if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_subcol))) & call endrun('pbuf_set_field: sub-column field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_subcol),field,start,kount) ! Set grid field else if (col_type_use == col_type_grid) then if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_grid))) & call endrun('pbuf_set_field: grid field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_grid),field,start,kount) else call endrun('pbuf_set_field: Trying to set '//trim(pbuf(index)%hdr%name)//& ' but col_type is neigher col_type_grid nor col_type_subcol ') end if else ! Set sub-column field if (col_type_use == col_type_subcol) then if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_subcol))) & call endrun('pbuf_set_field: sub-column field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_subcol),field) ! Set grid field else if (col_type_use == col_type_grid) then if (.not. buffer_field_is_alloc(pbuf(index)%bfg(col_type_grid))) & call endrun('pbuf_set_field: grid field not allocated for '//trim(pbuf(index)%hdr%name)) call buffer_set_field(pbuf(index)%bfg(col_type_grid),field) else call endrun('pbuf_set_field: Trying to set '//trim(pbuf(index)%hdr%name)//& ' but col_type is neither col_type_grid nor col_type_subcol ') end if endif end subroutine set_pbuf1d_field_by_index_{DIMS}d_{TYPE} function pbuftype2piotype(pbuftype) result(piotype) use pio, only : pio_double, pio_int, pio_real integer, intent(in) :: pbuftype integer :: piotype select case(pbuftype) case (TYPEDOUBLE) piotype = pio_double case(TYPEINT) piotype = pio_int case(TYPEREAL) piotype = pio_real ! case(TYPELONG) ! piotype = pio_int case default write(iulog, *) 'Dtype = ', pbuftype call endrun('No restart support for dtype') end select end function pbuftype2piotype ! ! Initialize a restart file to write - all additional dims in a field are ! bundled into a single dimension for output and a dim pbuf_xxxxx is declared ! in the file if it does not already exist. ! subroutine pbuf_init_restart(File, pbuf2d) use pio, only: file_desc_t, pio_int use cam_pio_utils, only: cam_pio_def_dim, cam_pio_def_var use phys_grid, only: phys_decomp use cam_grid_support, only: cam_grid_get_file_dimids, cam_grid_dimensions ! Dummy Variables type(file_desc_t), intent(inout) :: file type(physics_buffer_desc), pointer :: pbuf2d(:,:) ! Local Variables type(physics_buffer_desc), pointer :: pbuf integer :: i, grid_select integer :: adimid(3) ! PIO IDs integer :: hdimcnt ! # grid dims integer :: dimcnt ! # array dims integer :: mdimsize, piodtype character(len=10) :: dimname character(len=24) :: varname ! Use adimid as a temp to find number of horizontal dims call cam_grid_dimensions(phys_decomp, adimid(1:2), hdimcnt) call cam_grid_get_file_dimids(phys_decomp, File, adimid) do i = 1, currentpbufflds pbuf => pbuf2d(i,begchunk) ! Only save global pbufs for restart if(pbuf%hdr%persistence /= persistence_global) cycle piodtype = pbuftype2piotype(pbuf%hdr%dtype) do grid_select = 1, ngrid_types ! For subcol fields, mdimsize includes psubcols in size mdimsize = product(pbuf%hdr%dimsizes(:,grid_select))/pcols if(mdimsize > 1) then write(dimname,'(a,i5.5)') 'pbuf_',mdimsize call cam_pio_def_dim(File, dimname, mdimsize, adimid(hdimcnt+1), & existOK=.true.) dimcnt = hdimcnt + 1 else dimcnt = hdimcnt end if if (mdimsize > 0) then varname = trim(pbuf%hdr%name)//trim(field_grid_suff(grid_select)) call cam_pio_def_var(File, varname, piodtype, adimid(1:dimcnt), & pbuf%hdr%vardesc(grid_select), existOK=.false.) end if end do end do call cam_pio_def_var(File, 'pbuf_time_idx', pio_int, timeidx_desc, & existOK=.false.) end subroutine pbuf_init_restart subroutine pbuf_write_restart(File, pbuf2d) use pio, only: file_desc_t, pio_put_var use cam_grid_support, only: cam_grid_dimensions use phys_grid, only: phys_decomp ! Dummy Variables type(file_desc_t), intent(inout) :: file type(physics_buffer_desc), pointer :: pbuf2d(:,:) ! Local Variables type(physics_buffer_desc), pointer :: pbufhdr(:) integer :: index, dtype, ierr integer :: gdimlens(2) integer :: grank pbufhdr => pbuf_get_chunk(pbuf2d, begchunk) gdimlens = 1 call cam_grid_dimensions(phys_decomp, gdimlens(1:2), grank) do index = 1, currentpbufflds if(pbufhdr(index)%hdr%persistence == persistence_global) then dtype = pbufhdr(index)%hdr%dtype select case(dtype) case (TYPEDOUBLE) call write_restart_field_double(File, pbuf2d, index, gdimlens, grank) case (TYPEREAL) call write_restart_field_real(File, pbuf2d, index, gdimlens, grank) case (TYPEINT) call write_restart_field_int(File, pbuf2d, index, gdimlens, grank) end select end if end do ierr = pio_put_var(File, timeidx_desc, (/old_time_idx/)) end subroutine pbuf_write_restart subroutine pbuf_restart_dimsizes(pbuf, gdimlens, grank, gridnum, adimlens, & fdimlens, frank) ! Dummy arguments type(physics_buffer_desc), pointer :: pbuf integer, intent(in) :: gdimlens(2) ! Grid horiz. dim sizes integer, intent(in) :: grank ! Array rank on file integer, intent(in) :: gridnum ! pbuf grid selector integer, intent(out) :: adimlens(3) ! Array dims integer, intent(out) :: fdimlens(3) ! Array dims on file integer, intent(out) :: frank ! array rank on file ! Local variable integer :: mdimsize mdimsize = PRODUCT(pbuf%hdr%dimsizes(:,gridnum)) / pcols fdimlens(1:2) = gdimlens if (grank == 1) then ! Unstructured grid fdimlens(2) = mdimsize fdimlens(3) = 1 frank = 2 else fdimlens(3) = mdimsize frank = 3 end if ! Restart does not write a dimension if mdimsize == 1 if (mdimsize == 1) then frank = frank - 1 end if adimlens(1) = pcols adimlens(2) = mdimsize adimlens(3) = endchunk - begchunk + 1 end subroutine pbuf_restart_dimsizes ! TYPE int,double,real subroutine write_restart_field_{TYPE}(File, pbuf2d, index, gdimlens, grank) use cam_grid_support, only: cam_grid_write_dist_array use phys_grid, only: phys_decomp use pio, only: file_desc_t ! Dummy Variables type(file_desc_t), intent(inout) :: File type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: index integer, intent(in) :: gdimlens(2) integer, intent(in) :: grank ! Local Variables type(physics_buffer_desc), pointer :: pbuf {VTYPE}, allocatable :: field(:,:,:) integer :: grid_select ! 1=grid, 2=subcol integer :: c ! Chunk index integer :: fdimlens(3) ! Array dimensions on NetCDF file integer :: adimlens(3) ! Array dimensions integer :: frank pbuf => pbuf2d(index, begchunk) do grid_select = 1, ngrid_types call pbuf_restart_dimsizes(pbuf, gdimlens, grank, grid_select, & adimlens, fdimlens, frank) if ((PRODUCT(adimlens) > 0) .and. (.not. pbuf%hdr%is_copy(grid_select))) then allocate(field(adimlens(1), adimlens(2), begchunk:endchunk)) do c = begchunk, endchunk call pbuf_get_field_restart(pbuf2d, c, index, field(:,:,c), & adimlens(2), grid_select) end do if (size(field,2) == 1) then ! Special case for 2-D pbuf fields adimlens(2) = adimlens(3) call cam_grid_write_dist_array(File, phys_decomp, adimlens(1:2), & fdimlens(1:frank), field(:,1,:), pbuf%hdr%vardesc(grid_select)) else call cam_grid_write_dist_array(File, phys_decomp, adimlens, & fdimlens(1:frank), field, pbuf%hdr%vardesc(grid_select)) end if deallocate(field) end if end do end subroutine write_restart_field_{TYPE} subroutine pbuf_read_restart(File, pbuf2d) use cam_grid_support, only: cam_grid_dimensions use phys_grid, only: phys_decomp use pio, only: file_desc_t, pio_inq_varid, pio_get_var ! Dummy Variables type(File_desc_t), intent(inout) :: File type(physics_buffer_desc), pointer :: pbuf2d(:,:) ! Local Variables type(physics_buffer_desc), pointer :: pbufhdr(:) integer :: index, dtype, ierr integer :: gdimlens(2) ! Horizontal grid dimensions integer :: grank call pbuf_initialize(pbuf2d) pbufhdr => pbuf_get_chunk(pbuf2d, begchunk) gdimlens = 1 call cam_grid_dimensions(phys_decomp, gdimlens(1:2), grank) ierr = pio_inq_varid(File, 'pbuf_time_idx', timeidx_desc) ierr = pio_get_var(File, timeidx_desc, old_time_idx) do index = 1, currentpbufflds if(pbufhdr(index)%hdr%persistence == persistence_global) then dtype = pbufhdr(index)%hdr%dtype select case(dtype) case (TYPEDOUBLE) call read_restart_field_double(File, pbuf2d, index, gdimlens, grank) case (TYPEREAL) call read_restart_field_real(File, pbuf2d, index, gdimlens, grank) case (TYPEINT) call read_restart_field_int(File, pbuf2d, index, gdimlens, grank) end select end if end do end subroutine pbuf_read_restart ! TYPE int,double,real subroutine read_restart_field_{TYPE} (File, pbuf2d, index, gdimlens, grank) use pio, only: file_desc_t use pio, only: pio_inq_varid use cam_pio_utils, only: cam_pio_handle_error use cam_grid_support, only: cam_grid_read_dist_array use subcol_utils, only: subcol_pack use phys_grid, only: phys_decomp ! Dummy Variables type(file_desc_t), intent(inout) :: File type(physics_buffer_desc), pointer :: pbuf2d(:,:) integer, intent(in) :: index integer, intent(in) :: gdimlens(2) integer, intent(in) :: grank ! Local Variables type(physics_buffer_desc), pointer :: pbuf {VTYPE}, pointer :: fld6(:,:,:,:,:,:) {VTYPE}, allocatable :: fld5_pack(:,:,:,:,:) {VTYPE}, allocatable :: field(:,:,:) integer :: grid_select ! (1=grid, 2=subcol) integer :: ierr, c integer :: start(6) integer :: dimsizes(6) integer :: fdimlens(3) ! Array dimensions on NetCDF file integer :: frank ! Array rank on file integer :: adimlens(3) ! Array dimensions character(len=24) :: varname character(len=*), parameter :: subname = 'read_restart_field_{TYPE}' pbuf => pbuf2d(index, begchunk) start = 1 do grid_select = 1, ngrid_types dimsizes(:) = pbuf%hdr%dimsizes(:, grid_select) if(all(dimsizes(:) == 0)) then ! None of this grid type for this variable cycle end if call pbuf_restart_dimsizes(pbuf, gdimlens, grank, grid_select, & adimlens, fdimlens, frank) ! Fix up dimensions for subcolumn field if (grid_select == col_type_subcol) then ! Field stored as (pcols, psubcols*restDims, chunks) ! Read then pack to (pcols*psubcols,restDims, chunks) allocate(fld5_pack(dimsizes(1),dimsizes(2),dimsizes(3),dimsizes(4),dimsizes(5))) do c = 5, 2, -1 dimsizes(c + 1) = dimsizes(c) end do dimsizes(1) = pcols dimsizes(2) = psubcols end if varname = trim(pbuf%hdr%name)//trim(field_grid_suff(grid_select)) ierr = pio_inq_varid(File, varname, pbuf%hdr%vardesc(grid_select)) call cam_pio_handle_error(ierr, trim(subname)//': '//trim(varname)//' not found') allocate(field(adimlens(1), adimlens(2), begchunk:endchunk)) if (size(field,2) == 1) then ! Special case for 2-D pbuf fields adimlens(2) = adimlens(3) call cam_grid_read_dist_array(File, phys_decomp, adimlens(1:2), & fdimlens(1:frank), field(:,1,:), pbuf%hdr%vardesc(grid_select)) else call cam_grid_read_dist_array(File, phys_decomp, adimlens, & fdimlens(1:frank), field, pbuf%hdr%vardesc(grid_select)) end if do c = begchunk, endchunk if (grid_select == col_type_grid) then call buffer_get_field_ptr(pbuf2d(index,c)%bfg(grid_select), fld6, & start, dimsizes) fld6 = reshape(field(:,:,c), dimsizes) else if (grid_select == col_type_subcol) then call subcol_pack(c, reshape(field(:,:,c), dimsizes), fld5_pack) call buffer_set_field(pbuf2d(index,c)%bfg(grid_select), fld5_pack) else call endrun('read_restart_field_{TYPE}: invalid grid selector - must be either 1 or 2') end if end do nullify(fld6) deallocate(field) if (allocated(fld5_pack)) then deallocate(fld5_pack) end if end do end subroutine read_restart_field_{TYPE} subroutine pbuf_dump_pbuf(pbuf2d, name, num) use cam_pio_utils, only: cam_pio_dump_field use spmd_utils, only: masterproc ! Dummy arguments type(physics_buffer_desc), pointer :: pbuf2d(:,:) character(len=*), optional, intent(in) :: name integer, optional, intent(in) :: num ! Local variables integer, parameter :: max_name = 64 character(len=max_name) :: field_name integer :: index, grid_select, c integer :: namelen integer :: dimstart(6), dimend(6) integer :: ierr integer :: dtype type(physics_buffer_desc), pointer :: pbufhdr(:) real(r8), allocatable :: field(:,:,:,:,:,:) real(r8), pointer :: fld5(:,:,:,:,:) pbufhdr => pbuf_get_chunk(pbuf2d, begchunk) dimstart = 1 do index = 1, currentpbufflds dtype = pbufhdr(index)%hdr%dtype do grid_select = 1, ngrid_types namelen = len_trim(pbufhdr(index)%hdr%name) namelen = namelen + len_trim(field_grid_suff(grid_select)) if (present(name)) then namelen = namelen + len_trim(name) end if if (present(num)) then namelen = namelen + int(log10(real(num))) + 2 end if if (namelen > 64) then call endrun("PBUF_DUMP_PBUF: Name string too long") end if if (present(name)) then if (present(num)) then write(field_name, '(4a,i0)') trim(pbufhdr(index)%hdr%name), & trim(field_grid_suff(grid_select)), trim(name), '_', num else write(field_name, '(3a)') trim(pbufhdr(index)%hdr%name), & trim(field_grid_suff(grid_select)), trim(name) end if else if (present(num)) then write(field_name, '(3a,i0)') trim(pbufhdr(index)%hdr%name), & trim(field_grid_suff(grid_select)), '_', num else write(field_name, '(2a)') trim(pbufhdr(index)%hdr%name), & trim(field_grid_suff(grid_select)) end if dimend = pbufhdr(index)%hdr%dimsizes(:, grid_select) if (PRODUCT(dimend) > 0) then if (dimend(6) /= 1) then if (masterproc) then write(iulog, *) 'PBUF_DUMP_PBUF: ', trim(field_name), dimend end if call endrun("PBUF_DUMP_PBUF: No support for 6D pbuf field") end if dimend(6) = endchunk - begchunk + 1 select case(dtype) case (TYPEDOUBLE) allocate(field(dimend(1), dimend(2), dimend(3), dimend(4), & dimend(5), dimend(6))) do c = begchunk, endchunk call pbuf_get_field(pbuf2d, c, index, fld5, & col_type=grid_select, errcode=ierr) field(:,:,:,:,:,c-begchunk+1) = fld5(:,:,:,:,:) end do call cam_pio_dump_field(field_name, dimstart, dimend, & field, compute_maxdim_in=.false.) deallocate(field) case (TYPEREAL) if (masterproc) then write(iulog, *) 'PBUF_DUMP_PBUF: No support for real fields' end if case (TYPEINT) if (masterproc) then write(iulog, *) 'PBUF_DUMP_PBUF: No support for integer fields' end if end select end if end do end do end subroutine pbuf_dump_pbuf end module physics_buffer