diff --git a/.gitmodules b/.gitmodules index 1325b4228..387d23fcd 100644 --- a/.gitmodules +++ b/.gitmodules @@ -6,8 +6,8 @@ fxDONOTUSEurl = https://github.com/NCAR/ccpp-framework [submodule "history"] path = src/history/buffers - url = https://github.com/ESMCI/history_output - fxtag = history01_02 + url = https://github.com/peverwhee/history_output + fxtag = a8e4bf4ffd69e1332aebab7f5d34ddf2ffb228fd fxrequired = AlwaysRequired fxDONOTUSEurl = https://github.com/ESMCI/history_output [submodule "mpas"] diff --git a/cime_config/buildlib b/cime_config/buildlib index 4d9222c8b..64834b1ee 100755 --- a/cime_config/buildlib +++ b/cime_config/buildlib @@ -84,6 +84,7 @@ def _build_cam(): dycore = config.get_value('dyn') reg_dir = config.get_value('reg_dir') init_dir = config.get_value('init_dir') + restart_dir = config.get_value('restart_dir') phys_dirs_str = config.get_value('phys_dirs') #Convert the phys_dirs_str into a proper list: @@ -96,7 +97,7 @@ def _build_cam(): filepath_src = os.path.join(caseroot, "Buildconf", "camconf", "Filepath") filepath_dst = os.path.join(bldroot, "Filepath") - paths = [source_mods_dir, reg_dir, init_dir, + paths = [source_mods_dir, reg_dir, init_dir, restart_dir, os.path.join(atm_root, "src", "data"), os.path.join(atm_root, "src", "control"), os.path.join(atm_root, "src", "cpl", diff --git a/cime_config/cam_autogen.py b/cime_config/cam_autogen.py index 997e793ef..6711652fc 100644 --- a/cime_config/cam_autogen.py +++ b/cime_config/cam_autogen.py @@ -33,6 +33,7 @@ # Import needed registry and other src/data scripts: from generate_registry_data import gen_registry from write_init_files import write_init_files +from write_restart_physics import write_restart_physics ############################################################################### @@ -454,7 +455,7 @@ def generate_registry(data_search, build_cache, atm_root, bldroot, gen_fort_indent, source_mods_dir, atm_root, logger=_LOGGER, schema_paths=data_search, error_on_no_validate=True) - retcode, reg_file_list, ic_names, registry_constituents, vars_init_value = retvals + retcode, reg_file_list, ic_names, registry_constituents, restart_vars, vars_init_value = retvals # Raise error if gen_registry failed: if retcode != 0: emsg = "ERROR:Unable to generate CAM data structures from {}, err = {}" @@ -468,16 +469,17 @@ def generate_registry(data_search, build_cache, atm_root, bldroot, # Save build details in the build cache reg_file_paths = [x.file_path for x in reg_file_list if x.file_path] build_cache.update_registry(gen_reg_file, registry_files, dycore, - reg_file_paths, ic_names, registry_constituents, vars_init_value) + reg_file_paths, ic_names, registry_constituents, restart_vars, vars_init_value) else: # If we did not run the registry generator, retrieve info from cache reg_file_paths = build_cache.reg_file_list() ic_names = build_cache.ic_names() registry_constituents = build_cache.constituents() + restart_vars = build_cache.restart_vars() vars_init_value = build_cache.vars_init_value() # End if - return genreg_dir, do_gen_registry, reg_file_paths, ic_names, registry_constituents, vars_init_value + return genreg_dir, do_gen_registry, reg_file_paths, ic_names, registry_constituents, restart_vars, vars_init_value ############################################################################### def generate_physics_suites(build_cache, preproc_defs, host_name, @@ -808,6 +810,55 @@ def generate_init_routines(build_cache, bldroot, force_ccpp, force_init, return init_dir +############################################################################### +def generate_restart_routines(build_cache, bldroot, force_ccpp, force_restart, + source_mods_dir, gen_fort_indent, cap_database, + ic_names, registry_constituents, restart_vars): +############################################################################### + """ + Generate the physics restart routines (restart_physics.F90) using + both the registry and the CCPP physics suites if required + (new case or changes to registry or CCPP source(s), meta-data, + and/or script). + """ + + #Add new directory to build path: + restart_dir = os.path.join(bldroot, "phys_restart") + # Use this for cache check + gen_restart_file = os.path.join(_REG_GEN_DIR, "write_restart_physics.py") + + # Figure out if we need to generate new restart routines: + if os.path.exists(restart_dir): + # Check if registry and / or CCPP suites were modified: + if force_ccpp or force_restart: + do_gen_restart = True + else: + #If not, then check cache to see if actual + #"restart_physics.py" was modified: + do_gen_restart = build_cache.restart_write_mismatch(gen_restart_file) + else: + #If no directory exists, then one will need + # to create new routines: + os.mkdir(restart_dir) + do_gen_restart = True + # End if + if do_gen_restart: + source_paths = [source_mods_dir, _REG_GEN_DIR] + retmsg = write_restart_physics(cap_database, ic_names, registry_constituents, + restart_vars, restart_dir, _find_file, + source_paths, gen_fort_indent, _LOGGER) + + #Check that script ran properly + if retmsg: + emsg = f"ERROR: Unable to generate CAM restart source code, error message is:\n{retmsg}" + raise CamAutoGenError(emsg) + #----- + + # save build details in the build cache + build_cache.update_restart_gen(gen_restart_file) + + return restart_dir + ############# # End of file ############# diff --git a/cime_config/cam_build_cache.py b/cime_config/cam_build_cache.py index 3e23b74a8..94f2b0d51 100644 --- a/cime_config/cam_build_cache.py +++ b/cime_config/cam_build_cache.py @@ -213,6 +213,7 @@ def __init__(self, build_cache): # Set empty values sure to trigger processing self.__gen_reg_file = None self.__gen_init_file = None + self.__gen_restart_file = None self.__registry_files = {} self.__dycore = None self.__sdfs = {} @@ -227,6 +228,7 @@ def __init__(self, build_cache): self.__reg_gen_files = [] self.__ic_names = {} self.__constituents = [] + self.__restart_dict = {} self.__vars_init_value = [] if os.path.exists(build_cache): # Initialize build cache state @@ -240,6 +242,9 @@ def __init__(self, build_cache): elif item.tag == 'generate_init_file': new_entry = new_entry_from_xml(item) self.__gen_init_file = new_entry + elif item.tag == 'generate_restart_file': + new_entry = new_entry_from_xml(item) + self.__gen_restart_file = new_entry elif item.tag == 'registry_file': new_entry = new_entry_from_xml(item) self.__registry_files[new_entry.key] = new_entry @@ -258,6 +263,10 @@ def __init__(self, build_cache): stdname = item.get('standard_name') itext = clean_xml_text(item) self.__constituents.append(itext) + elif item.tag == 'restart_entry': + stdname = item.get('standard_name') + itext = clean_xml_text(item) + self.__restart_dict[stdname] = itext elif item.tag == 'vars_init_value_entry': itext = clean_xml_text(item) self.__vars_init_value.append(itext) @@ -322,7 +331,8 @@ def __init__(self, build_cache): # end if def update_registry(self, gen_reg_file, registry_source_files, - dycore, reg_file_list, ic_names, constituents, vars_init_value): + dycore, reg_file_list, ic_names, constituents, + restart_dict, vars_init_value): """Replace the registry cache data with input data """ self.__dycore = dycore @@ -338,6 +348,7 @@ def update_registry(self, gen_reg_file, registry_source_files, # and should already be of type dict: self.__ic_names = ic_names self.__constituents = constituents + self.__restart_dict = restart_dict self.__vars_init_value = vars_init_value def update_ccpp(self, suite_definition_files, scheme_files, host_files, @@ -383,6 +394,14 @@ def update_init_gen(self, gen_init_file): """ self.__gen_init_file = FileStatus(gen_init_file, 'generate_init_file') + def update_restart_gen(self, gen_restart_file): + """ + Replace the restart writer + (write_restart_physics.py) cache + data with input data + """ + self.__gen_restart_file = FileStatus(gen_restart_file, 'generate_restart_file') + def write(self): """Write out the current cache state""" new_cache = ET.Element("CAMBuildCache") @@ -394,6 +413,9 @@ def write(self): new_xml_entry(registry, 'generate_registry_file', self.__gen_reg_file.file_path, self.__gen_reg_file.file_hash) + new_xml_entry(registry, 'generate_restart_file', + self.__gen_restart_file.file_path, + self.__gen_restart_file.file_hash) for rfile in self.__registry_files.values(): new_xml_entry(registry, 'registry_file', rfile.file_path, rfile.file_hash) @@ -415,6 +437,11 @@ def write(self): const_entry = ET.SubElement(registry, 'constituent_entry') const_entry.text = stdname # end for + for stdname, restart_name in self.__restart_dict.items(): + restart_entry = ET.SubElement(registry, 'restart_entry') + restart_entry.set('standard_name', stdname) + restart_entry.text = restart_name + # end for for stdname in self.__vars_init_value: var_entry = ET.SubElement(registry, 'vars_init_value_entry') var_entry.text = stdname @@ -603,6 +630,22 @@ def init_write_mismatch(self, gen_init_file): #Return mismatch logical: return mismatch + def restart_write_mismatch(self, gen_restart_file): + """ + Determine if the restart_files writer (write_restart_files.py) + differs from the data stored in our cache. Return True + if the data differs. + """ + + # Initialize variable + mismatch = False + + # Check file hash to see if mis-match exists: + mismatch = self.__gen_restart_file.hash_mismatch(gen_restart_file) + + # Return mismatch logical: + return mismatch + def scheme_nl_metadata(self): """Return the stored list of scheme namelist metadata files""" return self.__scheme_nl_metadata @@ -625,6 +668,10 @@ def constituents(self): """Return a copy of the registry constituents list""" return list(self.__constituents) + def restart_vars(self): + """Return a copy of the registry's list of variables for the restart file""" + return dict(self.__restart_dict) + def vars_init_value(self): """Return a copy of the list of variables with initial_value""" return list(self.__vars_init_value) diff --git a/cime_config/cam_config.py b/cime_config/cam_config.py index 2224ded85..85752c1f5 100644 --- a/cime_config/cam_config.py +++ b/cime_config/cam_config.py @@ -32,6 +32,7 @@ # Import fortran auto-generation routines: from cam_autogen import generate_registry, generate_physics_suites from cam_autogen import generate_init_routines +from cam_autogen import generate_restart_routines ############################################################################### #HELPER FUNCTIONS @@ -867,7 +868,7 @@ def generate_cam_src(self, gen_fort_indent): retvals = generate_registry(data_search, build_cache, self.__atm_root, self.__bldroot, source_mods_dir, dyn, gen_fort_indent) - reg_dir, force_ccpp, reg_files, ic_names, registry_constituents, vars_init_value = retvals + reg_dir, force_ccpp, reg_files, ic_names, registry_constituents, restart_vars, vars_init_value = retvals #Add registry path to config object: reg_dir_desc = "Location of auto-generated registry code." @@ -907,6 +908,19 @@ def generate_cam_src(self, gen_fort_indent): init_dir_desc = "Location of auto-generated physics initialization code." self.create_config("init_dir", init_dir_desc, init_dir) + #--------------------------------------------------------- + # Create restart init/write/read routines + #--------------------------------------------------------- + restart_dir = generate_restart_routines(build_cache, self.__bldroot, + force_ccpp, force_init, + source_mods_dir, gen_fort_indent, + capgen_db, ic_names, + registry_constituents, restart_vars) + + #Add registry path to config object: + restart_dir_desc = "Location of auto-generated physics restart code." + self.create_config("restart_dir", restart_dir_desc, restart_dir) + #-------------------------------------------------------------- # write out the cache here as we have completed pre-processing #-------------------------------------------------------------- diff --git a/src/control/cam_comp.F90 b/src/control/cam_comp.F90 index 186299fc1..b8f9ca25d 100644 --- a/src/control/cam_comp.F90 +++ b/src/control/cam_comp.F90 @@ -467,51 +467,24 @@ end subroutine cam_run3 !----------------------------------------------------------------------- ! - subroutine cam_run4(rstwr, nlend, & - yr_spec, mon_spec, day_spec, sec_spec) + subroutine cam_run4(rstwr, nlend) + logical, intent(in) :: rstwr ! write restart file + logical, intent(in) :: nlend ! this is final timestep !----------------------------------------------------------------------- ! - ! Purpose: Final phase of atmosphere model run method. This consists - ! of all the restart output, history writes, and other - ! file output. + ! Purpose: ! !----------------------------------------------------------------------- -! use cam_restart, only: cam_write_restart ! use qneg_module, only: qneg_print_summary - logical, intent(in) :: rstwr ! write restart file - logical, intent(in) :: nlend ! this is final timestep - integer, intent(in), optional :: yr_spec ! Simulation year - integer, intent(in), optional :: mon_spec ! Simulation month - integer, intent(in), optional :: day_spec ! Simulation day - integer, intent(in), optional :: sec_spec ! Secs in current simulation day - - ! - ! Write restart files - ! - if (rstwr) then - call t_startf('cam_write_restart') - if (present(yr_spec) .and. present(mon_spec) .and. & - present(day_spec).and.present(sec_spec)) then -!!XXgoldyXX: v need to import this -! call cam_write_restart(cam_in, cam_out, dyn_out, yr_spec=yr_spec, & -! mon_spec=mon_spec, day_spec=day_spec, sec_spec= sec_spec) -!!XXgoldyXX: ^ need to import this - else -!!XXgoldyXX: v need to import this -! call cam_write_restart(cam_in, cam_out, dyn_out) -!!XXgoldyXX: ^ need to import this - end if - call t_stopf('cam_write_restart') - end if - end subroutine cam_run4 ! !----------------------------------------------------------------------- ! - subroutine cam_timestep_final(rstwr, nlend, do_ncdata_check, do_history_write) + subroutine cam_timestep_final(rstwr, nlend, do_ncdata_check, & + yr_spec, mon_spec, day_spec, sec_spec, do_history_write) !----------------------------------------------------------------------- ! ! Purpose: Timestep final runs at the end of each timestep @@ -521,15 +494,20 @@ subroutine cam_timestep_final(rstwr, nlend, do_ncdata_check, do_history_write) use phys_comp, only: phys_timestep_final use cam_history, only: history_write_files use cam_history, only: history_wrap_up + use cam_restart, only: cam_write_restart logical, intent(in) :: rstwr ! write restart file logical, intent(in) :: nlend ! this is final timestep !Flag for whether a snapshot (ncdata) check should be run or not ! - flag is true if this is not the first or last step logical, intent(in) :: do_ncdata_check + integer, intent(in), optional :: yr_spec ! Simulation year + integer, intent(in), optional :: mon_spec ! Simulation month + integer, intent(in), optional :: day_spec ! Simulation day + integer, intent(in), optional :: sec_spec ! Secs in current simulation day !Flag for whether to perform the history write logical, optional, intent(in) :: do_history_write - logical :: history_write_loc + logical :: history_write_loc if (present(do_history_write)) then history_write_loc = do_history_write @@ -540,7 +518,20 @@ subroutine cam_timestep_final(rstwr, nlend, do_ncdata_check, do_history_write) if (history_write_loc) then call history_write_files() end if - ! peverwhee - todo: handle restarts + ! + ! Write restart files + ! + if (rstwr) then + call t_startf('cam_write_restart') + if (present(yr_spec) .and. present(mon_spec) .and. & + present(day_spec).and.present(sec_spec)) then + call cam_write_restart(dyn_out, yr_spec=yr_spec, & + mon_spec=mon_spec, day_spec=day_spec, sec_spec= sec_spec) + else + call cam_write_restart(dyn_out) + end if + call t_stopf('cam_write_restart') + end if call history_wrap_up(rstwr, nlend) ! diff --git a/src/control/cam_restart.F90 b/src/control/cam_restart.F90 new file mode 100644 index 000000000..359b730d6 --- /dev/null +++ b/src/control/cam_restart.F90 @@ -0,0 +1,149 @@ +module cam_restart + + implicit none + private + + public :: cam_write_restart + +CONTAINS + subroutine cam_write_restart(dyn_out, yr_spec, mon_spec, day_spec, sec_spec) + use cam_filenames, only: interpret_filename_spec + use cam_pio_utils, only: cam_pio_createfile, cam_pio_set_fill + use restart_dynamics, only: write_restart_dynamics, init_restart_dynamics + use restart_physics, only: restart_physics_write, restart_physics_init + use cam_history, only: history_restart_init, history_restart_write + use cam_instance, only: inst_suffix + use pio, only: file_desc_t, io_desc_t, pio_double, pio_global + use pio, only: pio_put_att, pio_enddef, pio_closefile + use cam_grid_support, only: cam_grid_write_attr, cam_grid_id + use cam_grid_support, only: cam_grid_header_info_t, cam_grid_write_var + use cam_grid_support, only: cam_grid_dimensions, cam_grid_get_decomp + use physics_grid, only: phys_decomp, num_global_phys_cols + use dyn_comp, only: dyn_export_t + use cam_control_mod, only: caseid + use cam_abortutils, only: endrun + use shr_kind_mod, only: cl=>shr_kind_cl + use runtime_obj, only: cam_runtime_opts + + ! Arguments + type(dyn_export_t), intent(in) :: dyn_out + integer, optional, intent(in) :: yr_spec ! Simulation year + integer, optional, intent(in) :: mon_spec ! Simulation month + integer, optional, intent(in) :: day_spec ! Simulation day + integer, optional, intent(in) :: sec_spec ! Seconds into current simulation day + + ! Local variables + character(len=cl) :: rfilename_spec ! filename specifier for primary restart file + character(len=cl) :: fname ! Restart filename + type(file_desc_t) :: fh + integer :: ierr, i, errflg + character(len=512) :: errmsg + integer :: grid_id + integer :: dims(3), gdims(3) + integer :: nhdims, ndims + type(io_desc_t), pointer :: iodesc + type(cam_grid_header_info_t) :: info + + ! Set template for primary restart filename based on instance suffix + ! (%c = caseid, $y = year, $m = month, $d = day, $s = seconds in day, %t = number) + rfilename_spec = '%c.cam' // trim(inst_suffix) //'.r.%y-%m-%d-%s.nc' + fname = interpret_filename_spec(rfilename_spec, yr_spec=yr_spec, & + mon_spec=mon_spec, day_spec=day_spec, sec_spec= sec_spec) + + call cam_pio_createfile(fh, trim(fname), 0) + ierr = cam_pio_set_fill(fh) + + call init_restart_dynamics(fh, dyn_out) + + ! TODO: initialize ionosphere restart + + ! Initialize physics grid variable + grid_id = cam_grid_id('physgrid') + call cam_grid_write_attr(fh, grid_id, info) + + call restart_physics_init(fh, errmsg, errflg) + if (errflg /= 0) then + call endrun(errmsg) + end if + call history_restart_init(fh) + + ierr = pio_put_att(fh, pio_global, 'caseid', caseid) + + ! The MPAS dycore will end the definition phase for us + if (cam_runtime_opts%get_dycore() /= 'mpas') then + ierr = pio_enddef(fh) + end if + + ! Restart dynamics must come first in write phase + ! because MPAS restart write includes the definitions + call write_restart_dynamics(fh, dyn_out) + + ! TODO: write ionosphere restart + + ! Write physics grid info + call cam_grid_write_var(fh, phys_decomp) + + ! Write physics restart variables + call restart_physics_write(fh, grid_id, errmsg, errflg) + if (errflg /= 0) then + call endrun(errmsg) + end if + + call history_restart_write(fh) + + ! Close the primary restart file + call pio_closefile(fh) + + ! Update the restart pointer file + call write_rest_pfile(fname, yr_spec=yr_spec, mon_spec=mon_spec, & + day_spec=day_spec, sec_spec= sec_spec ) + + + end subroutine cam_write_restart + +!======================================================================================== + + subroutine write_rest_pfile(restart_file, yr_spec, mon_spec, day_spec, sec_spec) + + ! Write the restart pointer file + use cam_instance, only: inst_suffix + use cam_filenames, only: interpret_filename_spec + use shr_kind_mod, only: cl=>shr_kind_cl + use spmd_utils, only: masterproc + use cam_logfile, only: iulog + use cam_abortutils, only: endrun + use ioFileMod, only: cam_open_file + + character(len=*), intent(in) :: restart_file + integer, optional, intent(in) :: yr_spec ! Simulation year + integer, optional, intent(in) :: mon_spec ! Simulation month + integer, optional, intent(in) :: day_spec ! Simulation day + integer, optional, intent(in) :: sec_spec ! Seconds into current simulation day + + integer :: ierr, unitn + character(len=CL) :: rest_pfile + character(len=*), parameter :: sub='write_rest_pfile' + !--------------------------------------------------------------------------- + + if (masterproc) then + rest_pfile = interpret_filename_spec('rpointer.cam'//trim(inst_suffix)//'.'//'%y-%m-%d-%s',& + yr_spec=yr_spec, mon_spec=mon_spec, day_spec=day_spec, sec_spec= sec_spec ) + call cam_open_file(rest_pfile, unitn, 'f', status="unknown") + rewind unitn + write(unitn, '(a)', iostat=ierr) trim(restart_file) + if (ierr /= 0) then + call endrun(sub//': ERROR: writing rpointer file') + end if + close(unitn) + + write(iulog,*)'(WRITE_REST_PFILE): successfully wrote local restart pointer file ',& + trim(rest_pfile) + write(iulog,'("---------------------------------------")') + end if + + end subroutine write_rest_pfile + +!======================================================================================== + + +end module cam_restart diff --git a/src/cpl/nuopc/atm_comp_nuopc.F90 b/src/cpl/nuopc/atm_comp_nuopc.F90 index d67b08c4a..6535aaddc 100644 --- a/src/cpl/nuopc/atm_comp_nuopc.F90 +++ b/src/cpl/nuopc/atm_comp_nuopc.F90 @@ -1183,10 +1183,10 @@ subroutine ModelAdvance(gcomp, rc) call t_stopf ('CAM_run3') call t_startf ('CAM_run4') - call cam_run4( rstwr, nlend, & - yr_spec=yr_sync, mon_spec=mon_sync, day_spec=day_sync, sec_spec=tod_sync) + call cam_run4(rstwr, nlend) call t_stopf ('CAM_run4') - call cam_timestep_final(rstwr, nlend, do_ncdata_check=do_ncdata_check) + call cam_timestep_final(rstwr, nlend, do_ncdata_check=do_ncdata_check, & + yr_spec=yr_sync, mon_spec=mon_sync, day_spec=day_sync, sec_spec=tod_sync) ! Advance cam time step diff --git a/src/data/generate_registry_data.py b/src/data/generate_registry_data.py index a738fc78c..ce888ed81 100755 --- a/src/data/generate_registry_data.py +++ b/src/data/generate_registry_data.py @@ -166,6 +166,8 @@ def __init__(self, elem_node, local_name, dimensions, allocation_dimensions, self.__allocatable = elem_node.get('allocatable', default=alloc_default) self.__constituent = elem_node.get("constituent", default=False) self.__advected = elem_node.get("advected", default=False) + self.__restart = elem_node.get("restart", default=False) + self.__diag_name = elem_node.get("diagnostic_name", default=local_name) self.__tstep_init = elem_node.get("phys_timestep_init_zero", default=tstep_init_default) if self.__allocatable == "none": @@ -437,6 +439,11 @@ def is_advected(self): """Return True if this variable is advected""" return self.__advected + @property + def is_restart(self): + """Return True if this variable should be included on the physics restart file""" + return self.__restart + @property def tstep_init(self): """Return True if variable will be set to zero every physics timestep.""" @@ -544,7 +551,8 @@ class Variable(VarBase): "constituent", "dycore", "extends", "kind", "local_name", "name", "phys_timestep_init_zero", "standard_name", - "type", "units", "version"] + "type", "units", "version", "restart", + "diagnostic_name"] def __init__(self, var_node, known_types, vdict, dycore, logger): # pylint: disable=too-many-locals @@ -1826,6 +1834,35 @@ def _create_constituent_list(registry): # end for return constituent_list +############################################################################### +def _create_restart_dict(registry): +############################################################################### + """ + Create a list of all registry variables that need to be included in + the physics restart file. + To be used by write_physics_restart.py + """ + restart_list = {} + for section in registry: + if section.tag == 'file': + for obj in section: + if obj.tag == 'variable': + if obj.get('restart'): + stdname = obj.get('standard_name') + diagnostic_name = obj.get('diagnostic_name') + if diagnostic_name: + restart_list[stdname] = diagnostic_name + else: + local_name = obj.get('local_name') + restart_list[stdname] = local_name + # end if + # end if (ignore non-restart variables) + # end if (ignore other node types) + # end for + # end if (ignore other node types) + # end for + return restart_list + ############################################################################### def _create_variables_with_initial_value_list(registry): ############################################################################### @@ -1922,10 +1959,11 @@ def gen_registry(registry_file, dycore, outdir, indent, # See comment in _create_ic_name_dict ic_names = _create_ic_name_dict(registry) registry_constituents = _create_constituent_list(registry) + restart_vars = _create_restart_dict(registry) vars_init_value = _create_variables_with_initial_value_list(registry) retcode = 0 # Throw exception on error # end if - return retcode, files, ic_names, registry_constituents, vars_init_value + return retcode, files, ic_names, registry_constituents, restart_vars, vars_init_value def main(): """Function to execute when module called as a script""" diff --git a/src/data/registry.xml b/src/data/registry.xml index fd9a8aa0c..dc458a99e 100644 --- a/src/data/registry.xml +++ b/src/data/registry.xml @@ -139,14 +139,16 @@ + allocatable="pointer" restart="true" + diagnostic_name="FRONTGF"> horizontal_dimension vertical_layer_dimension frontgf pbuf_FRONTGF + allocatable="pointer" restart="true" + diagnostic_name="FRONTGA"> horizontal_dimension vertical_layer_dimension frontga pbuf_FRONTGA @@ -161,7 +163,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="QPERT"> horizontal_dimension 0.0_kind_phys qpert pbuf_qpert @@ -177,7 +180,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="CMFMC_DP"> horizontal_dimension vertical_interface_dimension 0.0_kind_phys cmfmc_dp pbuf_CMFMC_DP @@ -185,7 +189,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="CMFMC_SH"> horizontal_dimension vertical_interface_dimension 0.0_kind_phys cmfmc_sh pbuf_CMFMC_SH @@ -193,7 +198,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="DLFZM"> horizontal_dimension vertical_layer_dimension 0.0_kind_phys dlfzm pbuf_DLFZM @@ -217,7 +223,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="RPRDDP"> horizontal_dimension vertical_layer_dimension 0.0_kind_phys pbuf_RPRDDP @@ -605,7 +612,7 @@ + units="kg kg-1" type="real" constituent="true" restart="true" diagnostic_name="cldliq"> Cloud water mass mixing ratio with respect to moist air plus all airborne condensates horizontal_dimension vertical_layer_dimension CLDLIQ cnst_CLDLIQ @@ -681,7 +688,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="LANDM"> horizontal_dimension 0.0_kind_phys @@ -690,7 +698,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="AST"> horizontal_dimension vertical_layer_dimension 0.0_kind_phys pbuf_AST @@ -732,7 +741,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="CFLX"> horizontal_dimension number_of_ccpp_constituents 0.0_kind_phys cflx cam_in_cflx @@ -922,7 +932,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="SHF"> horizontal_dimension 0.0_kind_phys shf cam_in_shf @@ -1027,7 +1038,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="wsx"> horizontal_dimension 0.0_kind_phys wsx cam_in_wsx @@ -1035,7 +1047,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="wsy"> horizontal_dimension 0.0_kind_phys wsy cam_in_wsy @@ -1045,7 +1058,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="BCPHIDRY"> horizontal_dimension 0.0_kind_phys bcphidry cam_out_bcphidry @@ -1061,7 +1075,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="BCPHODRY"> horizontal_dimension 0.0_kind_phys bcphodry cam_out_bcphodry @@ -1069,7 +1084,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="CO2DIAG"> horizontal_dimension 0.0_kind_phys co2diag cam_out_co2diag @@ -1077,7 +1093,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="CO2PROG"> horizontal_dimension 0.0_kind_phys co2prog cam_out_co2prog @@ -1085,7 +1102,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="DSTDRY1"> horizontal_dimension 0.0_kind_phys dstdry1 cam_out_dstdry1 @@ -1093,7 +1111,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="DSTDRY2"> horizontal_dimension 0.0_kind_phys dstdry2 cam_out_dstdry2 @@ -1101,7 +1120,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="DSTDRY3"> horizontal_dimension 0.0_kind_phys dstdry3 cam_out_dstdry3 @@ -1109,7 +1129,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="DSTDRY4"> horizontal_dimension 0.0_kind_phys dstdry4 cam_out_dstdry4 @@ -1149,7 +1170,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="FLWDS"> horizontal_dimension 0.0_kind_phys flwds cam_out_flwds @@ -1189,7 +1211,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="OCPHIDRY"> horizontal_dimension 0.0_kind_phys ocphidry cam_out_ocphidry @@ -1205,7 +1228,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="OCPHODRY"> horizontal_dimension 0.0_kind_phys ocphodry cam_out_ocphodry @@ -1286,7 +1310,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="SOLL"> horizontal_dimension 0.0_kind_phys soll cam_out_soll @@ -1294,7 +1319,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="SOLLD"> horizontal_dimension 0.0_kind_phys solld cam_out_solld @@ -1302,7 +1328,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="SOLS"> horizontal_dimension 0.0_kind_phys sols cam_out_sols @@ -1310,7 +1337,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="SOLSD"> horizontal_dimension 0.0_kind_phys solsd cam_out_solsd @@ -1462,7 +1490,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="PBLH"> horizontal_dimension @@ -1472,7 +1501,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="TPERT"> horizontal_dimension 273.15_kind_phys tpert pbuf_tpert @@ -1558,7 +1588,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="QCWAT"> horizontal_dimension vertical_layer_dimension -1.0_kind_phys @@ -1567,7 +1598,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="TCWAT"> horizontal_dimension vertical_layer_dimension -1.0_kind_phys @@ -1576,7 +1608,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="LCWAT"> horizontal_dimension vertical_layer_dimension -1.0_kind_phys @@ -1604,7 +1637,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="TAURESX"> horizontal_dimension 0.0_kind_phys pbuf_tauresx @@ -1612,7 +1646,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="TAURESY"> horizontal_dimension 0.0_kind_phys pbuf_tauresy @@ -1620,7 +1655,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="DRAGBLJ"> horizontal_dimension vertical_layer_dimension 0.0_kind_phys pbuf_dragblj @@ -1633,7 +1669,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="TAUBLJX"> horizontal_dimension 0.0_kind_phys pbuf_taubljx @@ -1641,7 +1678,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="TAUBLJY"> horizontal_dimension 0.0_kind_phys pbuf_taubljy @@ -1686,7 +1724,8 @@ + allocatable="allocatable" restart="true" + diagnostic_name="CLUBBTOP"> vertical layer index for top of CLUBB for running HB above horizontal_dimension @@ -1792,14 +1831,16 @@ + allocatable="allocatable" restart="true" + diagnostic_name="QL"> horizontal_dimension vertical_layer_dimension ICWMRDP pbuf_ICWMRDP + allocatable="allocatable" restart="true" + diagnostic_name="CLDFRC"> horizontal_dimension vertical_layer_dimension 0.7_kind_phys CLD pbuf_CLD @@ -1843,21 +1884,24 @@ + allocatable="allocatable" restart="true" + diagnostic_name="QINI"> horizontal_dimension vertical_layer_dimension QINI pbuf_QINI + allocatable="allocatable" restart="true" + diagnostic_name="LIQINI"> horizontal_dimension vertical_layer_dimension TOTLIQINI pbuf_TOTLIQINI + allocatable="allocatable" restart="true" + diagnostic_name="ICEINI"> horizontal_dimension vertical_layer_dimension TOTICEINI pbuf_TOTICEINI diff --git a/src/data/registry_v1_0.xsd b/src/data/registry_v1_0.xsd index a4772445d..8cd8bef61 100644 --- a/src/data/registry_v1_0.xsd +++ b/src/data/registry_v1_0.xsd @@ -110,10 +110,14 @@ + + @@ -147,8 +151,10 @@ + + @@ -171,8 +177,10 @@ + + diff --git a/src/data/write_restart_physics.py b/src/data/write_restart_physics.py new file mode 100644 index 000000000..7b2c4a772 --- /dev/null +++ b/src/data/write_restart_physics.py @@ -0,0 +1,556 @@ +""" +Use variable meta-data from "generate_registry_data.py" +to generate a CAM fortran file to manage physics restart +information and writing (restart_physics.F90). +""" + +# Python library import statements: +from collections import OrderedDict +import os.path + +# CCPP Framework import statements +from ccpp_state_machine import CCPP_STATE_MACH +from fortran_tools import FortranWriter +# Exclude these standard names from init processing +# Some are internal names (e.g., suite_name) +# Some are from the CCPP framework (e.g., ccpp_num_constituents) +# Some are for efficiency and to avoid dependency loops (e.g., log_output_unit) +_EXCLUDED_STDNAMES = {'suite_name', 'suite_part', + 'number_of_ccpp_constituents', + 'number_of_ccpp_advected_constituents', + 'ccpp_constituents', + 'ccpp_constituent_tendencies', + 'ccpp_constituent_properties', + 'ccpp_constituent_minimum_values', + 'ccpp_error_message', + 'ccpp_error_code', + 'log_output_unit', 'do_log_output', + 'mpi_communicator', 'mpi_root', 'mpi_rank', + 'number_of_mpi_tasks'} +# Variable input types +_INPUT_TYPES = set(['in', 'inout']) +_OUTPUT_TYPES = set(['out', 'inout']) + +# Increase allowed line lengths needed to fit extra-long CCPP standard names: +_LINE_FILL_LEN = 150 +_MAX_LINE_LEN = 200 + +############## +#Main function +############## + +def write_restart_physics(cap_database, ic_names, registry_constituents, + restart_vars, outdir, file_find_func, source_paths, indent, logger, + phys_restart_filename=None): + + """ + Create restart_physics.F90 using a database created + by the CCPP Framework generator (capgen). + """ + + #Initialize return message: + retmsg = "" + + # Gather all the host model variables that are required by + # any of the compiled CCPP physics suites. + in_vars, out_vars, constituent_set, retmsg = gather_ccpp_req_vars(cap_database, registry_constituents) + + # Quit now if there are missing variables + if retmsg: + return retmsg + # end if + + # ----------------------------------------- + # Generate "restart_physics.F90" file: + # ----------------------------------------- + + # Open new file: + if phys_restart_filename: + ofilename = os.path.join(outdir, phys_restart_filename) + # Get file name, ignoring file type: + phys_restart_fname_str = os.path.splitext(phys_restart_filename)[0] + else: + ofilename = os.path.join(outdir, "restart_physics.F90") + phys_restart_fname_str = "restart_physics" + # end if + + # Log file creation: + logger.info(f"Writing restart physics source file, {ofilename}") + + # Open file using CCPP's FortranWriter: + file_desc = "physics restart source file" + with FortranWriter(ofilename, "w", file_desc, + phys_restart_fname_str, + line_fill=_LINE_FILL_LEN, + line_max=_MAX_LINE_LEN, + indent=indent) as outfile: + + # Add use statements + outfile.write("use pio, only: var_desc_t", 1) + outfile.blank_line() + # Add boilerplate code: + outfile.write_preamble() + + # Add public function declarations: + outfile.write("!! public interfaces", 0) + outfile.write("public :: restart_physics_init", 1) + outfile.write("public :: restart_physics_write", 1) + outfile.write("public :: restart_physics_read", 1) + + all_req_vars = in_vars + out_vars + # Grab the required variables that are also restart variables + host_dict = cap_database.host_model_dict() + required_restart_vars, constituent_dimmed_vars, used_vars = gather_required_restart_variables(all_req_vars, constituent_set, restart_vars, host_dict) + + outfile.blank_line() + outfile.comment("Private module data", 0) + for key, value in required_restart_vars.items(): + outfile.write(f"type(var_desc_t) :: {value['diag_name'].lower()}_desc", 1) + # end for + + for key, value in constituent_dimmed_vars.items(): + outfile.write(f"type(var_desc_t), allocatable :: {value['diag_name'].lower()}_desc(:)", 1) + # end for + + # Add "contains" statement: + outfile.end_module_header() + + # Write the restart init subroutine + outfile.blank_line() + dim_use_stmts, dims = write_restart_physics_init(outfile, required_restart_vars, constituent_dimmed_vars, host_dict) + + # Write the restart write subroutine + outfile.blank_line() + write_restart_physics_write(outfile, required_restart_vars, constituent_dimmed_vars, used_vars, dim_use_stmts, dims) + + # Write the restart read subroutine + outfile.blank_line() + write_restart_physics_read(outfile, required_restart_vars, constituent_dimmed_vars, used_vars) + + # end with + +def write_restart_physics_init(outfile, required_vars, constituent_dimmed_vars, host_dict): + """ + Write the init routine for the physics restart variables. This + routine initializes the physics fields on the restart (cam.r) file + """ + + outfile.write("subroutine restart_physics_init(file, errmsg, errflg)", 1) + + dimensions_dict = {} + num_dimensions = 0 + dim_use_stmt_dict = {} + + # Find all unique dimensions + for _, value in required_vars.items(): + for dimension in value['dims']: + dim_name = dimension.split(':')[1] + if dim_name not in dimensions_dict.keys(): + var = host_dict.find_variable(dim_name) + dimensions_dict[dim_name] = {'local_name': var.get_prop_value('local_name'), 'index': -1} + num_dimensions = num_dimensions + 1 + # Add dimension to use statement dictionary + if var.source.name not in dim_use_stmt_dict.keys(): + dim_use_stmt_dict[var.source.name] = set() + # end if + dim_use_stmt_dict[var.source.name].add(var.get_prop_value('local_name')) + # end if + # end for + # end for + for _, value in constituent_dimmed_vars.items(): + for dimension in value['dims']: + dim_name = dimension.split(':')[1] + if dim_name not in dimensions_dict.keys() and dim_name != 'number_of_ccpp_constituents': + var = host_dict.find_variable(dim_name) + dimensions_dict[dim_name] = {'local_name': var.get_prop_value('local_name'), 'index': -1} + num_dimensions = num_dimensions + 1 + # Add dimension to use statement dictionary + if var.source.name not in dim_use_stmt_dict.keys(): + dim_use_stmt_dict[var.source.name] = set() + # end if + dim_use_stmt_dict[var.source.name].add(var.get_prop_value('local_name')) + # end if + # end for + # end for + + static_use_stmts = [["pio", ["file_desc_t", "pio_double"]], + ["cam_pio_utils", ["cam_pio_def_dim", "cam_pio_def_var"]], + ["cam_ccpp_cap", ["cam_model_const_properties"]], + ["physics_grid", ["num_global_phys_cols"]], + ["ccpp_constituent_prop_mod", ["ccpp_constituent_prop_ptr_t"]]] + + # Gather up dimension imports + dim_use_stmts = [] + dims = [] + for key, value in dim_use_stmt_dict.items(): + imports = [] + dim_use_stmt = [] + for var_import in value: + imports.append(f"{var_import}") + dims.append(f"{var_import}") + # end for + dim_use_stmt.append(key) + dim_use_stmt.append(imports) + dim_use_stmts.append(dim_use_stmt) + # end for + + # Output required, registered Fortran module use statements: + write_use_statements(outfile, static_use_stmts, 2) + write_use_statements(outfile, dim_use_stmts, 2) + + outfile.write("type(file_desc_t), intent(inout) :: file", 2) + outfile.write("character(len=512),intent(out) :: errmsg", 2) + outfile.write("integer, intent(out) :: errflg", 2) + outfile.blank_line() + + outfile.comment("Local variables", 2) + outfile.write("integer, allocatable :: dimids(:)", 2) + outfile.write("integer :: constituent_idx", 2) + outfile.write("type(ccpp_constituent_prop_ptr_t), pointer :: const_props(:)", 2) + outfile.write("character(len=256) :: const_diag_name", 2) + + outfile.blank_line() + + outfile.blank_line() + outfile.comment("Allocate dimids to the number of unique dimensions", 2) + outfile.write(f"allocate(dimids({num_dimensions}), stat=errflg, errmsg=errmsg)", 2) + outfile.write("if (errflg /= 0) then", 2) + outfile.write("return", 3) + outfile.write("end if", 2) + + dim_index = 1 + + outfile.comment("Define required restart variables on the restart file", 2) + for key, value in required_vars.items(): + hdimids = [] + for dimension in value['dims']: + dimname = dimension.split(':')[1] + if dimensions_dict[dimname]['index'] < 0: + outfile.comment(f"Define potentially new dimension '{dimname}'", 2) + if dimname == 'horizontal_dimension': + dim_loc_name = 'ncol' + dimsize = 'num_global_phys_cols' + elif dimname == 'vertical_layer_dimension': + dim_loc_name = 'lev' + dimsize = 'pver' + elif dimname == 'vertical_interface_dimension': + dim_loc_name = 'ilev' + dimsize = 'pverp' + else: + dim_loc_name = dimensions_dict[dimname]['local_name'] + outfile.write(f"call cam_pio_def_dim(file, '{dim_loc_name}', {dimsize}, dimids({dim_index}), existOK=.true.)", 2) # grab use statements for phys variables for nonstandard dimensions! + dimensions_dict[dimname]['index'] = dim_index + dim_index = dim_index + 1 + # end if + hdimids.append(dimensions_dict[dimname]['index']) + # end for + desc_name = f"{value['diag_name'].lower()}_desc" + dimids_array = ", ".join([f"dimids({i})" for i in hdimids]) + outfile.write(f"call cam_pio_def_var(file, '{value['diag_name']}', pio_double, (/{dimids_array}/), {desc_name}, existOK=.false.)", 2) + outfile.write("if (errflg /= 0) then", 2) + outfile.write(f"write(errmsg,*) 'restart_physics_init: error defining variable {value['diag_name']}'", 3) + outfile.write("return", 3) + outfile.write("end if", 2) + outfile.blank_line() + # end for + + # Handle constituent-dimensioned variables + if constituent_dimmed_vars: + outfile.write("const_props => cam_model_const_properties()", 2) + outfile.blank_line() + for key, value in constituent_dimmed_vars.items(): + hdimids = [] + for dimension in value['dims']: + dimname = dimension.split(':')[1] + if dimname != 'number_of_ccpp_constituents': + if dimensions_dict[dimname]['index'] < 0: + outfile.comment(f"Define potentially new dimension '{dimname}'", 2) + if dimname == 'horizontal_dimension': + dim_loc_name = 'ncol' + dimsize = 'num_global_phys_cols' + elif dimname == 'vertical_layer_dimension': + dim_loc_name = 'lev' + dimsize = 'pver' + elif dimname == 'vertical_interface_dimension': + dim_loc_name = 'ilev' + dimsize = 'pverp' + else: + dim_loc_name = dimensions_dict[dimname]['local_name'] + # end if + outfile.write(f"call cam_pio_def_dim(file, '{dim_loc_name}', {dimsize}, dimids({dim_index}), existOK=.true.)", 2) # grab use statements for phys variables for nonstandard dimensions! + dimensions_dict[dimname]['index'] = dim_index + dim_index = dim_index + 1 + # end if + hdimids.append(dimensions_dict[dimname]['index']) + # end if (ignore number of constituents dimensions) + # end for + outfile.comment(f"Handling for constituent-dimensioned variable '{value['diag_name']}'", 2) + outfile.write(f"allocate({value['diag_name'].lower()}_desc(size(const_props)), stat=errflg, errmsg=errmsg)", 2) + outfile.write("if (errflg /= 0) then", 2) + outfile.write("return", 3) + outfile.write("end if", 2) + outfile.write("do constituent_idx = 1, size(const_props)", 2) + outfile.comment("Grab constituent diagnostic name:", 3) + outfile.write("call const_props(constituent_idx)%diagnostic_name(const_diag_name)", 3) + desc_name = f"{value['diag_name'].lower()}_desc" + dimids_array = ", ".join([f"dimids({i})" for i in hdimids]) + outfile.write(f"call cam_pio_def_var(file, '{value['diag_name']}_'//trim(const_diag_name), pio_double, (/{dimids_array}/), {desc_name}(constituent_idx), existOK=.false.)", 3) + outfile.write("end do", 2) + outfile.blank_line() + # end for + # end if + + outfile.write("end subroutine restart_physics_init", 1) + return dim_use_stmts, dims + +def write_restart_physics_write(outfile, required_vars, constituent_dimmed_vars, used_vars, dim_use_stmts, dims): + """ + Write the 'write' routine for the physics restart variables. This + routine writes the physics fields to the restart (cam.r) file + """ + outfile.write("subroutine restart_physics_write(file, grid_id, errmsg, errflg)", 1) + + use_stmts = [["pio", ["file_desc_t", "io_desc_t", "pio_write_darray", "pio_double"]], + ["cam_ccpp_cap", ["cam_model_const_properties","cam_constituents_array"]], + ["ccpp_kinds", ["kind_phys"]], + ["ccpp_constituent_prop_mod", ["ccpp_constituent_prop_ptr_t"]], + ["physics_grid", ["num_global_phys_cols"]], + ["cam_grid_support", ["cam_grid_id", "cam_grid_write_dist_array"]]] + + write_use_statements(outfile, use_stmts, 2) + write_use_statements(outfile, dim_use_stmts, 2) + + for var in used_vars: + outfile.write(f"use physics_types, only: {var}", 2) + # end for + + outfile.blank_line() + outfile.write("type(file_desc_t), intent(inout) :: file", 2) + outfile.write("integer, intent(in) :: grid_id", 2) + outfile.write("character(len=512),intent(out) :: errmsg", 2) + outfile.write("integer, intent(out) :: errflg", 2) + + outfile.blank_line() + + outfile.comment("Local variables", 2) + outfile.write(f"integer :: dims({len(dims)})", 2) + outfile.write("integer :: grid_decomp", 2) + outfile.write("integer :: grid_dims(2)", 2) + outfile.write("integer :: field_shape(2)", 2) + outfile.write("integer :: constituent_idx", 2) + outfile.write("type(ccpp_constituent_prop_ptr_t), pointer :: const_props(:)", 2) + + # Just exit if we don't have any variables! + if len(required_vars) == 0 and len(constituent_dimmed_vars) == 0: + outfile.write("end subroutine restart_physics_write", 1) + return + # end if + + outfile.comment("Grab physics grid", 2) + outfile.write("grid_decomp = cam_grid_id('physgrid')", 2) + dim_index = 1 + for dim in dims: + outfile.write(f"dims({dim_index}) = {dim}", 2) + dim_index = dim_index + 1 + # end if + outfile.comment("Write required restart variables to the restart file", 2) + for key, value in required_vars.items(): + desc_name = f"{value['diag_name'].lower()}_desc" + if len(value['dims']) == 1 and 'horizontal_dimension' in value['dims'][0]: + outfile.comment("Handle horizontal-only field", 2) + outfile.write("field_shape(1) = num_global_phys_cols", 2) + outfile.write(f"call cam_grid_write_dist_array(file, grid_decomp, (/dims(1)/), (/field_shape(1)/), {key}, {desc_name})", 2) + elif len(value["dims"]) > 1: + for index, dim in enumerate(value["dims"]): + if 'horizontal_dimension' in dim: + outfile.write("field_shape(1) = num_global_phys_cols", 2) + else: + outfile.write(f"field_shape({index + 1}) = size({key}, {index + 1})", 2) + # end if + # end if + outfile.write(f"call cam_grid_write_dist_array(file, grid_decomp, dims, field_shape, {key}, {desc_name})", 2) + # end if + # end for + + # Handle constituent-dimensioned variables + if constituent_dimmed_vars: + outfile.write("const_props => cam_model_const_properties()", 2) + outfile.blank_line() + for key, value in constituent_dimmed_vars.items(): + outfile.comment(f"Handling for constituent-dimensioned variable '{value['diag_name']}'", 2) + outfile.write("do constituent_idx = 1, size(const_props)", 2) + desc_name = f"{value['diag_name'].lower()}_desc" + outfile.write("field_shape(1) = num_global_phys_cols", 3) + outfile.write(f"call cam_grid_write_dist_array(file, grid_decomp, (/dims(1)/), (/field_shape(1)/), {key}(:,constituent_idx), {desc_name}(constituent_idx))", 3) + outfile.write("end do", 2) + outfile.blank_line() + # end for + # end if + outfile.write("end subroutine restart_physics_write", 1) + +def write_restart_physics_read(outfile, required_vars, constituent_dimmed_vars, used_vars): + outfile.write("subroutine restart_physics_read()", 1) + outfile.write("end subroutine restart_physics_read", 1) + +################# +#HELPER FUNCTIONS +################# + +############################################################################## +def _find_and_add_host_variable(stdname, host_dict, var_dict): + """Find in and add it to if found and + not of type, 'host'. + If not found, add to . + If found and added to , also process the standard names of + any intrinsic sub-elements of . + Return the list of (if any). + Note: This function has a side effect (adding to ). + """ + missing_vars = [] + hvar = host_dict.find_variable(stdname) + if hvar and (hvar.source.ptype != 'host'): + var_dict[stdname] = hvar + # Process elements (if any) + ielem = hvar.intrinsic_elements() + # List elements are the only ones we care about + if isinstance(ielem, list): + for sname in ielem: + smissing = _find_and_add_host_variable(sname, host_dict, + var_dict) + missing_vars.extend(smissing) + # end for + # end if + # end if + if not hvar: + missing_vars.append(stdname) + # end if + return missing_vars + +############################################################################## +def gather_ccpp_req_vars(cap_database, registry_constituents): + """ + Generate a list of host-model and constituent variables + required by the CCPP physics suites potentially being used + in this model run. + is the database object returned by capgen. + It is an error if any physics suite variable is not accessible in + the host model. + Return several values: + - A list of host model variables + - An error message (blank for no error) + """ + + # Dictionary of all 'in' and 'inout' suite variables. + # Key is standard name, value is host-model or constituent variable + in_vars = {} + out_vars = {} + missing_vars = set() + constituent_vars = set() + retmsg = "" + # Host model dictionary + host_dict = cap_database.host_model_dict() + + # Create CCPP datatable required variables-listing object: + # XXgoldyXX: Choose only some phases here? + for phase in CCPP_STATE_MACH.transitions(): + for cvar in cap_database.call_list(phase).variable_list(): + stdname = cvar.get_prop_value('standard_name') + intent = cvar.get_prop_value('intent') + is_const = cvar.get_prop_value('advected') or cvar.get_prop_value('constituent') + if ((intent in _INPUT_TYPES) and + (stdname not in in_vars) and + (stdname not in _EXCLUDED_STDNAMES)): + if is_const: + #Add variable to constituent set: + constituent_vars.add(stdname) + #Add variable to required variable list if it's not a registry constituent + if stdname not in registry_constituents: + in_vars[stdname] = cvar + # end if + else: + # We need to work with the host model version of this variable + missing = _find_and_add_host_variable(stdname, host_dict, + in_vars) + missing_vars.update(missing) + # end if + # end if (only input variables) + if ((intent in _OUTPUT_TYPES) and + (stdname not in out_vars) and + (stdname not in _EXCLUDED_STDNAMES)): + if not is_const: + missing = _find_and_add_host_variable(stdname, host_dict, + out_vars) + # do nothing with missing variables + # end if + # end if (only output variables) + # end for (loop over call list) + # end for (loop over phases) + + if missing_vars: + mvlist = ', '.join(sorted(missing_vars)) + retmsg = f"Error: Missing required host variables: {mvlist}" + # end if + # Return the required variables as a list + return list(in_vars.values()), list(out_vars.values()), constituent_vars, retmsg + +############################################################################## +def gather_required_restart_variables(all_req_vars, registry_constituents, restart_vars, host_dict): + """ + Return lists of the required restart non-constituent variables, and the required + restart constituent variables + """ + + required_restart_vars = {} + required_constituent_dimensioned_vars = {} + used_vars = set() + + for required_var in all_req_vars: + stdname = required_var.get_prop_value('standard_name') + local_name = required_var.call_string(host_dict) + used_var = required_var.var.get_prop_value('local_name') + if stdname in restart_vars.keys(): + diagnostic_name = restart_vars[stdname] + used_vars.add(used_var) + dimensions = required_var.get_dimensions() + if 'ccpp_constant_one:number_of_ccpp_constituents' in dimensions: + required_constituent_dimensioned_vars[local_name] = {'diag_name': diagnostic_name, 'dims': dimensions} + else: + required_restart_vars[local_name] = {'diag_name': diagnostic_name, 'dims': dimensions} + # end if + # end if (ignore all non-restart variables) + # end for + + return required_restart_vars, required_constituent_dimensioned_vars, used_vars + +##### + +def write_use_statements(outfile, use_stmts, indent): + """Output Fortran module use (import) statements listed in . + """ + # Don't do anything if we don't have use statements! + if len(use_stmts) == 0: + return + # end if + # The plus one is for a comma + max_modname = max(len(x[0]) for x in use_stmts) + 1 + # max_modspace is the max chars of the module plus other 'use' statement + # syntax (e.g., 'only:') + max_modspace = (outfile.indent_size * indent) + max_modname + 10 + mod_space = outfile.line_fill - max_modspace + for use_item in use_stmts: + # Break up imported interfaces to clean up use statements + larg = 0 + num_imports = len(use_item[1]) + while larg < num_imports: + int_str = use_item[1][larg] + larg = larg + 1 + while ((larg < num_imports) and + ((len(int_str) + len(use_item[1][larg]) + 2) < mod_space)): + int_str += f", {use_item[1][larg]}" + larg = larg + 1 + # end while + modname = use_item[0] + ',' + outfile.write(f"use {modname: <{max_modname}} only: {int_str}", + indent) + # end while + # end for diff --git a/src/dynamics/mpas/restart_dynamics.F90 b/src/dynamics/mpas/restart_dynamics.F90 new file mode 100644 index 000000000..1babfac8a --- /dev/null +++ b/src/dynamics/mpas/restart_dynamics.F90 @@ -0,0 +1,50 @@ +module restart_dynamics + +! Writing and reading grid and dynamics state information to/from restart files is +! delegated to MPAS utility code. This module provides the CAM interfaces for the +! restart functionality. CAM just provides MPAS with the PIO filehandle to the +! restart file. + +use dyn_comp, only: dyn_export_t, mpas_dynamical_core +use pio, only: file_desc_t + +implicit none +private +save + +public :: & + init_restart_dynamics, & + write_restart_dynamics + +!========================================================================================= +contains +!========================================================================================= + +subroutine init_restart_dynamics(file, dyn_out) + + ! arguments + type(file_desc_t), target :: File + type(dyn_export_t), intent(in) :: dyn_out + !---------------------------------------------------------------------------- + + ! STUB-ROUTINE + +end subroutine init_restart_dynamics + +!========================================================================================= + +subroutine write_restart_dynamics(File, dyn_out) + + ! arguments + type(File_desc_t), target :: File + type(dyn_export_t), intent(in) :: dyn_out + type(file_desc_t), pointer :: file_ptr + !---------------------------------------------------------------------------- + + file_ptr => File + + call mpas_dynamical_core % read_write_stream(File, 'w', 'invariant+restart+input') + +end subroutine write_restart_dynamics + +end module restart_dynamics diff --git a/src/dynamics/none/restart_dynamics.F90 b/src/dynamics/none/restart_dynamics.F90 new file mode 100644 index 000000000..db2a6c19b --- /dev/null +++ b/src/dynamics/none/restart_dynamics.F90 @@ -0,0 +1,41 @@ +module restart_dynamics + +use dyn_comp, only: dyn_export_t +use pio, only: file_desc_t + +implicit none +private +save + +public :: & + init_restart_dynamics, & + write_restart_dynamics + +!========================================================================================= +contains +!========================================================================================= + +subroutine init_restart_dynamics(file, dyn_out) + + ! arguments + type(file_desc_t), target :: File + type(dyn_export_t), intent(in) :: dyn_out + !---------------------------------------------------------------------------- + + ! STUB-ROUTINE + +end subroutine init_restart_dynamics + +!========================================================================================= + +subroutine write_restart_dynamics(File, dyn_out) + + ! arguments + type(File_desc_t), target :: File + type(dyn_export_t), intent(in) :: dyn_out + !---------------------------------------------------------------------------- + + ! STUB-ROUTINE +end subroutine write_restart_dynamics + +end module restart_dynamics diff --git a/src/dynamics/se/restart_dynamics.F90 b/src/dynamics/se/restart_dynamics.F90 new file mode 100644 index 000000000..3704c1a54 --- /dev/null +++ b/src/dynamics/se/restart_dynamics.F90 @@ -0,0 +1,563 @@ +module restart_dynamics + +! Write and read dynamics fields from the restart file. For exact restart +! it is necessary to write all element data, including duplicate columns, +! to the file. The namelist option, se_write_restart_unstruct, is +! available to write just the unique columns to the restart file using the +! same unstructured grid used by the history and initial files. This +! results in the introduction of a roundoff size difference on restart, but +! writes the fields in the unstructured grid format which is easier to +! modify if the user desires to introduce perturbations or other +! adjustments into the run. The restart file containing the unstructured +! grid format may also be used for an initial run. + +use pio, only: var_desc_t + +implicit none +private +save + +public :: init_restart_dynamics +public :: write_restart_dynamics +!public :: read_restart_dynamics + +! these variables are module data so they can be shared between the +! file definition and write phases +type(var_desc_t) :: psdry_desc, udesc, vdesc, tdesc +type(var_desc_t), allocatable :: qdesc_dp(:) +type(var_desc_t) :: dp_fvm_desc +type(var_desc_t), pointer :: c_fvm_desc(:) + +integer, private :: nelem_tot = -1 ! Correct total number of elements + +!========================================================================================= +CONTAINS +!========================================================================================= + +subroutine init_nelem_tot() + use spmd_utils, only: mpicom + use mpi, only: mpi_integer, mpi_sum + use dimensions_mod, only: nelemd + + integer :: ierr + + if (nelem_tot < 0) then + call MPI_Allreduce(nelemd, nelem_tot, 1, MPI_INTEGER, MPI_SUM, mpicom, ierr) + end if +end subroutine init_nelem_tot + +subroutine init_restart_dynamics(file, dyn_out) + use hycoef, only: init_restart_hycoef + use cam_ccpp_cap, only: cam_model_const_properties + use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t + use pio, only: file_desc_t, pio_global, pio_unlimited, pio_double, pio_seterrorhandling + use pio, only: pio_bcast_error, pio_def_dim, pio_def_var, pio_put_att + use dyn_comp, only: dyn_export_t, write_restart_unstruct + use cam_grid_support, only: cam_grid_header_info_t, cam_grid_id + use cam_grid_support, only: cam_grid_write_attr + use dimensions_mod, only: np, npsq, ne, nelemd, fv_nphys + ! Define dimensions, variables, attributes for restart file. + + ! This is not really an "init" routine. It is called before + ! write_restart_dynamics every time an restart is written. + + ! arguments + type(file_desc_t), intent(inout) :: file + type(dyn_export_t), intent(in) :: dyn_out + + ! local variables + integer :: i + integer :: vdimids(2) + integer :: nlev_dimid + integer :: ncol_dimid + integer :: ncol_fvm_dimid + integer :: time_dimid + + integer :: ierr, err_handling + + integer :: grid_id + type(cam_grid_header_info_t) :: info + integer :: constituent_idx + type(ccpp_constituent_prop_ptr_t), pointer :: const_props(:) + character(len=256) :: const_diag_name + + !---------------------------------------------------------------------------- + + call init_nelem_tot() + call init_restart_hycoef(file, vdimids) + nlev_dimid = vdimids(1) + + call pio_seterrorhandling(File, pio_bcast_error, err_handling) + + ierr = PIO_Def_Dim(File, 'time', PIO_UNLIMITED, time_dimid) + + ! GLL restart fields + + ! number of columns written to restart depends on whether all columns in the + ! element structures are written, or just the unique columns (unstructured grid) + if (write_restart_unstruct) then + grid_id = cam_grid_id('GLL') + call cam_grid_write_attr(File, grid_id, info) + ncol_dimid = info%get_hdimid(1) + else + ierr = PIO_Def_Dim(File,'nenpnp', nelem_tot*np*np, ncol_dimid) + ierr = PIO_Put_Att(File, PIO_GLOBAL, 'ne', ne) + ierr = PIO_Put_Att(File, PIO_GLOBAL, 'np', np) + end if + + ierr = PIO_Def_Var(File, 'PSDRY', pio_double, (/ncol_dimid, time_dimid/), psdry_desc) + ierr = PIO_Def_Var(File, 'U', pio_double, (/ncol_dimid, nlev_dimid, time_dimid/), Udesc) + ierr = PIO_Def_Var(File, 'V', pio_double, (/ncol_dimid, nlev_dimid, time_dimid/), Vdesc) + ierr = PIO_Def_Var(File, 'T', pio_double, (/ncol_dimid, nlev_dimid, time_dimid/), Tdesc) + + const_props => cam_model_const_properties() + allocate(qdesc_dp(size(const_props))) + do constituent_idx = 1, size(const_props) + ! Grab constituent diagnostic name: + call const_props(constituent_idx)%diagnostic_name(const_diag_name) + ierr = PIO_Def_Var(File,"dp"//trim(const_diag_name), pio_double, & + (/ncol_dimid, nlev_dimid, time_dimid/), Qdesc_dp(constituent_idx)) + end do + + ! CSLAM restart fields + + if (fv_nphys > 0) then + + grid_id = cam_grid_id('FVM') + call cam_grid_write_attr(File, grid_id, info) + ncol_fvm_dimid = info%get_hdimid(1) + + ierr = PIO_Def_Var(File, 'dp_fvm', pio_double, & + (/ncol_fvm_dimid, nlev_dimid, time_dimid/), dp_fvm_desc) + + allocate(c_fvm_desc(size(const_props))) + do constituent_idx = 1, size(const_props) + call const_props(constituent_idx)%diagnostic_name(const_diag_name) + ierr = PIO_Def_Var(File, trim(const_diag_name)//"_fvm", pio_double, & + (/ncol_fvm_dimid, nlev_dimid, time_dimid/), c_fvm_desc(constituent_idx)) + end do + + end if + + call pio_seterrorhandling(File, err_handling) + +end subroutine init_restart_dynamics + +!========================================================================================= + +subroutine write_restart_dynamics(File, dyn_out) + use control_mod, only: qsplit + use pio, only: file_desc_t, pio_offset_kind, io_desc_t, pio_double + use pio, only: pio_initdecomp, pio_freedecomp, pio_setframe, pio_write_darray + use dyn_comp, only: dyn_export_t, write_restart_unstruct + use dyn_grid, only: timelevel + use cam_grid_support, only: cam_grid_id, cam_grid_write_var, cam_grid_get_decomp, cam_grid_dimensions + use element_mod, only: element_t + use fvm_control_volume_mod, only: fvm_struct + use hycoef, only: write_restart_hycoef + use time_mod, only: TimeLevel_Qdp + use parallel_mod, only: par + use dimensions_mod, only: nc, np, npsq, ne, nelemd, fv_nphys, nlev + use cam_ccpp_cap, only: cam_model_const_properties + use ccpp_constituent_prop_mod, only: ccpp_constituent_prop_ptr_t + use cam_pio_utils, only: pio_subsystem + use thread_mod, only: horz_num_threads + use spmd_utils, only: iam + use shr_kind_mod, only: r8 => shr_kind_r8 + + type(file_desc_t), intent(inout) :: File + type(dyn_export_t), intent(in) :: dyn_out + + ! local variables + integer(pio_offset_kind), parameter :: t_idx = 1 + + type(element_t), pointer :: elem(:) + type(fvm_struct), pointer :: fvm(:) + + integer :: tl, tlqdp + integer :: i, ie, ii, j, k, m + integer :: ierr + + integer :: grid_id + integer :: grid_dimlens(2) + + + + integer :: array_lens(3) + integer :: file_lens(2) + type(io_desc_t), pointer :: iodesc3d_fvm + real(r8), allocatable :: buf3d(:,:,:) + type(ccpp_constituent_prop_ptr_t), pointer :: const_props(:) + + + + character(len=*), parameter :: sub = 'write_restart_dynamics' + !---------------------------------------------------------------------------- + + call write_restart_hycoef(File) + + tl = timelevel%n0 + call TimeLevel_Qdp(timelevel, qsplit, tlQdp) + const_props => cam_model_const_properties() + + if (iam .lt. par%nprocs) then + elem => dyn_out%elem + fvm => dyn_out%fvm + else + allocate (elem(0), fvm(0)) + endif + + ! write fields on GLL grid + + if (write_restart_unstruct) then + call write_unstruct() + else + call write_elem() + end if + + ! write CSLAM fields + + if (fv_nphys > 0) then + + grid_id = cam_grid_id('FVM') + + ! write coords for FVM grid + call cam_grid_write_var(File, grid_id) + + call cam_grid_dimensions(grid_id, grid_dimlens) + allocate(buf3d(nc*nc,nlev,nelemd)) + array_lens = (/nc*nc, nlev, nelemd/) + file_lens = (/grid_dimlens(1), nlev/) + call cam_grid_get_decomp(grid_id, array_lens, file_lens, pio_double, iodesc3d_fvm) + + do ie = 1, nelemd + do k = 1, nlev + ii = 1 + do j = 1, nc + do i = 1, nc + buf3d(ii,k,ie) = fvm(ie)%dp_fvm(i,j,k) + ii = ii + 1 + end do + end do + end do + end do + call PIO_Setframe(file, dp_fvm_desc, t_idx) + call PIO_Write_Darray(file, dp_fvm_desc, iodesc3d_fvm, buf3d, ierr) + + do m = 1, size(const_props) + do ie = 1, nelemd + do k = 1, nlev + ii = 1 + do j = 1, nc + do i = 1, nc + buf3d(ii,k,ie) = fvm(ie)%c(i,j,k,m) + ii = ii + 1 + end do + end do + end do + end do + call PIO_Setframe(file, c_fvm_desc(m), t_idx) + call PIO_Write_Darray(file, c_fvm_desc(m), iodesc3d_fvm, buf3d, ierr) + end do + + deallocate(c_fvm_desc) + deallocate(buf3d) + ! should this call be made on a pointer? + !call pio_freedecomp(File, iodesc3d_fvm) + + end if + + if (iam >= par%nprocs) then + deallocate(elem, fvm) + endif + +!------------------------------------------------------------------------------- +contains +!------------------------------------------------------------------------------- + +subroutine write_elem() + + ! local variables + integer :: i, ie, j, k + integer :: ierr + integer, pointer :: ldof(:) + + type(io_desc_t) :: iodesc2d, iodesc3d + + real(kind=r8), pointer :: var3d(:,:,:,:), var2d(:,:,:) + !---------------------------------------------------------------------------- + + ldof => get_restart_decomp(elem, 1) + call PIO_InitDecomp(pio_subsystem, pio_double, (/nelem_tot*np*np/), ldof, iodesc2d) + deallocate(ldof) + + ldof => get_restart_decomp(elem, nlev) + call PIO_InitDecomp(pio_subsystem, pio_double, (/nelem_tot*np*np,nlev/), ldof, iodesc3d) + deallocate(ldof) + + allocate(var2d(np,np,nelemd)) + allocate(var3d(np,np,nelemd,nlev)) + + !$omp parallel do num_threads(horz_num_threads) private(ie, j, i) + do ie = 1, nelemd + do j = 1, np + do i = 1, np + var2d(i,j,ie) = elem(ie)%state%psdry(i,j) + end do + end do + end do + call PIO_Setframe(File, psdry_desc, t_idx) + call PIO_Write_Darray(File, psdry_desc, iodesc2d, var2d, ierr) + + !$omp parallel do num_threads(horz_num_threads) private(ie, k, j, i) + do ie = 1, nelemd + do k = 1, nlev + do j = 1, np + do i = 1, np + var3d(i,j,ie,k) = elem(ie)%state%V(i,j,1,k,tl) + end do + end do + end do + end do + call PIO_Setframe(File, Udesc, t_idx) + call PIO_Write_Darray(File, Udesc, iodesc3d, var3d, ierr) + + !$omp parallel do num_threads(horz_num_threads) private(ie, k, j, i) + do ie = 1, nelemd + do k = 1, nlev + do j = 1, np + do i = 1, np + var3d(i,j,ie,k) = elem(ie)%state%V(i,j,2,k,tl) + end do + end do + end do + end do + call PIO_Setframe(File, Vdesc, t_idx) + call PIO_Write_Darray(File, Vdesc, iodesc3d, var3d, ierr) + + !$omp parallel do num_threads(horz_num_threads) private(ie, k, j, i) + do ie = 1, nelemd + do k = 1, nlev + do j = 1, np + do i = 1, np + var3d(i,j,ie,k) = elem(ie)%state%T(i,j,k,tl) + end do + end do + end do + end do + call PIO_Setframe(File, Tdesc, t_idx) + call PIO_Write_Darray(File, Tdesc, iodesc3d, var3d, ierr) + + do m = 1, size(const_props) + + !$omp parallel do num_threads(horz_num_threads) private(ie, k, j, i) + do ie = 1, nelemd + do k = 1, nlev + do j = 1, np + do i = 1, np + var3d(i,j,ie,k) = elem(ie)%state%Qdp(i,j,k,m,tlQdp) + end do + end do + end do + end do + call PIO_Setframe(File, Qdesc_dp(m), t_idx) + call PIO_Write_Darray(File, Qdesc_dp(m), iodesc3d, var3d, ierr) + + end do + + deallocate(var2d) + deallocate(var3d) + deallocate(qdesc_dp) + + call pio_freedecomp(File, iodesc2d) + call pio_freedecomp(File, iodesc3d) + +end subroutine write_elem + +!------------------------------------------------------------------------------- + +subroutine write_unstruct() + + ! local variables + integer :: i, ie, ii, j, k + integer :: ierr + + integer :: array_lens_3d(3), array_lens_2d(2) + integer :: file_lens_2d(2), file_lens_1d(1) + + type(io_desc_t), pointer :: iodesc + real(r8), allocatable :: var2d(:,:), var3d(:,:,:) + !---------------------------------------------------------------------------- + + grid_id = cam_grid_id('GLL') + + ! write coordinate variables for unstructured GLL grid + call cam_grid_write_var(File, grid_id) + + ! create map for distributed write + call cam_grid_dimensions(grid_id, grid_dimlens) + + ! create map for distributed write of 2D fields + array_lens_2d = (/npsq, nelemd/) + file_lens_1d = (/grid_dimlens(1)/) + call cam_grid_get_decomp(grid_id, array_lens_2d, file_lens_1d, pio_double, iodesc) + + allocate(var2d(npsq,nelemd)) + + do ie = 1, nelemd + ii = 1 + do j = 1, np + do i = 1, np + var2d(ii,ie) = elem(ie)%state%psdry(i,j) + ii = ii + 1 + end do + end do + end do + call PIO_Setframe(File, psdry_desc, t_idx) + call PIO_Write_Darray(File, psdry_desc, iodesc, var2d, ierr) + + nullify(iodesc) + deallocate(var2d) + + ! create map for distributed write of 3D fields + array_lens_3d = (/npsq, nlev, nelemd/) + file_lens_2d = (/grid_dimlens(1), nlev/) + call cam_grid_get_decomp(grid_id, array_lens_3d, file_lens_2d, pio_double, iodesc) + + allocate(var3d(npsq,nlev,nelemd)) + + do ie = 1, nelemd + do k = 1, nlev + ii = 1 + do j = 1, np + do i = 1, np + var3d(ii,k,ie) = elem(ie)%state%V(i,j,1,k,tl) + ii = ii + 1 + end do + end do + end do + end do + call PIO_Setframe(File, Udesc, t_idx) + call PIO_Write_Darray(File, Udesc, iodesc, var3d, ierr) + + do ie = 1, nelemd + do k = 1, nlev + ii = 1 + do j = 1, np + do i = 1, np + var3d(ii,k,ie) = elem(ie)%state%V(i,j,2,k,tl) + ii = ii + 1 + end do + end do + end do + end do + call PIO_Setframe(File, Vdesc, t_idx) + call PIO_Write_Darray(File, Vdesc, iodesc, var3d, ierr) + + do ie = 1, nelemd + do k = 1, nlev + ii = 1 + do j = 1, np + do i = 1, np + var3d(ii,k,ie) = elem(ie)%state%T(i,j,k,tl) + ii = ii + 1 + end do + end do + end do + end do + call PIO_Setframe(File, Tdesc, t_idx) + call PIO_Write_Darray(File, Tdesc, iodesc, var3d, ierr) + + do m = 1, size(const_props) + + !$omp parallel do num_threads(horz_num_threads) private(ie, k, j, i) + do ie = 1, nelemd + do k = 1, nlev + ii = 1 + do j = 1, np + do i = 1, np + var3d(ii,k,ie) = elem(ie)%state%Qdp(i,j,k,m,tlQdp) + ii = ii + 1 + end do + end do + end do + end do + call PIO_Setframe(File, Qdesc_dp(m), t_idx) + call PIO_Write_Darray(File, Qdesc_dp(m), iodesc, var3d, ierr) + + end do + + deallocate(var3d) + deallocate(qdesc_dp) + +end subroutine write_unstruct + +!------------------------------------------------------------------------------- + +end subroutine write_restart_dynamics + +!========================================================================================= + +!========================================================================================= +! Private +!========================================================================================= + +function get_restart_decomp(elem, lev) result(ldof) + use dimensions_mod, only: nelemd, np + use element_mod, only: element_t + + ! Get the integer mapping of a variable in the dynamics decomp in memory. + ! The canonical ordering is as on the file. A 0 value indicates that the + ! variable is not on the file (eg halo or boundary values) + + type(element_t), intent(in) :: elem(:) + integer, intent(in) :: lev + integer, pointer :: ldof(:) + + integer :: i, j, k, ie + !---------------------------------------------------------------------------- + + allocate(ldof(nelemd*np*np*lev)) + + j = 1 + do k = 1, lev + do ie = 1, nelemd + do i = 1, np*np + ldof(j) = (elem(ie)%GlobalID-1)*np*np + (k-1)*nelem_tot*np*np + i + j = j + 1 + end do + end do + end do + +end function get_restart_decomp + +!========================================================================================= + +function get_restart_decomp_fvm(elem, lev) result(ldof) + use dimensions_mod, only: nc, nelemd + use element_mod, only: element_t + + type(element_t), intent(in) :: elem(:) + integer, intent(in) :: lev + integer, pointer :: ldof(:) + + integer :: i, j, k, ie + !---------------------------------------------------------------------------- + + allocate(ldof(nelemd*nc*nc*lev)) + + j = 1 + do k = 1, lev + do ie = 1, nelemd + do i = 1, nc*nc + ldof(j) = (elem(ie)%GlobalID-1)*nc*nc + (k-1)*nelem_tot*nc*nc + i + j = j + 1 + end do + end do + end do + +end function get_restart_decomp_fvm + +!========================================================================================= + +end module restart_dynamics diff --git a/src/history/buffers b/src/history/buffers index 170bcb959..a8e4bf4ff 160000 --- a/src/history/buffers +++ b/src/history/buffers @@ -1 +1 @@ -Subproject commit 170bcb9593822b21cdf639c3afcf1477808a4b7a +Subproject commit a8e4bf4ffd69e1332aebab7f5d34ddf2ffb228fd diff --git a/src/history/cam_hist_file.F90 b/src/history/cam_hist_file.F90 index dc4c243ea..8f03456df 100644 --- a/src/history/cam_hist_file.F90 +++ b/src/history/cam_hist_file.F90 @@ -14,6 +14,7 @@ module cam_hist_file use runtime_obj, only: UNSET_I => unset_int use runtime_obj, only: UNSET_C => unset_str use runtime_obj, only: UNSET_R8 => unset_real + use cam_logfile, only: iulog implicit none private @@ -27,6 +28,7 @@ module cam_hist_file integer, parameter :: nl_gname_len = len(hist_nl_group_name) integer, public, parameter :: instantaneous_file_index = 1 integer, public, parameter :: accumulated_file_index = 2 + integer, public, parameter :: restart_file_index = 3 character(len=fieldname_suffix_len ) :: fieldname_suffix = '&IC' ! Suffix appended to field names for IC file logical, parameter :: PATCH_DEF = .true. @@ -41,10 +43,13 @@ module cam_hist_file integer, parameter :: hfile_type_sat_track = 3 integer, parameter :: hfile_type_restart = 4 + character(len=*), parameter :: rh_filename_spec = '%c.cam.r%u.%y-%m-%d-%s.nc' + type :: hist_file_t ! History file configuration information character(len=vlen), private :: volume = UNSET_C type(file_desc_t), private :: hist_files(max_split_files) ! PIO file ids + type(file_desc_t), private :: restart_file ! PIO file id for restart file integer, private :: rl_kind = OUTPUT_DEF integer, private :: max_frames = UNSET_I integer, private :: output_freq_mult = UNSET_I @@ -56,6 +61,8 @@ module cam_hist_file character(len=max_fldlen), allocatable, private :: field_names(:) character(len=3), allocatable, private :: accumulate_types(:) type(var_desc_t), allocatable, private :: file_varids(:,:) + type(var_desc_t), allocatable, private :: nacs_varids(:,:) + type(var_desc_t), allocatable, private :: varbuff_varids(:,:) integer, allocatable, private :: grids(:) integer, private :: hfile_type = hfile_type_default logical, private :: collect_patch_output = PATCH_DEF @@ -67,6 +74,7 @@ module cam_hist_file logical, private :: files_open = .false. type(interp_info_t), pointer, private :: interp_info => NULL() character(len=CL), allocatable, private :: file_names(:) + character(len=CL), allocatable, private :: restart_file_name ! PIO IDs type(var_desc_t), private :: timeid type(var_desc_t), private :: dateid @@ -95,12 +103,25 @@ module cam_hist_file procedure :: get_volume => config_volume procedure :: get_filenames => config_get_filenames procedure :: get_filename_spec => config_get_filename_spec + procedure :: get_restart_filename => config_get_restart_filename procedure :: get_last_month_written => config_get_last_month_written procedure :: get_last_year_written => config_get_last_year_written procedure :: precision => config_precision procedure :: max_frame => config_max_frame procedure :: get_num_samples => config_get_num_samples procedure :: get_beg_time => config_get_beg_time + procedure :: get_field_list => config_get_field_list + procedure :: get_averaging_flags => config_get_averaging_flags + procedure :: get_decompositions => config_get_decompositions + procedure :: get_num_levels => config_get_num_levels + procedure :: get_cell_methods => config_get_cell_methods + procedure :: get_long_names => config_get_long_names + procedure :: get_units => config_get_units + procedure :: get_fill_flags => config_get_fill_flags + procedure :: get_fill_values => config_get_fill_values + procedure :: get_dimension_indices => config_get_dimension_indices + procedure :: get_num_dimensions => config_get_num_dimensions + procedure :: get_num_fields => config_get_num_fields procedure :: output_freq => config_output_freq procedure :: output_freq_separate => config_output_freq_separate procedure :: is_history_file => config_history_file @@ -110,9 +131,18 @@ module cam_hist_file procedure :: do_write_nstep0 => config_do_write_nstep0 procedure :: file_is_setup => config_file_is_setup procedure :: are_files_open => config_files_open + procedure :: has_accumulated_fields => config_has_accumulated_fields ! Actions procedure :: reset => config_reset procedure :: configure => config_configure + procedure :: set_up_dimensions => config_set_up_dimensions + procedure :: define_dimensions => config_define_dimensions + procedure :: define_header_info => config_define_header_info + procedure :: define_instantaneous_header_info => config_define_instantaneous_header_info + procedure :: define_additional_header_info => config_define_additional_header_info + procedure :: define_fields => config_define_fields + procedure :: define_restart_fields => config_define_restart_fields + procedure :: write_invariant_header => config_write_invariant_header procedure :: print_config => config_print_config procedure :: set_beg_time => config_set_beg_time procedure :: set_end_time => config_set_end_time @@ -122,9 +152,13 @@ module cam_hist_file procedure :: set_up_fields => config_set_up_fields procedure :: find_in_field_list => config_find_in_field_list procedure :: define_file => config_define_file + procedure :: define_restart_file => config_define_restart_file procedure :: write_time_dependent_variables => config_write_time_dependent_variables + procedure :: write_restart_file => config_write_restart_file procedure :: write_field => config_write_field + procedure :: write_restart_fields => config_write_restart_fields procedure :: close_files => config_close_files + procedure :: close_restart_file => config_close_restart_file procedure :: clear_buffers => config_clear_buffers procedure :: reset_samples => config_reset_samples end type hist_file_t @@ -187,6 +221,9 @@ subroutine config_set_filenames(this) incomplete_ok=.false.) end do + this%restart_file_name = interpret_filename_spec(rh_filename_spec, & + unit=this%volume, incomplete_ok=.false.) + end subroutine config_set_filenames ! ======================================================================== @@ -224,6 +261,18 @@ end function config_get_filenames ! ======================================================================== + function config_get_restart_filename(this) result(rhfile) + use cam_filenames, only: interpret_filename_spec + ! Dummy arguments + class(hist_file_t), intent(in) :: this + character(len=CL) :: rhfile + + rhfile = this%restart_file_name + + end function config_get_restart_filename + + ! ======================================================================== + pure function config_get_filename_spec(this) result(filename_spec) ! Dummy arguments class(hist_file_t), intent(in) :: this @@ -315,6 +364,238 @@ end function config_get_beg_time ! ======================================================================== + function config_get_field_list(this) result(field_list) + use cam_abortutils, only: endrun + ! Dummy arguments + class(hist_file_t), intent(in) :: this + character(len=max_fldlen), allocatable :: field_list(:) + ! Local variables + character(len=512) :: errmsg + integer :: ierr, idx + + allocate(field_list(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_field_list: failed to allocate field_list; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + field_list(idx) = this%field_list(idx)%diag_name() + end do + + end function config_get_field_list + + ! ======================================================================== + + function config_get_averaging_flags(this) result(avgflags) + use cam_abortutils, only: endrun + use cam_history_support, only: max_chars + ! Dummy arguments + class(hist_file_t), intent(in) :: this + character(len=max_chars), allocatable :: avgflags(:) + ! Local variables + character(len=512) :: errmsg + integer :: ierr, idx + + allocate(avgflags(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_averaging_flags: failed to allocate avgflags; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + avgflags(idx) = this%field_list(idx)%accumulate_type() + end do + end function config_get_averaging_flags + + ! ======================================================================== + + function config_get_decompositions(this) result(decomps) + use cam_abortutils, only: endrun + ! Dummy arguments + class(hist_file_t), intent(in) :: this + integer, allocatable :: decomps(:) + ! Local variables + character(len=512) :: errmsg + integer :: ierr, idx + allocate(decomps(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_decompositions: failed to allocate decomps; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + decomps(idx) = this%field_list(idx)%decomp() + end do + end function config_get_decompositions + + ! ======================================================================== + + function config_get_num_levels(this) result(numlevs) + use cam_abortutils, only: endrun + ! Dummy arguments + class(hist_file_t), intent(in) :: this + integer, allocatable :: numlevs(:) + ! Local variables + character(len=512) :: errmsg + integer :: ierr, idx + allocate(numlevs(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_num_levels: failed to allocate numlevs; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + numlevs(idx) = this%field_list(idx)%num_levels() + end do + end function config_get_num_levels + + ! ======================================================================== + + function config_get_cell_methods(this) result(cell_methods) + use cam_abortutils, only: endrun + use cam_history_support, only: max_chars + ! Dummy arguments + class(hist_file_t), intent(in) :: this + character(len=max_chars), allocatable :: cell_methods(:) + ! Local variables + character(len=512) :: errmsg + integer :: ierr, idx + allocate(cell_methods(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_cell_methods: failed to allocate cell_methods; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + cell_methods(idx) = this%field_list(idx)%cell_methods() + end do + end function config_get_cell_methods + + ! ======================================================================== + + function config_get_long_names(this) result(long_names) + use cam_abortutils, only: endrun + use cam_history_support, only: max_chars + ! Dummy arguments + class(hist_file_t), intent(in) :: this + character(len=max_chars), allocatable :: long_names(:) + ! Local variables + character(len=512) :: errmsg + integer :: idx, ierr + allocate(long_names(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_long_names: failed to allocate long_names; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + long_names(idx) = this%field_list(idx)%long_name() + end do + end function config_get_long_names + + ! ======================================================================== + + function config_get_units(this) result(units) + use cam_abortutils, only: endrun + use cam_history_support, only: max_chars + ! Dummy arguments + class(hist_file_t), intent(in) :: this + character(len=max_chars), allocatable :: units(:) + ! Local variables + character(len=512) :: errmsg + integer :: idx, ierr + allocate(units(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_units: failed to allocate units; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + units(idx) = this%field_list(idx)%units() + end do + end function config_get_units + + ! ======================================================================== + + function config_get_fill_flags(this) result(flags) + use cam_abortutils, only: endrun + ! Dummy arguments + class(hist_file_t), intent(in) :: this + integer, allocatable :: flags(:) + ! Local variables + character(len=512) :: errmsg + integer :: idx, ierr + allocate(flags(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_fill_flags: failed to allocate flags; errmsg = '//errmsg) + end if + flags = 0 + do idx = 1, size(this%field_list) + if (this%field_list(idx)%flag_xyfill()) then + flags(idx) = 1 + end if + end do + end function config_get_fill_flags + + ! ======================================================================== + + function config_get_fill_values(this) result(fill_values) + use cam_abortutils, only: endrun + ! Dummy arguments + class(hist_file_t), intent(in) :: this + real(r8), allocatable :: fill_values(:) + ! Local variables + character(len=512) :: errmsg + integer :: idx, ierr + allocate(fill_values(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_fill_values: failed to allocate fill_values; errmsg = '//errmsg) + end if +! do idx = 1, size(this%field_list) +! fill_values(idx) = this%field_list(idx)%fill_value() +! end do + fill_values = 0._r8 + end function config_get_fill_values + + ! ======================================================================== + + function config_get_dimension_indices(this) result(dimids) + use cam_abortutils, only: endrun + ! Dummy arguments + class(hist_file_t), intent(in) :: this + integer, allocatable :: dimids(:,:) + ! Local variables + character(len=512) :: errmsg + integer :: idx, ierr + allocate(dimids(4, size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_dimension_indices: failed to allocate dimids; errmsg = '//errmsg) + end if + dimids = 0 + do idx = 1, size(this%field_list) + if (size(this%field_list(idx)%dimensions()) > 0) then + dimids(1:size(this%field_list(idx)%dimensions()), idx) = this%field_list(idx)%dimensions() + end if + end do + end function config_get_dimension_indices + + ! ======================================================================== + + function config_get_num_dimensions(this) result (ndims) + use cam_abortutils, only: endrun + class(hist_file_t), intent(in) :: this + integer, allocatable :: ndims(:) + ! Local variables + character(len=512) :: errmsg + integer :: ierr, idx + + allocate(ndims(size(this%field_list)), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_get_num_dimensions: failed to allocate ndims; errmsg = '//errmsg) + end if + do idx = 1, size(this%field_list) + ndims(idx) = size(this%field_list(idx)%dimensions()) + end do + end function config_get_num_dimensions + + ! ======================================================================== + + pure function config_get_num_fields(this) result(nflds) + class(hist_file_t), intent(in) :: this + integer :: nflds + + nflds = size(this%field_list) + + end function config_get_num_fields + ! ======================================================================== + function config_output_freq(this) result(out_freq) use shr_kind_mod, only: CS => SHR_KIND_CS use shr_string_mod, only: to_lower => shr_string_toLower @@ -433,6 +714,16 @@ end function config_files_open ! ======================================================================== + pure function config_has_accumulated_fields(this) result(has_accum) + ! Dummy arguments + class(hist_file_t), intent(in) :: this + logical :: has_accum + + has_accum = this%has_accumulated + end function config_has_accumulated_fields + + ! ======================================================================== + subroutine config_reset(this) ! Dummy argument class(hist_file_t), intent(inout) :: this @@ -901,24 +1192,16 @@ end subroutine AvgflagToString ! ======================================================================== - subroutine config_define_file(this, restart, logname, host, model_doi_url) - use pio, only: PIO_CLOBBER, pio_file_is_open, pio_unlimited - use pio, only: pio_double, pio_def_var, pio_put_att, pio_int - use pio, only: PIO_GLOBAL, pio_char, pio_real, PIO_NOERR, pio_enddef - use pio, only: pio_put_var - use cam_pio_utils, only: cam_pio_createfile, cam_pio_def_var - use cam_pio_utils, only: cam_pio_def_dim, cam_pio_handle_error + subroutine config_define_file(this, logname, host, model_doi_url) + use pio, only: PIO_CLOBBER, pio_file_is_open + use pio, only: PIO_NOERR, pio_enddef + use cam_pio_utils, only: cam_pio_createfile use cam_grid_support, only: cam_grid_header_info_t, cam_grid_write_attr - use cam_grid_support, only: cam_grid_write_var - use cam_history_support, only: write_hist_coord_attrs use cam_history_support, only: write_hist_coord_vars use cam_history_support, only: max_chars use shr_kind_mod, only: r4 => shr_kind_r4 use time_manager, only: get_ref_date, timemgr_get_calendar_cf use time_manager, only: get_step_size - use string_utils, only: date2yyyymmdd, sec2hms, stringify - use cam_control_mod, only: caseid - use cam_initfiles, only: ncdata, bnd_topo use cam_abortutils, only: check_allocate, endrun use cam_logfile, only: iulog use spmd_utils, only: masterproc @@ -926,55 +1209,54 @@ subroutine config_define_file(this, restart, logname, host, model_doi_url) ! Define the metadata for the file(s) for this volume ! Dummy arguments class(hist_file_t), intent(inout) :: this - logical, intent(in) :: restart character(len=*), intent(in) :: logname character(len=*), intent(in) :: host character(len=*), intent(in) :: model_doi_url ! Local variables integer :: amode, ierr - integer :: grid_index, split_file_index, field_index, idx, jdx - integer :: dtime ! timestep size + integer :: split_file_index, idx integer :: yr, mon, day ! year, month, day components of a date integer :: nbsec ! time of day component of base date [seconds] integer :: nbdate ! base date in yyyymmdd format integer :: nsbase = 0 ! seconds component of base time integer :: ndbase = 0 ! days component of base time - integer :: ncreal ! real data type for output - integer :: grd - integer :: mdimsize, num_hdims, fdims + integer :: num_hdims + integer, allocatable :: fdims(:) integer :: num_patches integer, allocatable :: mdims(:) logical :: is_satfile logical :: is_initfile logical :: varid_set - character(len=16) :: time_per_freq - character(len=max_chars) :: str ! character temporary + logical :: restart, split_file character(len=max_chars) :: calendar ! Calendar type - character(len=max_chars) :: cell_methods ! For cell_methods attribute - character(len=max_chars) :: fname_tmp ! local copy of field name - character(len=CM) :: errmsg - type(var_desc_t) :: varid ! NetCDF dimensions integer :: timdim ! unlimited dimension id integer :: bnddim ! bounds dimension id integer :: chardim ! character dimension id - integer :: dimenchar(2) ! character dimension ids - integer :: nacsdims(2) ! dimension ids for nacs (used in restart file) - integer :: max_mdims ! maximum number of middle dimensions - integer :: max_hdims ! maximum number of grid dimensions integer, allocatable :: dimindex(:) ! dimension ids for variable declaration ! A structure to hold the horizontal dimension and coordinate info type(cam_grid_header_info_t), allocatable :: header_info(:) - integer, allocatable :: mdimids(:) - integer, parameter :: max_netcdf_len = 256 + integer, allocatable :: mdimids(:) character(len=*), parameter :: subname = 'config_define_file: ' is_initfile = (this%hfile_type == hfile_type_init_value) is_satfile = (this%hfile_type == hfile_type_sat_track) + + restart = .false. + split_file = .true. + + ! Get date and calendar info + call get_ref_date(yr, mon, day, nbsec) + nbdate = yr*10000 + mon*100 + day + calendar = timemgr_get_calendar_cf() + ! v peverwhee - remove when patch output is set up + num_patches = 1 + ! ^ peverwhee - remove when patch output is set up + ! Log what we're doing if (masterproc) then if (this%has_accumulated) then @@ -1001,339 +1283,32 @@ subroutine config_define_file(this, restart, logname, host, model_doi_url) this%files_open = .true. - allocate(header_info(size(this%grids)), stat=ierr) - call check_allocate(ierr, subname, 'header_info', & - file=__FILE__, line=__LINE__-1) - - max_hdims = 0 - do grid_index = 1, size(this%grids) - do split_file_index = 1, max_split_files - if (pio_file_is_open(this%hist_files(split_file_index))) then - call cam_grid_write_attr(this%hist_files(split_file_index), & - this%grids(grid_index), header_info(grid_index), & - file_index=split_file_index) - max_hdims = max(max_hdims, header_info(grid_index)%num_hdims()) - end if - end do - end do - ! Determine the maximum number of dimensions - max_mdims = 0 - do field_index = 1, size(this%field_list) - max_mdims = max(max_mdims, size(this%field_list(field_index)%dimensions())) - end do - ! Allocate dimindex to the maximum possible dimensions (plus 1 for time) - allocate(dimindex(max_hdims + max_mdims + 1), stat=ierr) - call check_allocate(ierr, subname, 'dimindex', file=__FILE__, line=__LINE__-1) - - call get_ref_date(yr, mon, day, nbsec) - nbdate = yr*10000 + mon*100 + day - calendar = timemgr_get_calendar_cf() - dtime = get_step_size() - ! v peverwhee - remove when patch output is set up - num_patches = 1 - ! ^ peverwhee - remove when patch output is set up - varid_set = .true. - ! Allocate the varid array - if (.not. allocated(this%file_varids)) then - allocate(this%file_varids(size(this%field_list), num_patches), stat=ierr) - call check_allocate(ierr, subname, 'this%file_varids', & - file=__FILE__, line=__LINE__-1) - varid_set = .false. - end if - - ! Format frequency - write(time_per_freq,999) trim(this%output_freq_type), '_', this%output_freq_mult -999 format(2a,i0) + call this%set_up_dimensions(num_patches, restart, header_info, dimindex, varid_set) + allocate(fdims(size(this%field_list))) do split_file_index = 1, max_split_files if (.not. pio_file_is_open(this%hist_files(split_file_index))) then cycle end if - ! Define the unlimited time dim - call cam_pio_def_dim(this%hist_files(split_file_index), 'time', pio_unlimited, timdim) - call cam_pio_def_dim(this%hist_files(split_file_index), 'nbnd', 2, bnddim, existOK=.true.) - call cam_pio_def_dim(this%hist_files(split_file_index), 'chars', 8, chardim) - ! Populate the history coordinate (well, mdims anyway) attributes - ! This routine also allocates the mdimids array - call write_hist_coord_attrs(this%hist_files(split_file_index), bnddim, mdimids, restart) - ! Define time variable - ierr=pio_def_var (this%hist_files(split_file_index),'time',pio_double,(/timdim/),this%timeid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "time" variable') - ierr=pio_put_att (this%hist_files(split_file_index), this%timeid, 'long_name', 'time') - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "time" variable') - str = 'days since ' // date2yyyymmdd(nbdate) // ' ' // sec2hms(nbsec) - ierr=pio_put_att (this%hist_files(split_file_index), this%timeid, 'units', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "units" attribtue to "time" variable') - ierr=pio_put_att (this%hist_files(split_file_index), this%timeid, 'calendar', trim(calendar)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "calendar" attribute to "time" variable') - - ! Define date variable - ierr=pio_def_var (this%hist_files(split_file_index),'date ',pio_int,(/timdim/),this%dateid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "date" variable') - str = 'current date (YYYYMMDD)' - ierr=pio_put_att (this%hist_files(split_file_index), this%dateid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "date" variable') - - ! Define datesec variable - ierr=pio_def_var (this%hist_files(split_file_index),'datesec ',pio_int,(/timdim/), this%datesecid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "datesec" variable') - str = 'current seconds of current date' - ierr=pio_put_att (this%hist_files(split_file_index), this%datesecid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "datesec" variable') - - ! - ! Character header information - ! - str = 'CF-1.0' - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'Conventions', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "Conventions" attribtue to file') - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'source', 'CAM-SIMA') - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "source" attribute to file') - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'case', caseid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "case" attribute to file') - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'logname',trim(logname)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "logname" attribute to file') - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'host', trim(host)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "host" attribute to file') - - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'initial_file', ncdata) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "initial_file" attribute to file') - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'topography_file', bnd_topo) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "topography_file" attribute to file') - if (len_trim(model_doi_url) > 0) then - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'model_doi_url', model_doi_url) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "model_doi_url" attribute to file') - end if - - ierr=pio_put_att (this%hist_files(split_file_index), PIO_GLOBAL, 'time_period_freq', trim(time_per_freq)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "time_period_freq" attribute to file') + ! Define the file dimensions + call this%define_dimensions(this%hist_files(split_file_index), restart, timdim, bnddim, chardim, mdimids) + ! Define the header information + call this%define_header_info(this%hist_files(split_file_index), timdim, bnddim, chardim, nbsec, nbdate, & + calendar, logname, host, model_doi_url) if (.not. is_satfile) then - ! Define time_bounds variable - ierr=pio_put_att (this%hist_files(split_file_index), this%timeid, 'bounds', 'time_bounds') - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "bounds" attribute to file') - ierr=pio_def_var (this%hist_files(split_file_index),'time_bounds',pio_double,(/bnddim,timdim/),this%tbndid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "time_bounds" variable') - ierr=pio_put_att (this%hist_files(split_file_index), this%tbndid, 'long_name', 'time interval endpoints') - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "time_bounds" variable') - str = 'days since ' // date2yyyymmdd(nbdate) // ' ' // sec2hms(nbsec) - ierr=pio_put_att (this%hist_files(split_file_index), this%tbndid, 'units', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "units" attribute to "time_bounds" variable') - ierr=pio_put_att (this%hist_files(split_file_index), this%tbndid, 'calendar', trim(calendar)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "calendar" attribute to "time_bounds" variable') - - ! - ! Character - ! - dimenchar(1) = chardim - dimenchar(2) = timdim - ierr=pio_def_var (this%hist_files(split_file_index),'date_written',pio_char,dimenchar,this%date_writtenid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "date_written" variable') - ierr=pio_def_var (this%hist_files(split_file_index),'time_written',pio_char,dimenchar,this%time_writtenid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "time_written" variable') - - ! - ! Integer header - ! - ! Define base day variables - ierr=pio_def_var (this%hist_files(split_file_index),'ndbase',PIO_INT,this%ndbaseid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "ndbase" variable') - str = 'base day' - ierr=pio_put_att (this%hist_files(split_file_index), this%ndbaseid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "ndbase" variable') - - ierr=pio_def_var (this%hist_files(split_file_index),'nsbase',PIO_INT,this%nsbaseid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "nsbase" variable') - str = 'seconds of base day' - ierr=pio_put_att (this%hist_files(split_file_index), this%nsbaseid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "nsbase" variable') - - ierr=pio_def_var (this%hist_files(split_file_index),'nbdate',PIO_INT,this%nbdateid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "nbdate" variable') - str = 'base date (YYYYMMDD)' - ierr=pio_put_att (this%hist_files(split_file_index), this%nbdateid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribtue to "nbdate" variable') - - ierr=pio_def_var (this%hist_files(split_file_index),'nbsec',PIO_INT,this%nbsecid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "nbsec" variable') - str = 'seconds of base date' - ierr=pio_put_att (this%hist_files(split_file_index), this%nbsecid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "nbsec" variable') - - ierr=pio_def_var (this%hist_files(split_file_index),'mdt',PIO_INT,this%mdtid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "mdt" variable') - ierr=pio_put_att (this%hist_files(split_file_index), this%mdtid, 'long_name', 'timestep') - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "mdt" variable') - ierr=pio_put_att (this%hist_files(split_file_index), this%mdtid, 'units', 's') - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "units" attribute to "mdt" variable') - - ! - ! Create variables for model timing and header information - ! + call this%define_additional_header_info(this%hist_files(split_file_index), bnddim, timdim, chardim, & + nbsec, nbdate, calendar) if (split_file_index == instantaneous_file_index) then - ierr=pio_def_var (this%hist_files(split_file_index),'ndcur ',pio_int,(/timdim/),this%ndcurid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "ndcur" variable') - str = 'current day (from base day)' - ierr=pio_put_att (this%hist_files(split_file_index), this%ndcurid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "ndcur" variable') - - ierr=pio_def_var (this%hist_files(split_file_index),'nscur ',pio_int,(/timdim/),this%nscurid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "nscur" variable') - str = 'current seconds of current day' - ierr=pio_put_att (this%hist_files(split_file_index), this%nscurid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "nscur" variable') - ierr=pio_def_var (this%hist_files(split_file_index),'nsteph',pio_int,(/timdim/),this%nstephid) - call cam_pio_handle_error(ierr, 'config_define_file: failed to define "nsteph" variable') - str = 'current timestep' - ierr=pio_put_att (this%hist_files(split_file_index), this%nstephid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: failed to add "long_name" attribute to "nsteph" variable') + call this%define_instantaneous_header_info(this%hist_files(split_file_index), timdim) end if end if ! .not. satfile - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - ! ! Create variables and attributes for field list - ! - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - - do field_index = 1, size(this%field_list) - if (.not. is_satfile .and. .not. restart .and. .not. is_initfile) then - if (split_file_index == accumulated_file_index) then - ! this is the accumulated file of a potentially split - ! history tape - skip instantaneous fields - if (this%field_list(field_index)%accumulate_type() == 'lst') then - cycle - end if - else - ! this is the instantaneous file of a potentially split - ! history tape - skip accumulated fields - if (this%field_list(field_index)%accumulate_type() /= 'lst') then - cycle - end if - end if - end if + call this%define_fields(this%hist_files(split_file_index), split_file_index, split_file, restart, mdimids, & + dimindex, timdim, num_patches, header_info, varid_set, num_hdims, fdims) - if (this%precision() == 'REAL32') then - ncreal = pio_real - else if (this%precision() == 'REAL64') then - ncreal = pio_double - end if - mdims = this%field_list(field_index)%dimensions() - mdimsize = size(mdims) - fname_tmp = strip_suffix(this%field_list(field_index)%diag_name()) - ! Ensure that fname_tmp is not longer than the maximum length for a - ! netcdf file - if (len_trim(fname_tmp) > max_netcdf_len) then - ! Endrun if the name is too long - write(errmsg, *) 'config_define_file: variable name ', trim(fname_tmp), & - ' too long for NetCDF file (len=', stringify((/len(trim(fname_tmp))/)), ' > ', & - stringify((/max_netcdf_len/)), ')' - call endrun(errmsg, file=__FILE__, line=__LINE__) - end if - ! - ! Create variables and atributes for fields written out as columns - ! - ! Find appropriate grid in header_info - if (.not. allocated(header_info)) then - ! Safety check - call endrun('config_define_file: header_info not allocated', file=__FILE__, line=__LINE__) - end if - grd = -1 - do idx = 1, size(header_info) - if (header_info(idx)%get_gridid() == this%field_list(field_index)%decomp()) then - grd = idx - exit - end if - end do - if (grd < 0) then - write(errmsg, '(a,i0,2a)') 'grid, ',this%field_list(field_index)%decomp(),', not found for ', & - trim(fname_tmp) - call endrun('config_define_file: '//errmsg, file=__FILE__, line=__LINE__) - end if - num_hdims = header_info(grd)%num_hdims() - do idx = 1, num_hdims - dimindex(idx) = header_info(grd)%get_hdimid(idx) - nacsdims(idx) = header_info(grd)%get_hdimid(idx) - end do - do idx = 1, num_patches - varid = this%file_varids(field_index, idx) - ! Figure the dimension ID array for this field - ! We have defined the horizontal grid dimensions in dimindex - fdims = num_hdims - do jdx = 1, mdimsize - fdims = fdims + 1 - dimindex(fdims) = mdimids(mdims(jdx)) - end do - if(.not. restart) then - ! Only add time dimension if this is not a restart history tape - fdims = fdims + 1 - dimindex(fdims) = timdim - end if - ! peverwhee - TODO: enable patch output - ! Define the variable - call cam_pio_def_var(this%hist_files(split_file_index), trim(fname_tmp), ncreal, & - dimindex(1:fdims), varid) - if (.not. varid_set) then - this%file_varids(field_index, idx) = varid - end if - if (mdimsize > 0) then - ierr = pio_put_att(this%hist_files(split_file_index), varid, 'mdims', mdims(1:mdimsize)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot define mdims for '//trim(fname_tmp)) - end if - str = this%field_list(field_index)%sampling_sequence() - if (len_trim(str) > 0) then - ierr = pio_put_att(this%hist_files(split_file_index), varid, 'Sampling_Sequence', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot define Sampling_Sequence for '//trim(fname_tmp)) - end if - if (this%field_list(field_index)%flag_xyfill()) then - ! Write _FillValue and missing_value attributes so that - ! netCDF-aware tools (ncview, ncdump, etc.) recognise fill - ! values. The attribute type must match the variable type. - if (ncreal == pio_double) then - ierr = pio_put_att(this%hist_files(split_file_index), varid, '_FillValue', this%field_list(field_index)%fill_value()) - call cam_pio_handle_error(ierr, subname//'cannot define _FillValue for '//trim(fname_tmp)) - ierr = pio_put_att(this%hist_files(split_file_index), varid, 'missing_value', this%field_list(field_index)%fill_value()) - call cam_pio_handle_error(ierr, subname//'cannot define missing_value for '//trim(fname_tmp)) - else - ierr = pio_put_att(this%hist_files(split_file_index), varid, '_FillValue', real(this%field_list(field_index)%fill_value(), r4)) - call cam_pio_handle_error(ierr, subname//'cannot define _FillValue for '//trim(fname_tmp)) - ierr = pio_put_att(this%hist_files(split_file_index), varid, 'missing_value', real(this%field_list(field_index)%fill_value(), r4)) - call cam_pio_handle_error(ierr, subname//'cannot define missing_value for '//trim(fname_tmp)) - end if - end if - str = this%field_list(field_index)%units() - if (len_trim(str) > 0) then - ierr=pio_put_att (this%hist_files(split_file_index), varid, 'units', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot define units for '//trim(fname_tmp)) - end if - str = this%field_list(field_index)%mixing_ratio() - if (len_trim(str) > 0) then - ierr=pio_put_att (this%hist_files(split_file_index), varid, 'mixing_ratio', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot define mixing_ratio for '//trim(fname_tmp)) - end if - str = this%field_list(field_index)%long_name() - ierr=pio_put_att (this%hist_files(split_file_index), varid, 'long_name', trim(str)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot define long_name for '//trim(fname_tmp)) - ! Assign field attributes defining valid levels and averaging info - cell_methods = '' - if (len_trim(this%field_list(field_index)%cell_methods()) > 0) then - if (len_trim(cell_methods) > 0) then - cell_methods = trim(cell_methods)//' '//trim(this%field_list(field_index)%cell_methods()) - else - cell_methods = trim(this%field_list(field_index)%cell_methods()) - end if - end if - ! Time cell methods is after field method because time averaging is - ! applied later (just before output) than field method which is applied - ! before outfld call. - call AvgflagToString(this%field_list(field_index)%accumulate_type(), str) - cell_methods = adjustl(trim(cell_methods)//' '//'time: '//str) - ierr = pio_put_att(this%hist_files(split_file_index), varid, 'cell_methods', trim(cell_methods)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot define cell_methods for '//trim(fname_tmp)) - end do ! end loop over patches - deallocate(mdims) - end do ! end loop over fields deallocate(mdimids) ierr = pio_enddef(this%hist_files(split_file_index)) if (ierr /= PIO_NOERR) then @@ -1343,34 +1318,12 @@ subroutine config_define_file(this, restart, logname, host, model_doi_url) write(iulog,*)'config_define_file: Successfully opened netcdf file '//trim(this%file_names(split_file_index)) end if - ! ! Write time-invariant portion of history header - if (.not. is_satfile) then - do idx = 1, size(this%grids) - call cam_grid_write_var(this%hist_files(split_file_index), this%grids(idx), & - file_index=split_file_index) - end do - ierr = pio_put_var(this%hist_files(split_file_index), this%mdtid, (/dtime/)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot put mdt') - - ! - ! Model date info - ! - ierr = pio_put_var(this%hist_files(split_file_index), this%ndbaseid, (/ndbase/)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot put ndbase') - ierr = pio_put_var(this%hist_files(split_file_index), this%nsbaseid, (/nsbase/)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot put nsbase') - - ierr = pio_put_var(this%hist_files(split_file_index), this%nbdateid, (/nbdate/)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot put nbdate') - ierr = pio_put_var(this%hist_files(split_file_index), this%nbsecid, (/nbsec/)) - call cam_pio_handle_error(ierr, 'config_define_file: cannot put nbsec') + call this%write_invariant_header(this%hist_files(split_file_index), split_file_index, & + ndbase, nsbase, nbdate, nbsec, restart) end if - ! Write the mdim variable data - call write_hist_coord_vars(this%hist_files(split_file_index), restart) - end do ! end loop over files if (allocated(header_info)) then @@ -1380,51 +1333,840 @@ subroutine config_define_file(this, restart, logname, host, model_doi_url) deallocate(header_info) end if + if (allocated(dimindex)) then + deallocate(dimindex) + end if + end subroutine config_define_file ! ======================================================================== - subroutine config_write_time_dependent_variables(this, restart) - use pio, only: pio_put_var, pio_file_is_open - use time_manager, only: get_nstep, get_curr_date, get_curr_time - use time_manager, only: set_date_from_time_float, get_step_size - use datetime_mod, only: datetime - use spmd_utils, only: masterproc - use cam_logfile, only: iulog - use perf_mod, only: t_startf, t_stopf - use cam_pio_utils, only: cam_pio_handle_error + subroutine config_define_restart_file(this, logname, host, model_doi_url) + use time_manager, only: timemgr_get_calendar_cf, get_ref_date + use spmd_utils, only: masterproc + use cam_logfile, only: iulog + use cam_history_support, only: max_chars, write_hist_coord_vars + use cam_pio_utils, only: cam_pio_createfile + use pio, only: pio_clobber, pio_noerr, pio_enddef + use cam_grid_support, only: cam_grid_header_info_t + use cam_abortutils, only: endrun + ! Define the metadata for the history restart file (rhX) ! Dummy arguments class(hist_file_t), intent(inout) :: this - logical, intent(in) :: restart - + character(len=*), intent(in) :: logname + character(len=*), intent(in) :: host + character(len=*), intent(in) :: model_doi_url ! Local variables - integer :: yr, mon, day ! year, month, and day components of a date - integer :: yr_mid, mon_mid, day_mid ! year, month, and day components of midpoint date - integer :: nstep ! current timestep number - integer :: ncdate(max_split_files) ! current (or midpoint) date in integer format [yyyymmdd] - integer :: ncsec(max_split_files) ! current (or midpoint) time of day [seconds] - integer :: ndcur ! day component of current time - integer :: nscur ! seconds component of current time - real(r8) :: time ! current (or midpoint) time - real(r8) :: time_interval(2) ! time interval boundaries - integer :: ierr - integer :: split_file_index, field_idx - integer :: start, count1 - integer :: startc(2) ! start values required by nf_put_vara (character) - integer :: countc(2) ! count values required by nf_put_vara (character) - logical :: is_initfile, is_satfile - character(len=8) :: cdate ! system date - character(len=8) :: ctime ! system time + integer :: amode, ierr, idx + integer :: yr, mon, day + integer :: nbsec, nbdate + integer :: ndbase = 0 + integer :: nsbase = 0 + integer :: num_patches + integer :: num_hdims + integer, allocatable :: fdims(:) + integer :: timdim, bnddim, chardim + logical :: varid_set + logical :: restart, split_file + character(len=max_chars) :: calendar ! Calendar type + integer, allocatable :: dimindex(:) ! dimension ids for variable declaration + ! A structure to hold the horizontal dimension and coordinate info + type(cam_grid_header_info_t), allocatable :: header_info(:) + integer, allocatable :: mdimids(:) - nstep = get_nstep() - call get_curr_date(yr, mon, day, ncsec(instantaneous_file_index)) - ncdate(instantaneous_file_index) = yr*10000 + mon*100 + day - call get_curr_time(ndcur, nscur) - time = ndcur + nscur/86400._r8 - time_interval(1) = this%beg_time - time_interval(2) = time - call set_date_from_time_float((time_interval(1) + time_interval(2)) / 2._r8, & - yr_mid, mon_mid, day_mid, ncsec(accumulated_file_index)) + restart = .true. + split_file = .false. + ! Get date and calendar info + call get_ref_date(yr, mon, day, nbsec) + nbdate = yr*10000 + mon*100 + day + calendar = timemgr_get_calendar_cf() + + ! v peverwhee - remove when patch output is set up + num_patches = 1 + ! ^ peverwhee - remove when patch output is set up + + ! Log what we're doing + if (masterproc) then + write(iulog,*) 'Opening history restart file to write, ', & + trim(this%restart_file_name) + end if + + amode = PIO_CLOBBER + + ! Create the restart history file + call cam_pio_createfile(this%restart_file, this%restart_file_name, amode) + + ! Deallocate file varids from history files (we don't need those anymore) + if (allocated(this%file_varids)) then + deallocate(this%file_varids) + end if + + ! Set up data structures for dimension handling + call this%set_up_dimensions(num_patches, restart, header_info, dimindex, varid_set) + + ! Begin header definition + !-- Define the file dimensions + call this%define_dimensions(this%restart_file, restart, timdim, bnddim, chardim, mdimids) + !-- Define the header information + call this%define_header_info(this%restart_file, timdim, bnddim, chardim, nbsec, nbdate, & + calendar, logname, host, model_doi_url) + call this%define_instantaneous_header_info(this%restart_file, timdim) + call this%define_additional_header_info(this%restart_file, bnddim, timdim, chardim, & + nbsec, nbdate, calendar) + ! End header definition + + ! Create variables and attributes for field list + allocate(fdims(size(this%field_list))) + call this%define_fields(this%restart_file, 1, split_file, restart, mdimids, & + dimindex, timdim, num_patches, header_info, varid_set, num_hdims, fdims) + + ! Define additional restart fields + call this%define_restart_fields(num_patches, dimindex, num_hdims, fdims) + + deallocate(mdimids) + ierr = pio_enddef(this%restart_file) + if (ierr /= PIO_NOERR) then + call endrun('config_define_restart_file: ERROR exiting define mode in PIO', file=__FILE__, line=__LINE__) + end if + if(masterproc) then + write(iulog,*)'config_define_restart_file: Successfully opened netcdf file '//trim(this%restart_file_name) + end if + + ! Write time-invariant header info + call this%write_invariant_header(this%restart_file, restart_file_index, ndbase, nsbase, nbdate, nbsec, restart) + + ! Clean up + if (allocated(header_info)) then + do idx = 1, size(header_info) + call header_info(idx)%deallocate() + end do + deallocate(header_info) + end if + + end subroutine config_define_restart_file + + ! ======================================================================== + + subroutine config_set_up_dimensions(this, num_patches, restart, header_info, dimindex, varid_set) + use cam_grid_support, only: cam_grid_header_info_t, cam_grid_write_attr + use cam_abortutils, only: check_allocate + use pio, only: pio_file_is_open + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + integer, intent(in) :: num_patches + logical, intent(in) :: restart + type(cam_grid_header_info_t), allocatable, intent(out) :: header_info(:) + integer, allocatable, intent(out) :: dimindex(:) + logical, intent(out) :: varid_set + ! Local variables + integer :: ierr + integer :: grid_index, field_index, split_file_index + integer :: max_hdims, max_mdims + character(len=512) :: errmsg + character(len=*), parameter :: subname = 'config_set_up_dimensions: ' + + allocate(header_info(size(this%grids)), stat=ierr, errmsg=errmsg) + call check_allocate(ierr, subname, 'header_info', & + file=__FILE__, line=__LINE__-1, errmsg=errmsg) + + max_hdims = 0 + do grid_index = 1, size(this%grids) + if (restart) then + if (pio_file_is_open(this%restart_file)) then + call cam_grid_write_attr(this%restart_file, this%grids(grid_index), & + header_info(grid_index), file_index=restart_file_index) + end if + else + do split_file_index = 1, max_split_files + if (pio_file_is_open(this%hist_files(split_file_index))) then + call cam_grid_write_attr(this%hist_files(split_file_index), & + this%grids(grid_index), header_info(grid_index), & + file_index=split_file_index) + max_hdims = max(max_hdims, header_info(grid_index)%num_hdims()) + end if + end do + end if + end do + ! Determine the maximum number of dimensions + max_mdims = 0 + do field_index = 1, size(this%field_list) + max_mdims = max(max_mdims, size(this%field_list(field_index)%dimensions())) + end do + ! Allocate dimindex to the maximum possible dimensions (plus 1 for time) + allocate(dimindex(max_hdims + max_mdims + 1), stat=ierr, errmsg=errmsg) + call check_allocate(ierr, subname, 'dimindex', file=__FILE__, line=__LINE__-1, & + errmsg=errmsg) + + varid_set = .true. + ! Allocate the varid array + if (.not. allocated(this%file_varids)) then + allocate(this%file_varids(size(this%field_list), num_patches), stat=ierr, errmsg=errmsg) + call check_allocate(ierr, subname, 'this%file_varids', & + file=__FILE__, line=__LINE__-1, errmsg=errmsg) + varid_set = .false. + end if + + end subroutine config_set_up_dimensions + + ! ======================================================================== + + subroutine config_define_dimensions(this, hist_file, restart, timdim, & + bnddim, chardim, mdimids) + use cam_pio_utils, only: cam_pio_def_dim + use pio, only: pio_unlimited + use cam_history_support, only: write_hist_coord_attrs + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + logical, intent(in) :: restart + integer, intent(out) :: timdim ! unlimited dimension id + integer, intent(out) :: bnddim ! bounds dimension id + integer, intent(out) :: chardim ! character dimension id + integer, allocatable, intent(out) :: mdimids(:) ! mdim ids + + ! Define the unlimited time dim + call cam_pio_def_dim(hist_file, 'time', pio_unlimited, timdim) + call cam_pio_def_dim(hist_file, 'nbnd', 2, bnddim, existOK=.true.) + call cam_pio_def_dim(hist_file, 'chars', 8, chardim) + ! Populate the history coordinate (well, mdims anyway) attributes + ! This routine also allocates the mdimids array + call write_hist_coord_attrs(hist_file, bnddim, mdimids, restart) + + end subroutine config_define_dimensions + + ! ======================================================================== + + subroutine config_define_header_info(this, hist_file, timdim, bnddim, chardim, & + nbsec, nbdate, calendar, logname, host, model_doi_url) + use pio, only: pio_def_var, pio_double, pio_int, pio_put_att, PIO_GLOBAL + use cam_pio_utils, only: cam_pio_handle_error + use string_utils, only: date2yyyymmdd, sec2hms + use cam_history_support, only: max_chars + use cam_control_mod, only: caseid + use cam_initfiles, only: ncdata, bnd_topo + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + integer, intent(in) :: timdim ! unlimited dimension id + integer, intent(in) :: bnddim ! bounds dimension id + integer, intent(in) :: chardim ! character dimension id + integer, intent(in) :: nbsec ! time of day component of base date [seconds] + integer, intent(in) :: nbdate ! base date in yyyymmdd format + character(len=*), intent(in) :: calendar ! Calendar type + character(len=*), intent(in) :: logname + character(len=*), intent(in) :: host + character(len=*), intent(in) :: model_doi_url + ! Local variables + integer :: ierr + character(len=max_chars) :: str ! character temporary + character(len=16) :: time_per_freq + + ! Format frequency + write(time_per_freq,999) trim(this%output_freq_type), '_', this%output_freq_mult +999 format(2a,i0) + + ! Define time variable + ierr=pio_def_var (hist_file,'time',pio_double,(/timdim/),this%timeid) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to define "time" variable') + ierr=pio_put_att (hist_file, this%timeid, 'long_name', 'time') + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "long_name" attribute to "time" variable') + str = 'days since ' // date2yyyymmdd(nbdate) // ' ' // sec2hms(nbsec) + ierr=pio_put_att (hist_file, this%timeid, 'units', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "units" attribtue to "time" variable') + ierr=pio_put_att (hist_file, this%timeid, 'calendar', trim(calendar)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "calendar" attribute to "time" variable') + + ! Define date variable + ierr=pio_def_var (hist_file,'date ',pio_int,(/timdim/),this%dateid) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to define "date" variable') + str = 'current date (YYYYMMDD)' + ierr=pio_put_att (hist_file, this%dateid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "long_name" attribute to "date" variable') + + ! Define datesec variable + ierr=pio_def_var (hist_file,'datesec ',pio_int,(/timdim/), this%datesecid) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to define "datesec" variable') + str = 'current seconds of current date' + ierr=pio_put_att (hist_file, this%datesecid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "long_name" attribute to "datesec" variable') + ! + ! Character header information + ! + str = 'CF-1.0' + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'Conventions', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "Conventions" attribtue to file') + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'source', 'CAM-SIMA') + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "source" attribute to file') + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'case', caseid) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "case" attribute to file') + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'logname',trim(logname)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "logname" attribute to file') + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'host', trim(host)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "host" attribute to file') + + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'initial_file', ncdata) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "initial_file" attribute to file') + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'topography_file', bnd_topo) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "topography_file" attribute to file') + if (len_trim(model_doi_url) > 0) then + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'model_doi_url', model_doi_url) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "model_doi_url" attribute to file') + end if + ierr=pio_put_att (hist_file, PIO_GLOBAL, 'time_period_freq', trim(time_per_freq)) + call cam_pio_handle_error(ierr, 'config_define_header_info: failed to add "time_period_freq" attribute to file') + end subroutine config_define_header_info + + ! ======================================================================== + + subroutine config_define_instantaneous_header_info(this, hist_file, timdim) + use pio, only: pio_def_var, pio_int, pio_put_att + use cam_pio_utils, only: cam_pio_handle_error + use cam_history_support, only: max_chars + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + integer, intent(in) :: timdim ! unlimited dimension id + ! Local variables + integer :: ierr + character(len=max_chars) :: str ! Temporary string + + ierr=pio_def_var (hist_file,'ndcur ',pio_int,(/timdim/),this%ndcurid) + call cam_pio_handle_error(ierr, 'config_define_instantaneous_header_info: failed to define "ndcur" variable') + str = 'current day (from base day)' + ierr=pio_put_att (hist_file, this%ndcurid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_instantaneous_header_info: failed to add "long_name" attribute to "ndcur" variable') + + ierr=pio_def_var (hist_file,'nscur ',pio_int,(/timdim/),this%nscurid) + call cam_pio_handle_error(ierr, 'config_define_instantaneous_header_info: failed to define "nscur" variable') + str = 'current seconds of current day' + ierr=pio_put_att (hist_file, this%nscurid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_instantaneous_header_info: failed to add "long_name" attribute to "nscur" variable') + ierr=pio_def_var (hist_file,'nsteph',pio_int,(/timdim/),this%nstephid) + call cam_pio_handle_error(ierr, 'config_define_instantaneous_header_info: failed to define "nsteph" variable') + str = 'current timestep' + ierr=pio_put_att (hist_file, this%nstephid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_instantaneous_header_info: failed to add "long_name" attribute to "nsteph" variable') + + end subroutine config_define_instantaneous_header_info + + ! ======================================================================== + + subroutine config_define_additional_header_info(this, hist_file, bnddim, timdim, chardim, & + nbsec, nbdate, calendar) + use pio, only: pio_def_var, pio_double, pio_put_att, pio_char, pio_int + use cam_pio_utils, only: cam_pio_handle_error + use string_utils, only: date2yyyymmdd, sec2hms + use cam_history_support, only: max_chars + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + integer, intent(in) :: timdim ! unlimited dimension id + integer, intent(in) :: bnddim ! bounds dimension id + integer, intent(in) :: chardim ! character dimension id + integer, intent(in) :: nbsec ! time of day component of base date [seconds] + integer, intent(in) :: nbdate ! base date in yyyymmdd format + character(len=*), intent(in) :: calendar ! Calendar type + ! Local variables + integer :: dimenchar(2), ierr + character(len=max_chars) :: str ! character temporary + + dimenchar(1) = chardim + dimenchar(2) = timdim + ! Define time_bounds variable + ierr=pio_put_att (hist_file, this%timeid, 'bounds', 'time_bounds') + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "bounds" attribute to file') + ierr=pio_def_var (hist_file,'time_bounds',pio_double,(/bnddim,timdim/),this%tbndid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "time_bounds" variable') + ierr=pio_put_att (hist_file, this%tbndid, 'long_name', 'time interval endpoints') + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "long_name" attribute to "time_bounds" variable') + str = 'days since ' // date2yyyymmdd(nbdate) // ' ' // sec2hms(nbsec) + ierr=pio_put_att (hist_file, this%tbndid, 'units', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "units" attribute to "time_bounds" variable') + ierr=pio_put_att (hist_file, this%tbndid, 'calendar', trim(calendar)) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "calendar" attribute to "time_bounds" variable') + + ! + ! Character + ! + ierr=pio_def_var (hist_file,'date_written',pio_char,dimenchar,this%date_writtenid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "date_written" variable') + ierr=pio_def_var (hist_file,'time_written',pio_char,dimenchar,this%time_writtenid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "time_written" variable') + ! + ! Integer header + ! + ! Define base day variables + ierr=pio_def_var (hist_file,'ndbase',PIO_INT,this%ndbaseid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "ndbase" variable') + str = 'base day' + ierr=pio_put_att (hist_file, this%ndbaseid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "long_name" attribute to "ndbase" variable') + + ierr=pio_def_var (hist_file,'nsbase',PIO_INT,this%nsbaseid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "nsbase" variable') + str = 'seconds of base day' + ierr=pio_put_att (hist_file, this%nsbaseid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "long_name" attribute to "nsbase" variable') + + ierr=pio_def_var (hist_file,'nbdate',PIO_INT,this%nbdateid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "nbdate" variable') + str = 'base date (YYYYMMDD)' + ierr=pio_put_att (hist_file, this%nbdateid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "long_name" attribtue to "nbdate" variable') + + ierr=pio_def_var (hist_file,'nbsec',PIO_INT,this%nbsecid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "nbsec" variable') + str = 'seconds of base date' + ierr=pio_put_att (hist_file, this%nbsecid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "long_name" attribute to "nbsec" variable') + + ierr=pio_def_var (hist_file,'mdt',PIO_INT,this%mdtid) + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to define "mdt" variable') + ierr=pio_put_att (hist_file, this%mdtid, 'long_name', 'timestep') + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "long_name" attribute to "mdt" variable') + ierr=pio_put_att (hist_file, this%mdtid, 'units', 's') + call cam_pio_handle_error(ierr, 'config_define_additional_header_info: failed to add "units" attribute to "mdt" variable') + + end subroutine config_define_additional_header_info + + ! ======================================================================== + + subroutine config_define_fields(this, hist_file, file_index, split_file, restart, mdimids, & + dimindex, timdim, num_patches, header_info, varid_set, num_hdims, fdims) + use string_utils, only: stringify + use pio, only: pio_def_var, pio_double, pio_real, pio_int, pio_put_att + use cam_pio_utils, only: cam_pio_handle_error, cam_pio_def_var + use cam_history_support, only: max_chars + use cam_abortutils, only: endrun + use cam_grid_support, only: cam_grid_header_info_t + use shr_kind_mod, only: r4 => shr_kind_r4 + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + integer, intent(in) :: file_index + logical, intent(in) :: split_file + logical, intent(in) :: restart + integer, intent(in) :: mdimids(:) + integer, intent(inout) :: dimindex(:) + integer, intent(in) :: timdim + integer, intent(in) :: num_patches + logical, intent(in) :: varid_set + type(cam_grid_header_info_t), allocatable, intent(in) :: header_info(:) + integer, intent(out) :: num_hdims + integer, intent(out) :: fdims(:) + + ! Local variables + type(var_desc_t) :: varid + integer :: dimenchar(2) + integer :: ierr + integer :: field_index, idx, jdx + integer :: ncreal + integer :: grd + integer :: mdimsize + integer, allocatable :: mdims(:) + integer, parameter :: max_netcdf_len = 256 + integer :: nacsdims(2) ! dimension ids for nacs (used in restart file) + character(len=max_chars) :: str ! character temporary + character(len=512) :: errmsg + character(len=max_chars) :: fname_tmp ! local copy of field name + character(len=max_chars) :: cell_methods + + do field_index = 1, size(this%field_list) + if (split_file) then + if (file_index == accumulated_file_index) then + ! this is the accumulated file of a potentially split + ! history tape - skip instantaneous fields + if (this%field_list(field_index)%accumulate_type() == 'lst') then + cycle + end if + else + ! this is the instantaneous file of a potentially split + ! history tape - skip accumulated fields + if (this%field_list(field_index)%accumulate_type() /= 'lst') then + cycle + end if + end if + end if + + ! Don't write instantaneous fields to the restart history file + if (restart .and. this%field_list(field_index)%accumulate_type() == 'lst') then + cycle + end if + + if (this%precision() == 'REAL32') then + ncreal = pio_real + else if (this%precision() == 'REAL64') then + ncreal = pio_double + end if + mdims = this%field_list(field_index)%dimensions() + mdimsize = size(mdims) + fname_tmp = strip_suffix(this%field_list(field_index)%diag_name()) + ! Ensure that fname_tmp is not longer than the maximum length for a + ! netcdf file + if (len_trim(fname_tmp) > max_netcdf_len) then + ! Endrun if the name is too long + write(errmsg, *) 'config_define_fields: variable name ', trim(fname_tmp), & + ' too long for NetCDF file (len=', stringify((/len(trim(fname_tmp))/)), ' > ', & + stringify((/max_netcdf_len/)), ')' + call endrun(errmsg, file=__FILE__, line=__LINE__) + end if + ! + ! Create variables and atributes for fields written out as columns + ! + ! Find appropriate grid in header_info + if (.not. allocated(header_info)) then + ! Safety check + call endrun('config_define_fields: header_info not allocated', file=__FILE__, line=__LINE__) + end if + grd = -1 + do idx = 1, size(header_info) + if (header_info(idx)%get_gridid() == this%field_list(field_index)%decomp()) then + grd = idx + exit + end if + end do + if (grd < 0) then + write(errmsg, '(a,i0,2a)') 'grid, ',this%field_list(field_index)%decomp(),', not found for ', & + trim(fname_tmp) + call endrun('config_define_fields: '//errmsg, file=__FILE__, line=__LINE__) + end if + num_hdims = header_info(grd)%num_hdims() + do idx = 1, num_hdims + dimindex(idx) = header_info(grd)%get_hdimid(idx) + if (restart) then + nacsdims(idx) = header_info(grd)%get_hdimid(idx) + end if + end do + do idx = 1, num_patches + varid = this%file_varids(field_index, idx) + ! Figure the dimension ID array for this field + ! We have defined the horizontal grid dimensions in dimindex + fdims(field_index) = num_hdims + do jdx = 1, mdimsize + fdims(field_index) = fdims(field_index) + 1 + dimindex(fdims(field_index)) = mdimids(mdims(jdx)) + end do + if(.not. restart) then + ! Only add time dimension if this is not a restart history tape + fdims(field_index) = fdims(field_index) + 1 + dimindex(fdims(field_index)) = timdim + end if + ! peverwhee - TODO: enable patch output + ! Define the variable + call cam_pio_def_var(hist_file, trim(fname_tmp), ncreal, & + dimindex(1:fdims(field_index)), varid) + if (.not. varid_set) then + this%file_varids(field_index, idx) = varid + end if + + if (mdimsize > 0) then + ierr = pio_put_att(hist_file, varid, 'mdims', mdims(1:mdimsize)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define mdims for '//trim(fname_tmp)) + end if + str = this%field_list(field_index)%sampling_sequence() + if (len_trim(str) > 0) then + ierr = pio_put_att(hist_file, varid, 'Sampling_Sequence', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define Sampling_Sequence for '//trim(fname_tmp)) + end if + if (this%field_list(field_index)%flag_xyfill()) then + ! Write _FillValue and missing_value attributes so that + ! netCDF-aware tools (ncview, ncdump, etc.) recognise fill + ! values. The attribute type must match the variable type. + if (ncreal == pio_double) then + ierr = pio_put_att(hist_file, varid, '_FillValue', this%field_list(field_index)%fill_value()) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define _FillValue for '//trim(fname_tmp)) + ierr = pio_put_att(hist_file, varid, 'missing_value', this%field_list(field_index)%fill_value()) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define missing_value for '//trim(fname_tmp)) + else + ierr = pio_put_att(hist_file, varid, '_FillValue', real(this%field_list(field_index)%fill_value(), r4)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define _FillValue for '//trim(fname_tmp)) + ierr = pio_put_att(hist_file, varid, 'missing_value', real(this%field_list(field_index)%fill_value(), r4)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define missing_value for '//trim(fname_tmp)) + end if + end if + str = this%field_list(field_index)%units() + if (len_trim(str) > 0) then + ierr=pio_put_att (hist_file, varid, 'units', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define units for '//trim(fname_tmp)) + end if + str = this%field_list(field_index)%mixing_ratio() + if (len_trim(str) > 0) then + ierr=pio_put_att (hist_file, varid, 'mixing_ratio', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define mixing_ratio for '//trim(fname_tmp)) + end if + str = this%field_list(field_index)%long_name() + ierr=pio_put_att (hist_file, varid, 'long_name', trim(str)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define long_name for '//trim(fname_tmp)) + ! Assign field attributes defining valid levels and averaging info + cell_methods = '' + if (len_trim(this%field_list(field_index)%cell_methods()) > 0) then + if (len_trim(cell_methods) > 0) then + cell_methods = trim(cell_methods)//' '//trim(this%field_list(field_index)%cell_methods()) + else + cell_methods = trim(this%field_list(field_index)%cell_methods()) + end if + end if + ! Time cell methods is after field method because time averaging is + ! applied later (just before output) than field method which is applied + ! before outfld call. + call AvgflagToString(this%field_list(field_index)%accumulate_type(), str) + cell_methods = adjustl(trim(cell_methods)//' '//'time: '//str) + ierr = pio_put_att(hist_file, varid, 'cell_methods', trim(cell_methods)) + call cam_pio_handle_error(ierr, 'config_define_fields: cannot define cell_methods for '//trim(fname_tmp)) + end do ! end loop over patches + deallocate(mdims) + end do ! end loop over fields + + end subroutine config_define_fields + + ! ======================================================================== + + subroutine config_define_restart_fields(this, num_patches, dimindex, num_hdims, fdims) + use cam_history_support, only: max_chars + use cam_pio_utils, only: cam_pio_def_var + use pio, only: pio_int, pio_double + use cam_abortutils, only: endrun + class(hist_file_t), intent(inout) :: this + integer, intent(in) :: num_patches + integer, intent(in) :: dimindex(:) + integer, intent(in) :: num_hdims + integer, intent(in) :: fdims(:) + ! Local variables + integer :: ierr, field_index, idx + integer :: nstdev, stdev_index + character(len=512) :: errmsg + character(len=max_chars) :: str + type(var_desc_t) :: varid + + ! Determine how many standard deviation fields we have + nstdev = 0 + stdev_index = 1 + do field_index = 1, size(this%field_list) + if (this%field_list(field_index)%accumulate_type() == 'var') then + nstdev = nstdev + 1 + end if + end do + + if (.not. allocated(this%nacs_varids)) then + allocate(this%nacs_varids(size(this%field_list), num_patches), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_define_restart_fields: failed to allocate this%nacs_varids; errmsg = '//trim(errmsg)) + end if + allocate(this%varbuff_varids(nstdev, num_patches), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('config_define_restart_fields: failed to allocate this%varbuff_varids; errmsg = '//trim(errmsg)) + end if + end if + do field_index = 1, size(this%field_list) + do idx = 1, num_patches + if (this%field_list(field_index)%accumulate_type() == 'lst') then + cycle + end if + ! Define "number of samples" variable per field + str = trim(this%field_list(field_index)%diag_name()) // '_nacs' + call cam_pio_def_var(this%restart_file, trim(str), pio_int, & + dimindex(1:num_hdims), varid) + this%nacs_varids(field_index, idx) = varid + ! Define standard deviation buffer variable where relevant + if (this%field_list(field_index)%accumulate_type() == 'var') then + str = trim(this%field_list(field_index)%diag_name()) // '_var_buffer' + call cam_pio_def_var(this%restart_file, trim(str), pio_double, & + dimindex(1:fdims(field_index)), varid) + this%varbuff_varids(stdev_index, idx) = varid + stdev_index = stdev_index + 1 + end if + end do + end do + + end subroutine config_define_restart_fields + + ! ======================================================================== + + subroutine config_write_invariant_header(this, hist_file, file_index, ndbase, & + nsbase, nbdate, nbsec, restart) + use pio, only: pio_put_var + use cam_pio_utils, only: cam_pio_handle_error + use cam_grid_support, only: cam_grid_write_var + use time_manager, only: get_step_size + use cam_history_support, only: write_hist_coord_vars + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + integer, intent(in) :: file_index + integer, intent(in) :: ndbase + integer, intent(in) :: nsbase + integer, intent(in) :: nbdate + integer, intent(in) :: nbsec + logical, intent(in) :: restart + ! Local variables + integer :: idx, ierr + integer :: dtime + + dtime = get_step_size() + + do idx = 1, size(this%grids) + call cam_grid_write_var(hist_file, this%grids(idx), & + file_index=file_index) + end do + ierr = pio_put_var(hist_file, this%mdtid, (/dtime/)) + call cam_pio_handle_error(ierr, 'config_write_invariant_header: cannot put mdt') + + ! + ! Model date info + ! + ierr = pio_put_var(hist_file, this%ndbaseid, (/ndbase/)) + call cam_pio_handle_error(ierr, 'config_write_invariant_header: cannot put ndbase') + ierr = pio_put_var(hist_file, this%nsbaseid, (/nsbase/)) + call cam_pio_handle_error(ierr, 'config_write_invariant_header: cannot put nsbase') + + ierr = pio_put_var(hist_file, this%nbdateid, (/nbdate/)) + call cam_pio_handle_error(ierr, 'config_write_invariant_header: cannot put nbdate') + ierr = pio_put_var(hist_file, this%nbsecid, (/nbsec/)) + call cam_pio_handle_error(ierr, 'config_write_invariant_header: cannot put nbsec') + + ! Write the mdim variable data + call write_hist_coord_vars(hist_file, restart) + + end subroutine config_write_invariant_header + + ! ======================================================================== + + subroutine config_write_restart_file(this) + use pio, only: pio_put_var, pio_file_is_open + use time_manager, only: get_nstep, get_curr_date, get_curr_time + use time_manager, only: set_date_from_time_float, get_step_size + use datetime_mod, only: datetime + use spmd_utils, only: masterproc + use cam_logfile, only: iulog + use perf_mod, only: t_startf, t_stopf + use cam_pio_utils, only: cam_pio_handle_error + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + ! Local variables + integer :: ierr + integer :: field_index, var_buffer_idx + integer :: yr, mon, day ! year, month, and day components of a date + integer :: nstep + integer :: ncdate + integer :: ncsec + integer :: ndcur + integer :: nscur + integer :: start, count1 + integer :: startc(2) ! start values required by nf_put_vara (character) + integer :: countc(2) ! count values required by nf_put_vara (character) + real(r8) :: time ! current (or midpoint) time + real(r8) :: time_interval(2) ! time interval boundaries + logical :: accumulate ! flag for write_field to indicate whether to do accumulation or not + logical :: write_var_buffer + character(len=8) :: cdate ! system date + character(len=8) :: ctime ! system time + + accumulate = .false. + nstep = get_nstep() + call get_curr_date(yr, mon, day, ncsec) + ncdate = yr*10000 + mon*100 + day + call get_curr_time(ndcur, nscur) + time = ndcur + nscur/86400._r8 + time_interval(1) = this%beg_time + time_interval(2) = time + + start = 1 + count1 = 1 + startc(1) = 1 + startc(2) = start + countc(2) = 1 + + call datetime (cdate, ctime) + + ierr = pio_put_var (this%restart_file,this%ndcurid,(/start/),(/count1/),(/ndcur/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "ndcur" variable') + ierr = pio_put_var (this%restart_file,this%nscurid,(/start/),(/count1/),(/nscur/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "nscur" variable') + ierr = pio_put_var (this%restart_file,this%nstephid,(/start/),(/count1/),(/nstep/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "nstephid" variable') + + ierr = pio_put_var (this%restart_file,this%dateid,(/start/),(/count1/),(/ncdate/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "ncdate" variable') + ierr = pio_put_var (this%restart_file,this%datesecid,(/start/),(/count1/),(/ncsec/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "ncsec" variable') + countc(1) = 2 + + ierr=pio_put_var (this%restart_file, this%timeid, (/start/),(/count1/),(/time/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write instantaneous "time" variable') + ierr=pio_put_var (this%restart_file, this%tbndid, startc, countc, time_interval) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "time_bounds" variable') + countc(1) = 8 + ierr = pio_put_var (this%restart_file, this%date_writtenid, startc, countc, (/cdate/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "cdate" variable') + ierr = pio_put_var (this%restart_file, this%time_writtenid, startc, countc, (/ctime/)) + call cam_pio_handle_error(ierr, 'config_write_restart_file: cannot write "ctime" variable') + + var_buffer_idx = 0 + + do field_index = 1, size(this%field_list) + if (this%field_list(field_index)%accumulate_type() == 'lst') then + cycle + end if + call this%write_field(this%restart_file, this%field_list(field_index), start, field_index, & + this%precision(), accumulate) + if (this%field_list(field_index)%accumulate_type() == 'var') then + write_var_buffer = .true. + var_buffer_idx = var_buffer_idx + 1 + else + write_var_buffer = .false. + end if + call this%write_restart_fields(this%restart_file, this%field_list(field_index), field_index, write_var_buffer, & + var_buffer_idx) + ! Clear the buffers + call this%field_list(field_index)%clear_buffers() + end do + + end subroutine config_write_restart_file + + ! ======================================================================== + + subroutine config_write_time_dependent_variables(this) + use pio, only: pio_put_var, pio_file_is_open + use time_manager, only: get_nstep, get_curr_date, get_curr_time + use time_manager, only: set_date_from_time_float, get_step_size + use datetime_mod, only: datetime + use spmd_utils, only: masterproc + use cam_logfile, only: iulog + use perf_mod, only: t_startf, t_stopf + use cam_pio_utils, only: cam_pio_handle_error + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + + ! Local variables + integer :: yr, mon, day ! year, month, and day components of a date + integer :: yr_mid, mon_mid, day_mid ! year, month, and day components of midpoint date + integer :: nstep ! current timestep number + integer :: ncdate(max_split_files) ! current (or midpoint) date in integer format [yyyymmdd] + integer :: ncsec(max_split_files) ! current (or midpoint) time of day [seconds] + integer :: ndcur ! day component of current time + integer :: nscur ! seconds component of current time + real(r8) :: time ! current (or midpoint) time + real(r8) :: time_interval(2) ! time interval boundaries + integer :: ierr + integer :: split_file_index, field_idx + integer :: start, count1 + integer :: startc(2) ! start values required by nf_put_vara (character) + integer :: countc(2) ! count values required by nf_put_vara (character) + logical :: accumulate ! flag for write_field to indicate whether to do accumulation or not + logical :: is_initfile, is_satfile + character(len=8) :: cdate ! system date + character(len=8) :: ctime ! system time + + accumulate = .true. + nstep = get_nstep() + call get_curr_date(yr, mon, day, ncsec(instantaneous_file_index)) + ncdate(instantaneous_file_index) = yr*10000 + mon*100 + day + call get_curr_time(ndcur, nscur) + time = ndcur + nscur/86400._r8 + time_interval(1) = this%beg_time + time_interval(2) = time + call set_date_from_time_float((time_interval(1) + time_interval(2)) / 2._r8, & + yr_mid, mon_mid, day_mid, ncsec(accumulated_file_index)) ncdate(accumulated_file_index) = yr_mid*10000 + mon_mid*100 + day_mid ! Increment samples @@ -1440,7 +2182,7 @@ subroutine config_write_time_dependent_variables(this, restart) startc(2) = start countc(2) = 1 - if(.not.restart) this%beg_time = time ! update beginning time of next interval + this%beg_time = time ! update beginning time of next interval call datetime (cdate, ctime) ! peverwhee - TODO handle composed fields @@ -1472,7 +2214,7 @@ subroutine config_write_time_dependent_variables(this, restart) ierr = pio_put_var (this%hist_files(split_file_index),this%datesecid,(/start/),(/count1/),(/ncsec(split_file_index)/)) call cam_pio_handle_error(ierr, 'config_write_time_dependent_variables: cannot write "ncsec" variable') countc(1) = 2 - if (split_file_index == accumulated_file_index .and. .not. restart .and. .not. is_initfile) then + if (split_file_index == accumulated_file_index .and. .not. is_initfile) then ! accumulated tape - time is midpoint of time_bounds ierr=pio_put_var (this%hist_files(split_file_index), this%timeid, (/start/),(/count1/), & @@ -1494,14 +2236,16 @@ subroutine config_write_time_dependent_variables(this, restart) ! we may have a history split, conditionally skip fields that are ! for the other file if ((this%field_list(field_idx)%accumulate_type() .eq. 'lst') .and. & - split_file_index == accumulated_file_index .and. .not. restart) then + split_file_index == accumulated_file_index) then cycle else if ((this%field_list(field_idx)%accumulate_type() .ne. 'lst') .and. & - split_file_index == instantaneous_file_index .and. .not. restart) then + split_file_index == instantaneous_file_index) then cycle end if - call this%write_field(this%field_list(field_idx), split_file_index, restart, start, field_idx, & - this%precision()) + call this%write_field(this%hist_files(split_file_index), this%field_list(field_idx), start, field_idx, & + this%precision(), accumulate) + ! Clear the buffers + call this%field_list(field_idx)%clear_buffers() end do end do call t_stopf ('write_field') @@ -1510,10 +2254,11 @@ end subroutine config_write_time_dependent_variables ! ======================================================================== - subroutine config_write_field(this, field, split_file_index, restart, & - sample_index, field_index, field_precision) + subroutine config_write_field(this, hist_file, field, sample_index, field_index, field_precision, & + accumulate) use pio, only: PIO_OFFSET_KIND, pio_setframe use hist_api, only: hist_field_norm_value + use hist_api, only: hist_field_value use cam_grid_support, only: cam_grid_write_dist_array use cam_abortutils, only: check_allocate, endrun use hist_field, only: hist_field_info_t @@ -1522,12 +2267,12 @@ subroutine config_write_field(this, field, split_file_index, restart, & use cam_logfile, only: iulog ! Dummy arguments class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file type(hist_field_info_t), intent(inout) :: field - integer, intent(in) :: split_file_index - logical, intent(in) :: restart integer, intent(in) :: sample_index integer, intent(in) :: field_index character(len=*), intent(in) :: field_precision + logical, intent(in) :: accumulate ! Local variables integer, allocatable :: field_shape(:) ! Field file dim sizes @@ -1573,45 +2318,133 @@ subroutine config_write_field(this, field, split_file_index, restart, & do patch_idx = 1, num_patches varid = this%file_varids(field_index, patch_idx) - call pio_setframe(this%hist_files(split_file_index), varid, int(sample_index,kind=PIO_OFFSET_KIND)) + call pio_setframe(hist_file, varid, int(sample_index,kind=PIO_OFFSET_KIND)) if (frank == 1) then - call hist_field_norm_value(field, field_data(:,1), logger=errors) + if (accumulate) then + call hist_field_norm_value(field, field_data(:,1), logger=errors) + else + call hist_field_value(field, field_data(:,1), logger=errors) + end if if (errors%num_errors() > 0) then call errors%output(iulog) write(errmsg, *) subname, 'ERROR writing field "', trim(field%diag_name()), '"' call endrun(errmsg) end if if (trim(field_precision) == 'REAL32') then - call cam_grid_write_dist_array(this%hist_files(split_file_index), field_decomp, (/dim_sizes(1)/), & + call cam_grid_write_dist_array(hist_file, field_decomp, (/dim_sizes(1)/), & field_shape, real(field_data(:,1), r4), varid) else - call cam_grid_write_dist_array(this%hist_files(split_file_index), field_decomp, (/dim_sizes(1)/), & + call cam_grid_write_dist_array(hist_file, field_decomp, (/dim_sizes(1)/), & field_shape, field_data(:,1), varid) end if else - call hist_field_norm_value(field, field_data, logger=errors) + if (accumulate) then + call hist_field_norm_value(field, field_data, logger=errors) + else + call hist_field_value(field, field_data, logger=errors) + end if if (errors%num_errors() > 0) then call errors%output(iulog) write(errmsg, *) subname, 'ERROR writing field "', trim(field%diag_name()), '"' call endrun(errmsg) end if if (trim(field_precision) == 'REAL32') then - call cam_grid_write_dist_array(this%hist_files(split_file_index), field_decomp, dim_sizes(1:frank), & + call cam_grid_write_dist_array(hist_file, field_decomp, dim_sizes(1:frank), & field_shape, real(field_data, r4), varid) else - call cam_grid_write_dist_array(this%hist_files(split_file_index), field_decomp, dim_sizes(1:frank), & + call cam_grid_write_dist_array(hist_file, field_decomp, dim_sizes(1:frank), & field_shape, field_data, varid) end if end if end do - ! Clear the buffers after writing - call field%clear_buffers(logger=errors) - end subroutine config_write_field ! ======================================================================== + subroutine config_write_restart_fields(this, hist_file, field, field_index, write_var_buffer, & + var_buff_idx) + use hist_field, only: hist_field_info_t + use hist_api, only: hist_field_samples + use hist_api, only: hist_field_var_buffer + use cam_abortutils, only: check_allocate + use pio, only: var_desc_t + use cam_grid_support, only: cam_grid_write_dist_array + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + type(file_desc_t), intent(inout) :: hist_file + type(hist_field_info_t), intent(inout) :: field + integer, intent(in) :: field_index + integer, intent(in) :: var_buff_idx + logical, intent(in) :: write_var_buffer + ! Local variables + integer, allocatable :: num_samples(:) + integer, allocatable :: beg_dims(:) + integer, allocatable :: end_dims(:) + integer, allocatable :: dimind(:) + integer, allocatable :: dim_sizes(:) + integer, allocatable :: field_shape(:) + real(r8), allocatable :: var_buffer(:,:) + type(var_desc_t) :: varid + integer :: idx, patch_index + integer :: field_decomp + integer :: num_patches, ierr, rank + character(len=512) :: errmsg + character(len=*), parameter :: subname = 'config_write_restart_fields: ' + + beg_dims = field%beg_dims() + end_dims = field%end_dims() + field_shape = field%shape() + rank = size(field_shape) + + allocate(num_samples(end_dims(1) - beg_dims(1) + 1), stat=ierr, errmsg=errmsg) + call check_allocate(ierr, subname, 'num_samples', file=__FILE__, line=__LINE__-1, errmsg=errmsg) + + if (rank == 1) then + allocate(var_buffer(end_dims(1) - beg_dims(1) + 1, 1), stat=ierr, errmsg=errmsg) + call check_allocate(ierr, subname, 'var_buffer', file=__FILE__, line=__LINE__-1, errmsg=errmsg) + else + allocate(var_buffer(end_dims(1) - beg_dims(1) + 1, field_shape(2)), stat=ierr, errmsg=errmsg) + call check_allocate(ierr, subname, 'var_buffer', file=__FILE__, line=__LINE__-1, errmsg=errmsg) + end if + + ! Shape of array + dimind = field%dimensions() + + allocate(dim_sizes(size(beg_dims)), stat=ierr) + call check_allocate(ierr, subname, 'dim_sizes', file=__FILE__, line=__LINE__-1) + do idx = 1, size(beg_dims) + dim_sizes(idx) = end_dims(idx) - beg_dims(idx) + 1 + end do + + field_decomp = field%decomp() + + num_patches = 1 + + do patch_index = 1, num_patches + varid = this%nacs_varids(field_index, patch_index) + call hist_field_samples(field, num_samples) + call cam_grid_write_dist_array(hist_file, field_decomp, (/dim_sizes(1)/), & + (/field_shape(1)/), num_samples, varid) + if (write_var_buffer) then + varid = this%varbuff_varids(var_buff_idx, patch_index) + if (rank == 1) then + call hist_field_var_buffer(field, var_buffer(:,1)) + call cam_grid_write_dist_array(hist_file, field_decomp, (/dim_sizes(1)/), & + (/field_shape(1)/), var_buffer(:,1), varid) + else + call hist_field_var_buffer(field, var_buffer) + call cam_grid_write_dist_array(hist_file, field_decomp, dim_sizes(1:rank), & + field_shape, var_buffer, varid) + end if + end if + end do + + + end subroutine config_write_restart_fields + + ! ======================================================================== + subroutine config_close_files(this) use pio, only: pio_file_is_open use cam_pio_utils, only: cam_pio_closefile @@ -1623,16 +2456,16 @@ subroutine config_close_files(this) ! Local variables integer :: split_file_index - if(pio_file_is_open(this%hist_files(accumulated_file_index)) .or. & - pio_file_is_open(this%hist_files(instantaneous_file_index))) then - deallocate(this%file_varids) - end if +! if(pio_file_is_open(this%hist_files(accumulated_file_index)) .or. & +! pio_file_is_open(this%hist_files(instantaneous_file_index))) then +! deallocate(this%file_varids) +! end if do split_file_index = 1, max_split_files if (pio_file_is_open(this%hist_files(split_file_index))) then if (masterproc) then write(iulog,*)'config_close_files: nf_close(', trim(this%volume),')=',& - this%file_names(split_file_index) + trim(this%file_names(split_file_index)) end if call cam_pio_closefile(this%hist_files(split_file_index)) end if @@ -1647,6 +2480,26 @@ end subroutine config_close_files ! ======================================================================== + subroutine config_close_restart_file(this) + use cam_pio_utils, only: cam_pio_closefile + use pio, only: pio_file_is_open + use spmd_utils, only: masterproc + use cam_logfile, only: iulog + ! Dummy arguments + class(hist_file_t), intent(inout) :: this + + if (pio_file_is_open(this%restart_file)) then + if (masterproc) then + write(iulog,*) 'config_close_restart_file: nf_close(', trim(this%volume), ')=', & + trim(this%restart_file_name) + end if + call cam_pio_closefile(this%restart_file) + end if + + end subroutine config_close_restart_file + + ! ======================================================================== + subroutine config_clear_buffers(this) use hist_msg_handler, only: hist_log_messages use spmd_utils, only: masterproc @@ -1698,7 +2551,8 @@ end function count_array ! ======================================================================== subroutine read_namelist_entry(unitn, hfile_config, hist_inst_fields, & - hist_avg_fields, hist_min_fields, hist_max_fields, hist_var_fields) + hist_avg_fields, hist_min_fields, hist_max_fields, hist_var_fields, & + num_fields) use mpi, only: MPI_CHARACTER, MPI_INTEGER, MPI_LOGICAL use string_utils, only: stringify use spmd_utils, only: masterproc, masterprocid, mpicom @@ -1720,6 +2574,7 @@ subroutine read_namelist_entry(unitn, hfile_config, hist_inst_fields, & character(len=max_fldlen), allocatable, intent(inout) :: hist_min_fields(:) character(len=max_fldlen), allocatable, intent(inout) :: hist_max_fields(:) character(len=max_fldlen), allocatable, intent(inout) :: hist_var_fields(:) + integer, intent(out) :: num_fields ! Local variables (namelist) character(len=vlen) :: hist_volume character(len=vlen) :: hist_precision @@ -1868,6 +2723,8 @@ subroutine read_namelist_entry(unitn, hfile_config, hist_inst_fields, & interp_nlat=hist_interp_nlat, interp_nlon=hist_interp_nlon, & interp_grid=hist_interp_grid, interp_type=hist_interp_type) call hfile_config%print_config() + num_fields = num_fields_inst + num_fields_avg + num_fields_min + & + num_fields_max + num_fields_var end subroutine read_namelist_entry @@ -1974,7 +2831,7 @@ end subroutine allocate_field_arrays ! ======================================================================== - subroutine hist_read_namelist_config(filename, config_arr) + subroutine hist_read_namelist_config(filename, config_arr, max_fields) use mpi, only: MPI_CHARACTER, MPI_INTEGER use shr_kind_mod, only: max_str =>SHR_KIND_CXX, CM => shr_kind_cm use shr_nl_mod, only: shr_nl_find_group_name @@ -1988,6 +2845,7 @@ subroutine hist_read_namelist_config(filename, config_arr) ! Dummy arguments character(len=*), intent(in) :: filename type(hist_file_t), allocatable, intent(inout) :: config_arr(:) + integer, intent(out) :: max_fields ! Local variables integer :: unitn integer :: read_status @@ -1995,6 +2853,7 @@ subroutine hist_read_namelist_config(filename, config_arr) integer :: line_num integer :: lindex integer :: num_configs + integer :: num_fields logical :: filefound character(len=max_fldlen), allocatable :: hist_inst_fields(:) character(len=max_fldlen), allocatable :: hist_avg_fields(:) @@ -2009,6 +2868,7 @@ subroutine hist_read_namelist_config(filename, config_arr) ! Variables for reading a namelist entry unitn = -1 ! Prevent reads on error or wrong tasks ierr = 0 + max_fields = 0 errmsg = '' if (masterproc) then @@ -2083,7 +2943,10 @@ subroutine hist_read_namelist_config(filename, config_arr) end if call read_namelist_entry(unitn, config_arr(lindex), & hist_inst_fields, hist_avg_fields, hist_min_fields, & - hist_max_fields, hist_var_fields) + hist_max_fields, hist_var_fields, num_fields) + if (num_fields > max_fields) then + max_fields = num_fields + end if end do ! ! Cleanup diff --git a/src/history/cam_hist_restart.F90 b/src/history/cam_hist_restart.F90 new file mode 100644 index 000000000..42057da28 --- /dev/null +++ b/src/history/cam_hist_restart.F90 @@ -0,0 +1,442 @@ +module cam_hist_restart + use pio, only: var_desc_t + use cam_history_support, only: max_fieldname_len + use shr_kind_mod, only: r4 => shr_kind_r4 + use shr_kind_mod, only: r8 => shr_kind_r8 + use cam_logfile, only: iulog + + implicit none + private + + integer, parameter :: max_dimensions = 4 + + type restart_variable_t + type(var_desc_t), pointer :: vdesc => null() + integer :: var_type + integer :: number_of_dimensions + integer :: dimension_ids(max_dimensions) + character(len=max_fieldname_len) :: var_name + logical :: fill_set = .false. + integer :: integer_fill + real(r4) :: real_fill + real(r8) :: double_fill + end type restart_variable_t + + type restart_dimension_t + integer :: dimension_length + integer :: dimension_id + character(len=max_fieldname_len) :: dimension_name + end type restart_dimension_t + + ! + ! The size of these parameters should match the assignments in restart_vars_setnames and restart_dims_setnames below + ! + integer, parameter :: num_restart_vars = 21 + integer, parameter :: num_restart_dims = 9 + type(restart_variable_t) :: restart_vars(num_restart_vars) + type(restart_dimension_t) :: restart_dims(num_restart_dims) + + integer, parameter :: max_num_fields = 1000 + + integer, parameter :: num_configs_dim_ind = 1 + integer, parameter :: max_string_len_dim_ind = 2 + integer, parameter :: max_fieldname_len_dim_ind = 3 + integer, parameter :: max_num_fields_dim_ind = 4 + integer, parameter :: max_chars_dim_ind = 5 + integer, parameter :: max_dims_dim_ind = 6 + integer, parameter :: registeredmdims_dim_ind = 7 + integer, parameter :: max_hcoordname_len_dim_ind = 8 + integer, parameter :: max_num_split_files_dim_ind = 9 + + private :: set_restart_variable_names + private :: set_restart_dimension_names + + public :: hist_restart_init + public :: hist_restart_write + +CONTAINS + + subroutine hist_restart_init(restart_file, num_hist_configs, max_fields) + use pio, only: file_desc_t, pio_def_var + use cam_pio_utils, only: cam_pio_handle_error, cam_pio_def_dim + type(file_desc_t), intent(inout) :: restart_file + integer, intent(in) :: num_hist_configs + integer, intent(in) :: max_fields + + ! Local variables + integer :: idx, ndims, kdx, ierr + integer :: dimids(4) + + ! Set the restart dimensions and variables for writing to the file + call set_restart_variable_names() + call set_restart_dimension_names(num_hist_configs, max_fields) + + do idx = 1, num_restart_dims + ! it's possible that one or more of these have been defined elsewhere + call cam_pio_def_dim(restart_file, restart_dims(idx)%dimension_name, restart_dims(idx)%dimension_length, & + restart_dims(idx)%dimension_id, existOK=.true.) + end do + + do idx = 1, num_restart_vars + ndims= restart_vars(idx)%number_of_dimensions + do kdx = 1, ndims + dimids(kdx)=restart_dims(restart_vars(idx)%dimension_ids(kdx))%dimension_id + end do + allocate(restart_vars(idx)%vdesc) + ierr = pio_def_var(restart_file, restart_vars(idx)%var_name, restart_vars(idx)%var_type, dimids(1:ndims), restart_vars(idx)%vdesc) + call cam_pio_handle_error(ierr, 'INIT_RESTART_HISTORY: Error defining '//trim(restart_vars(idx)%var_name)) + end do + + end subroutine hist_restart_init + + subroutine hist_restart_write(restart_file, hist_configs, max_num_fields, just_written) + use pio, only: file_desc_t, pio_put_var + use cam_hist_file, only: hist_file_t + use cam_history_support, only: max_chars, max_string_len, get_hist_coord_names, registeredmdims + use cam_grid_support, only: max_split_files, max_hcoordname_len + use cam_abortutils, only: endrun + type(file_desc_t), intent(inout) :: restart_file + type(hist_file_t), intent(in) :: hist_configs(:) + integer, intent(in) :: max_num_fields + logical, intent(in) :: just_written(:) + + ! Local variables + integer :: idx, jdx, ierr + integer :: has_rh_int(size(hist_configs)) + integer :: num_fields(size(hist_configs)) + integer :: num_frames(size(hist_configs)) + integer :: max_frames(size(hist_configs)) + integer :: ndims(max_num_fields) + integer :: decomp(max_num_fields, size(hist_configs)) + integer :: num_levels(max_num_fields, size(hist_configs)) + integer :: fill_flag(max_num_fields, size(hist_configs)) + integer :: dimensions(max_dimensions, max_num_fields, size(hist_configs)) + character(len=max_fieldname_len) :: field_list(max_num_fields, size(hist_configs)) + character(len=max_fieldname_len), allocatable :: field_list_config(:) + character(len=max_chars) :: output_freq(size(hist_configs)) + character(len=max_string_len) :: current_files(size(hist_configs), max_split_files) + character(len=max_chars) :: hist_precision(size(hist_configs)) + character(len=max_chars) :: avg_flag(max_num_fields, size(hist_configs)) + character(len=max_chars) :: long_name(max_num_fields, size(hist_configs)) + character(len=max_chars) :: cell_methods(max_num_fields, size(hist_configs)) + character(len=max_chars) :: units(max_num_fields, size(hist_configs)) + character(len=max_chars) :: volume(size(hist_configs)) + character(len=max_hcoordname_len) :: dim_names(registeredmdims) + character(len=max_string_len) :: restart_file_paths(size(hist_configs)) + real(r8) :: beg_time(size(hist_configs)) + real(r8) :: fill_value(max_num_fields, size(hist_configs)) + logical :: has_accum + + field_list = '' + has_rh_int = 0 + avg_flag = '' + decomp = 0 + num_levels = 0 + cell_methods = '' + long_name = '' + units = '' + fill_flag = 0 + fill_value = 0 + dimensions = 0 + restart_file_paths = '' + ! Compile all the necessary info for the restart file from the hist_configs array + do idx = 1, size(hist_configs) + ! Grab the field list + field_list_config = hist_configs(idx)%get_field_list() + num_fields(idx) = hist_configs(idx)%get_num_fields() + field_list(1:num_fields(idx), idx) = field_list_config(1:num_fields(idx)) + ! Determine whether or not there will be an rh file + if (hist_configs(idx)%has_accumulated_fields() .and. .not. just_written(idx)) then + has_rh_int(idx) = 1 + restart_file_paths(idx) = hist_configs(idx)%get_restart_filename() + end if + volume(idx) = hist_configs(idx)%get_volume() + ! Get the output frequency + output_freq(idx) = hist_configs(idx)%output_freq() + ! Get the number of samples written to the hist file + num_frames(idx) = hist_configs(idx)%get_num_samples() + ! Get the maximum number of samples written to the hist file + max_frames(idx) = hist_configs(idx)%max_frame() + ! Get the file names + current_files(idx,:) = hist_configs(idx)%get_filenames() + ! Get the precision + hist_precision(idx) = hist_configs(idx)%precision() + ! Get the interval start time + beg_time(idx) = hist_configs(idx)%get_beg_time() + ! Get field-specific info vv + ! Get accumulated flags for each field + avg_flag(1:num_fields(idx), idx) = hist_configs(idx)%get_averaging_flags() + ! Get field decompositions + decomp(1:num_fields(idx), idx) = hist_configs(idx)%get_decompositions() + ! Get num vertical levels + num_levels(1:num_fields(idx), idx) = hist_configs(idx)%get_num_levels() + ! Get cell methods + cell_methods(1:num_fields(idx), idx) = hist_configs(idx)%get_cell_methods() + ! Get long names + long_name(1:num_fields(idx), idx) = hist_configs(idx)%get_long_names() + ! Get units + units(1:num_fields(idx), idx) = hist_configs(idx)%get_units() + ! Get fill value flag and fill value + fill_flag(1:num_fields(idx), idx) = hist_configs(idx)%get_fill_flags() + fill_value(1:num_fields(idx), idx) = hist_configs(idx)%get_fill_values() + ! Get field dimension indices + ndims(1:num_fields(idx)) = hist_configs(idx)%get_num_dimensions() + do jdx = 1, num_fields(idx) + dimensions(1:ndims(jdx), 1:num_fields(idx), idx) = hist_configs(idx)%get_dimension_indices() + end do + ! End field specific info ^^ + end do + + ! Grab the dimension names + dim_names = get_hist_coord_names() + + ! Loop over the restart vars and write them to the file + do idx = 1, num_restart_vars + select case(trim(restart_vars(idx)%var_name)) + case ('has_rh_file') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, has_rh_int) + case ('volume') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, volume) + case ('output_frequency') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, output_freq) + case ('field_list') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, field_list) + case ('number_of_fields') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, num_fields) + case ('number_of_frames') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, num_frames) + case ('max_frames') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, max_frames) + case ('current_file_name') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, current_files) + case ('precision') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, hist_precision) + case ('interval_start_time') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, beg_time) + case ('field_average_flag') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, avg_flag) + case ('field_decomposition_type') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, decomp) + case ('field_num_vertical_levels') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, num_levels) + case ('field_cell_methods') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, cell_methods) + case ('field_long_name') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, long_name) + case ('field_units') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, units) + case ('field_fill_flag') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, fill_flag) + case ('field_fill_value') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, fill_value) + case ('field_dimensions') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, dimensions) + case ('dimension_names') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, dim_names) + case ('rh_file_path') + ierr = pio_put_var(restart_file, restart_vars(idx)%vdesc, restart_file_paths) + case default + end select + if (ierr /= 0) then + call endrun('hist_restart_write: failed to write variable '//restart_vars(idx)%var_name) + end if + end do + end subroutine hist_restart_write + + subroutine set_restart_variable_names() + use pio, only: pio_int, pio_double, pio_char + integer :: rvar_index + + rvar_index = 1 + restart_vars(rvar_index)%var_name = 'has_rh_file' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 1 + restart_vars(rvar_index)%dimension_ids(1) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'output_frequency' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'volume' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_list' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_fieldname_len_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'number_of_fields' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 1 + restart_vars(rvar_index)%dimension_ids(1) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'number_of_frames' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 1 + restart_vars(rvar_index)%dimension_ids(1) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'max_frames' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 1 + restart_vars(rvar_index)%dimension_ids(1) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'current_file_name' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_string_len_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = max_num_split_files_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'precision' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'interval_start_time' + restart_vars(rvar_index)%var_type = pio_double + restart_vars(rvar_index)%number_of_dimensions = 1 + restart_vars(rvar_index)%dimension_ids(1) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_average_flag' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_decomposition_type' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + restart_vars(rvar_index)%fill_set = .true. + restart_vars(rvar_index)%integer_fill = 0 + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_num_vertical_levels' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + restart_vars(rvar_index)%fill_set = .true. + restart_vars(rvar_index)%integer_fill = 0 + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_cell_methods' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_long_name' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_units' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_chars_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_fill_flag' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_fill_value' + restart_vars(rvar_index)%var_type = pio_double + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + restart_vars(rvar_index)%fill_set = .true. + restart_vars(rvar_index)%double_fill = 0.0_r8 + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'field_dimensions' + restart_vars(rvar_index)%var_type = pio_int + restart_vars(rvar_index)%number_of_dimensions = 3 + restart_vars(rvar_index)%dimension_ids(1) = max_dims_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = max_num_fields_dim_ind + restart_vars(rvar_index)%dimension_ids(3) = num_configs_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'dimension_names' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_hcoordname_len_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = registeredmdims_dim_ind + + rvar_index = rvar_index + 1 + restart_vars(rvar_index)%var_name = 'rh_file_path' + restart_vars(rvar_index)%var_type = pio_char + restart_vars(rvar_index)%number_of_dimensions = 2 + restart_vars(rvar_index)%dimension_ids(1) = max_string_len_dim_ind + restart_vars(rvar_index)%dimension_ids(2) = num_configs_dim_ind + + end subroutine set_restart_variable_names + + subroutine set_restart_dimension_names(num_hist_configs, max_fields) + use cam_history_support, only: max_string_len, max_chars, max_fieldname_len, registeredmdims + use cam_grid_support, only: max_hcoordname_len, max_split_files + integer, intent(in) :: num_hist_configs + integer, intent(in) :: max_fields + + restart_dims(num_configs_dim_ind)%dimension_name = 'number_of_hist_configs' + restart_dims(num_configs_dim_ind)%dimension_length = num_hist_configs + + restart_dims(max_string_len_dim_ind)%dimension_name = 'max_string_length' + restart_dims(max_string_len_dim_ind)%dimension_length = max_string_len + + restart_dims(max_fieldname_len_dim_ind)%dimension_name = 'max_fieldname_length' + restart_dims(max_fieldname_len_dim_ind)%dimension_length = max_fieldname_len + + restart_dims(max_num_fields_dim_ind)%dimension_name = 'max_fields_per_configuration' + restart_dims(max_num_fields_dim_ind)%dimension_length = max_fields + + restart_dims(max_chars_dim_ind)%dimension_name = 'max_chars' + restart_dims(max_chars_dim_ind)%dimension_length = max_chars + + restart_dims(max_dims_dim_ind)%dimension_name = 'max_variable_mdims' + restart_dims(max_dims_dim_ind)%dimension_length = max_dimensions + + restart_dims(registeredmdims_dim_ind)%dimension_name = 'registered_mdims' + restart_dims(registeredmdims_dim_ind)%dimension_length = registeredmdims + + restart_dims(max_hcoordname_len_dim_ind)%dimension_name = 'max_hcoordname_len' + restart_dims(max_hcoordname_len_dim_ind)%dimension_length = max_hcoordname_len + + restart_dims(max_num_split_files_dim_ind)%dimension_name = 'max_num_split_files' + restart_dims(max_num_split_files_dim_ind)%dimension_length = max_split_files + + end subroutine set_restart_dimension_names + +end module cam_hist_restart diff --git a/src/history/cam_history.F90 b/src/history/cam_history.F90 index a023b5574..3f8d9481b 100644 --- a/src/history/cam_history.F90 +++ b/src/history/cam_history.F90 @@ -13,6 +13,7 @@ module cam_history use cam_hist_file, only: hist_file_t use hist_field, only: hist_field_info_t use hist_hash_table, only: hist_hash_table_t + use cam_logfile, only: iulog implicit none private @@ -31,6 +32,8 @@ module cam_history public :: history_add_field ! Write to list of possible history fields for this run public :: history_out_field ! Accumulate field if its in use by one or more tapes public :: history_wrap_up ! Process history files at end of timestep or run + public :: history_restart_init ! Initialize history fields on restart file, if necessary + public :: history_restart_write ! Write restart files, if necessary interface history_out_field module procedure history_out_field_1d @@ -48,6 +51,8 @@ module cam_history type(hist_field_info_t), pointer :: possible_field_list_head type(hist_hash_table_t) :: possible_field_list integer :: num_possible_fields + logical, allocatable :: just_written(:) + integer :: max_num_fields CONTAINS @@ -57,13 +62,23 @@ subroutine history_readnl(nlfile) ! Purpose: Read in history namelist and set hist_configs ! !----------------------------------------------------------------------- - use cam_hist_file, only: hist_read_namelist_config + use cam_hist_file, only: hist_read_namelist_config + use cam_abortutils, only: endrun ! Dummy argument character(len=*), intent(in) :: nlfile ! path of namelist input file + integer :: ierr + character(len=256) :: alloc_errmsg + character(len=512) :: errmsg ! Read in CAM history configuration - call hist_read_namelist_config(nlfile, hist_configs) + call hist_read_namelist_config(nlfile, hist_configs, max_num_fields) + allocate(just_written(size(hist_configs)), stat=ierr, errmsg=alloc_errmsg) + if (ierr /= 0) then + write(errmsg,*) 'history_readnl: failed to allocate hist_configs; errmsg = ', alloc_errmsg + call endrun(errmsg) + end if + just_written = .false. end subroutine history_readnl @@ -111,6 +126,8 @@ subroutine history_write_files() ! peverwhee - TODO: remove when restarts are implemented restart = .false. + just_written = .false. + ! Loop over history volumes do file_idx = 1, size(hist_configs) ! Determine if it's time to write! @@ -210,13 +227,15 @@ subroutine history_write_files() end if end do end do - call hist_configs(file_idx)%define_file(restart, logname, host, model_doi_url) + call hist_configs(file_idx)%define_file(logname, host, model_doi_url) end if - call hist_configs(file_idx)%write_time_dependent_variables(restart) + call hist_configs(file_idx)%write_time_dependent_variables() if (nstep == 0) then ! Reset samples if nstep0 was written call hist_configs(file_idx)%reset_samples() end if + ! Flag that we've just written this volume + just_written(file_idx) = .true. end do end subroutine history_write_files @@ -880,6 +899,61 @@ subroutine history_wrap_up(restart_write, last_timestep) end subroutine history_wrap_up +!####################################################################### + + subroutine history_restart_init(restart_file) + use cam_hist_restart, only: hist_restart_init + use pio, only: file_desc_t + ! Dummy variables + type(file_desc_t), intent(inout) :: restart_file + ! Local variables + integer :: config_idx + + ! Set up history restart variables in the restart file + ! We can skip this if we have no output for this run + if (max_num_fields > 0) then + call hist_restart_init(restart_file, size(hist_configs), max_num_fields) + end if + + ! Create restart files for history configs if necessary (.rhX.) + do config_idx = 1, size(hist_configs) + ! Only write the rhX file if it has accumulated fields and has not just been written + if (.not. hist_configs(config_idx)%has_accumulated_fields() .or. just_written(config_idx)) then + cycle + end if + call hist_configs(config_idx)%define_restart_file(logname, host, model_doi_url) + end do + + end subroutine history_restart_init + +!####################################################################### + + subroutine history_restart_write(restart_file) + use pio, only: file_desc_t + use cam_hist_restart, only: hist_restart_write + ! Dummy variables + type(file_desc_t), intent(inout) :: restart_file + ! Local variables + integer :: config_idx + + if (max_num_fields == 0) then + ! Don't do anything if there aren't any history fields + return + end if + + call hist_restart_write(restart_file, hist_configs, max_num_fields, just_written) + + ! Write restart files for history configs if necessary (.rhX.) + do config_idx = 1, size(hist_configs) + ! Only write the rhX file if it has accumulated fields and has not just been written + if (.not. hist_configs(config_idx)%has_accumulated_fields() .or. just_written(config_idx)) then + cycle + end if + call hist_configs(config_idx)%write_restart_file() + call hist_configs(config_idx)%close_restart_file() + end do + + end subroutine history_restart_write !####################################################################### recursive function get_entry_by_name(listentry, name) result(entry) diff --git a/src/history/cam_history_support.F90 b/src/history/cam_history_support.F90 index c5a26d148..5aea98c70 100644 --- a/src/history/cam_history_support.F90 +++ b/src/history/cam_history_support.F90 @@ -9,6 +9,7 @@ module cam_history_support use shr_kind_mod, only: r8=>shr_kind_r8, shr_kind_cl, shr_kind_cxx use cam_grid_support, only: max_hcoordname_len + use cam_logfile, only: iulog implicit none private @@ -17,13 +18,13 @@ module cam_history_support integer, parameter, public :: fieldname_len = 32 ! max chars for field name integer, parameter, public :: fieldname_suffix_len = 3 ! length of field name suffix ("&IC") ! max_fieldname_len = max chars for field name (including suffix) - integer, parameter, public :: max_fieldname_len = fieldname_len + fieldname_suffix_len + integer, parameter, public :: max_fieldname_len = fieldname_len + fieldname_suffix_len ! default fill value for history NetCDF fields real(r8), parameter, public :: hist_default_fillvalue = 1.e36_r8 - integer, parameter, public :: pfiles = 12 ! max number of tapes + integer, parameter, public :: pfiles = 12 ! max number of history configurations integer, parameter, public :: max_chars = shr_kind_cl ! max chars for char variables integer, parameter, public :: max_string_len = shr_kind_cxx - real(r8), parameter, public :: fillvalue = 1.e36_r8 ! fill value for netcdf fields + real(r8), parameter, public :: fillvalue = 1.e36_r8 ! default fill value for netcdf fields ! A special symbol for declaring a field which has no vertical or ! non-grid dimensions. It is here (rather than cam_history) so that it ! can be checked by add_hist_coord @@ -118,6 +119,7 @@ module cam_history_support public :: lookup_hist_coord_indices public :: hist_coord_find_levels public :: get_hist_coord_index + public :: get_hist_coord_names public :: parse_multiplier ! Parse a repeat count and a token from input interface add_hist_coord @@ -162,6 +164,22 @@ pure integer function get_hist_coord_index(mdimname) end function get_hist_coord_index + function get_hist_coord_names() result(mdimnames) + use cam_abortutils, only: endrun + character(len=max_hcoordname_len), allocatable :: mdimnames(:) + character(len=512) :: errmsg + integer :: ierr, idx + + allocate(mdimnames(registeredmdims), stat=ierr, errmsg=errmsg) + if (ierr /= 0) then + call endrun('get_hist_coord_names: failed to allocate mdimnames; errmsg = '//trim(errmsg)) + end if + do idx = 1, registeredmdims + mdimnames(idx) = hist_coords(idx)%name + end do + + end function get_hist_coord_names + ! Functions to check consistent term definition for hist coords pure logical function check_hist_coord_char(defined, input) diff --git a/src/utils/cam_grid_support.F90 b/src/utils/cam_grid_support.F90 index 1e4fc9666..478d23a7a 100644 --- a/src/utils/cam_grid_support.F90 +++ b/src/utils/cam_grid_support.F90 @@ -17,6 +17,7 @@ module cam_grid_support integer, parameter, public :: max_hcoordname_len = 16 integer, parameter, public :: max_split_files = 2 + integer, parameter :: max_files = max_split_files + 1 ! Split files plus restart type, public :: vardesc_ptr_t type(var_desc_t), pointer :: p => NULL() @@ -41,8 +42,8 @@ module cam_grid_support integer(iMap), pointer :: map(:) => NULL() ! map (dof) for dist. coord logical :: latitude ! .false. means longitude real(r8), pointer :: bnds(:,:) => NULL() ! bounds, if present - type(vardesc_ptr_t) :: vardesc(max_split_files) ! If we are to write coord - type(vardesc_ptr_t) :: bndsvdesc(max_split_files) ! Set to write bounds + type(vardesc_ptr_t) :: vardesc(max_files) ! If we are to write coord + type(vardesc_ptr_t) :: bndsvdesc(max_files) ! Set to write bounds contains procedure :: get_coord_len => horiz_coord_len procedure :: num_elem => horiz_coord_num_elem @@ -63,7 +64,7 @@ module cam_grid_support type, abstract :: cam_grid_attribute_t character(len=max_hcoordname_len) :: name = '' ! attribute name character(len=max_chars) :: long_name = '' ! attr long_name - type(vardesc_ptr_t) :: vardesc(max_split_files) + type(vardesc_ptr_t) :: vardesc(max_files) ! We aren't going to use this until we sort out PGI issues class(cam_grid_attribute_t), pointer :: next => NULL() contains @@ -165,7 +166,7 @@ module cam_grid_support type(horiz_coord_t), pointer :: lon_coord => NULL() ! Longitude logical :: unstructured ! Is this needed? logical :: block_indexed ! .false. for lon/lat - logical :: attrs_defined(max_split_files) = .false. + logical :: attrs_defined(max_files) = .false. logical :: zonal_grid = .false. type(cam_filemap_t), pointer :: map => null() ! global dim map (dof) type(cam_grid_attr_ptr_t), pointer :: attributes => NULL() @@ -191,6 +192,7 @@ module cam_grid_support procedure :: read_darray_3d_double => cam_grid_read_darray_3d_double procedure :: read_darray_2d_real => cam_grid_read_darray_2d_real procedure :: read_darray_3d_real => cam_grid_read_darray_3d_real + procedure :: write_darray_1d_int => cam_grid_write_darray_1d_int procedure :: write_darray_2d_int => cam_grid_write_darray_2d_int procedure :: write_darray_3d_int => cam_grid_write_darray_3d_int procedure :: write_darray_1d_double => cam_grid_write_darray_1d_double @@ -352,6 +354,7 @@ end subroutine print_attr_spec end interface cam_grid_read_dist_array interface cam_grid_write_dist_array + module procedure cam_grid_write_dist_array_1d_int module procedure cam_grid_write_dist_array_2d_int module procedure cam_grid_write_dist_array_3d_int module procedure cam_grid_write_dist_array_1d_double @@ -1322,6 +1325,41 @@ subroutine cam_grid_read_dist_array_3d_real(File, id, adims, fdims, & end subroutine cam_grid_read_dist_array_3d_real + !------------------------------------------------------------------------ + ! + ! cam_grid_write_dist_array_1d_int + ! + ! Interface function for the grid%write_darray_1d_int method + ! + !------------------------------------------------------------------------ + subroutine cam_grid_write_dist_array_1d_int(File, id, adims, fdims, & + hbuf, varid) + use pio, only: file_desc_t + + ! Dummy arguments + type(file_desc_t), intent(inout) :: File ! PIO file handle + integer, intent(in) :: id + integer, intent(in) :: adims(:) + integer, intent(in) :: fdims(:) + integer, intent(in) :: hbuf(:) + type(var_desc_t), intent(inout) :: varid + + ! Local variable + integer :: gridid + character(len=shr_kind_cm) :: errormsg + + gridid = get_cam_grid_index(id) + if (gridid > 0) then + call cam_grids(gridid)%write_darray_1d_int(File, adims, fdims, & + hbuf, varid) + else + write(errormsg, *) & + 'cam_grid_write_dist_array_1d_int: Bad grid ID, ', id + call endrun(errormsg) + end if + + end subroutine cam_grid_write_dist_array_1d_int + !------------------------------------------------------------------------ ! ! cam_grid_write_dist_array_2d_int @@ -3589,6 +3627,36 @@ subroutine cam_grid_read_darray_3d_real(this, File, adims, fdims, & call cam_pio_handle_error(ierr, subname//': Error reading variable') end subroutine cam_grid_read_darray_3d_real + !------------------------------------------------------------------------ + ! + ! cam_grid_write_darray_1d_int: Write a variable defined on this grid + ! + !------------------------------------------------------------------------ + subroutine cam_grid_write_darray_1d_int(this, File, adims, fdims, & + hbuf, varid) + use pio, only: file_desc_t, io_desc_t + use pio, only: pio_write_darray, PIO_INT + + use cam_pio_utils, only: cam_pio_get_decomp + + ! Dummy arguments + class(cam_grid_t) :: this + type(file_desc_t), intent(inout) :: File ! PIO file handle + integer, intent(in) :: adims(:) + integer, intent(in) :: fdims(:) + integer, intent(in) :: hbuf(:) + type(var_desc_t), intent(inout) :: varid + + ! Local variables + type(io_desc_t), pointer :: iodesc + integer :: ierr + character(len=*), parameter :: subname = 'cam_grid_write_darray_1d_int' + + call cam_pio_get_decomp(iodesc, adims, fdims, PIO_INT, this%map) + call pio_write_darray(File, varid, iodesc, hbuf, ierr) + call cam_pio_handle_error(ierr, subname//': Error writing variable') + end subroutine cam_grid_write_darray_1d_int + !------------------------------------------------------------------------ ! ! cam_grid_write_darray_2d_int: Write a variable defined on this grid