diff --git a/src/Utilities/Idm/DefinitionSelect.f90 b/src/Utilities/Idm/DefinitionSelect.f90 index 656312e9cf8..9393373de24 100644 --- a/src/Utilities/Idm/DefinitionSelect.f90 +++ b/src/Utilities/Idm/DefinitionSelect.f90 @@ -16,6 +16,7 @@ module DefinitionSelectModule implicit none private public :: get_param_definition_type + public :: package_scoped_param_dfn public :: get_aggregate_definition_type public :: split_record_definition @@ -63,6 +64,53 @@ function get_param_definition_type(input_definition_types, & return end function get_param_definition_type + !> @brief Return parameter definition without checking blockname + !< + function package_scoped_param_dfn(input_definition_types, & + component_type, subcomponent_type, & + tagname, filename) & + result(idt) + type(InputParamDefinitionType), dimension(:), intent(in), target :: & + input_definition_types + character(len=*), intent(in) :: component_type !< component type, such as GWF or GWT + character(len=*), intent(in) :: subcomponent_type !< subcomponent type, such as DIS or NPF + character(len=*), intent(in) :: tagname !< name of the input tag + character(len=*), intent(in) :: filename !< input filename + type(InputParamDefinitionType), pointer :: idt !< corresponding InputParameterDefinitionType for this tag + type(InputParamDefinitionType), pointer :: tmp_ptr + integer(I4B) :: i + ! + idt => null() + ! + do i = 1, size(input_definition_types) + tmp_ptr => input_definition_types(i) + if (tmp_ptr%component_type == component_type .and. & + tmp_ptr%subcomponent_type == subcomponent_type .and. & + tmp_ptr%tagname == tagname) then + if (associated(idt)) then + write (errmsg, '(a,a,a)') & + 'Input file tag name "', trim(tagname), & + '" is not unique (file scope) in package definition set.' + call store_error(errmsg) + call store_error_filename(filename) + else + idt => input_definition_types(i) + end if + end if + end do + ! + if (.not. associated(idt)) then + write (errmsg, '(a,a,a)') & + 'Input file tag not found: "', trim(tagname), & + '".' + call store_error(errmsg) + call store_error_filename(filename) + end if + ! + ! -- return + return + end function package_scoped_param_dfn + !> @brief Return aggregate definition !< function get_aggregate_definition_type(input_definition_types, component_type, & diff --git a/src/Utilities/Idm/IdmLoad.f90 b/src/Utilities/Idm/IdmLoad.f90 index fe6d5b263e6..1f61462990d 100644 --- a/src/Utilities/Idm/IdmLoad.f90 +++ b/src/Utilities/Idm/IdmLoad.f90 @@ -132,10 +132,13 @@ end subroutine idm_ad !> @brief idm deallocate routine !< subroutine idm_da(iout) + use InputModelContextModule, only: ModelContextDestroy integer(I4B), intent(in) :: iout ! call dynamic_da(iout) ! + call ModelContextDestroy() + ! ! -- return return end subroutine idm_da @@ -189,12 +192,19 @@ end subroutine model_pkg_load !< subroutine load_model_pkgs(model_pkg_inputs, iout) use ModelPackageInputsModule, only: ModelPackageInputsType - use SourceLoadModule, only: open_source_file + use SourceLoadModule, only: open_source_file, create_context use IdmDfnSelectorModule, only: idm_integrated type(ModelPackageInputsType), intent(inout) :: model_pkg_inputs integer(i4B), intent(in) :: iout integer(I4B) :: itype, ipkg ! + ! -- create package context object + call create_context(model_pkg_inputs%modeltype, & + model_pkg_inputs%component_type, & + model_pkg_inputs%modelname, & + model_pkg_inputs%modelfname, & + model_pkg_inputs%pkglist, iout) + ! ! -- load package instances by type do itype = 1, size(model_pkg_inputs%pkglist) ! diff --git a/src/Utilities/Idm/ModflowInput.f90 b/src/Utilities/Idm/ModflowInput.f90 index 9500fc5dcb6..984f5df08d5 100644 --- a/src/Utilities/Idm/ModflowInput.f90 +++ b/src/Utilities/Idm/ModflowInput.f90 @@ -10,7 +10,8 @@ module ModflowInputModule use KindModule, only: I4B, LGP use ConstantsModule, only: LENMEMPATH, LENCOMPONENTNAME, & - LENPACKAGETYPE + LENPACKAGETYPE, LENPACKAGENAME, & + LINELENGTH use MemoryHelperModule, only: create_mem_path use InputDefinitionModule, only: InputParamDefinitionType, & InputBlockDefinitionType @@ -23,9 +24,9 @@ module ModflowInputModule private public :: ModflowInputType, getModflowInput - !> @brief derived type for storing input definition for a file + !> @brief type for storing input definition for a file !! - !! This derived type contains the information needed to read + !! This type contains the information needed to read !! a specific modflow input file, including block definitions, !! aggregate definitions (structarrays), and individual !! parameter definitions. @@ -63,7 +64,8 @@ function getModflowInput(pkgtype, component_type, subcomponent_type, & ! -- set subcomponent type if (present(filename)) then dfn_subcomponent_type = update_sc_type(pkgtype, filename, component_type, & - subcomponent_type) + subcomponent_type, component_name, & + subcomponent_name) else dfn_subcomponent_type = trim(subcomponent_type) end if @@ -90,23 +92,45 @@ function getModflowInput(pkgtype, component_type, subcomponent_type, & mf6_input%subcomponent_type) end function getModflowInput - function update_sc_type(filetype, filename, component_type, subcomponent_type) & + function update_sc_type(filetype, filename, component_type, subcomponent_type, & + component_name, subcomponent_name) & result(sc_type) - character(len=*), intent(in) :: component_type - character(len=*), intent(in) :: subcomponent_type + use SourceCommonModule, only: package_source_type + use InputModelContextModule, only: GetModelNC4Context + use NC4ModelInputsModule, only: NC4ModelInputsType, NC4ModelPackageInputType character(len=*), intent(in) :: filetype character(len=*), intent(in) :: filename + character(len=*), intent(in) :: component_type + character(len=*), intent(in) :: subcomponent_type + character(len=*), intent(in) :: component_name + character(len=*), intent(in) :: subcomponent_name ! -- result character(len=LENPACKAGETYPE) :: sc_type + ! -- local + character(len=LENPACKAGENAME) :: source_type + type(NC4ModelInputsType), pointer :: nc4_context + type(NC4ModelPackageInputType), pointer :: ncpkg ! sc_type = subcomponent_type + source_type = package_source_type(filename) ! - select case (subcomponent_type) - case ('RCH', 'EVT', 'SCP') - sc_type = read_as_arrays(filetype, filename, component_type, & - subcomponent_type) - case default - end select + if (source_type == 'MF6FILE') then + select case (subcomponent_type) + case ('RCH', 'EVT', 'SCP') + sc_type = read_as_arrays(filetype, filename, component_type, & + subcomponent_type) + case default + end select + else if (source_type == 'NETCDF4') then + select case (subcomponent_type) + case ('RCH', 'EVT', 'SCP') + nc4_context => GetModelNC4Context(component_name) + ncpkg => nc4_context%get_package(component_name, & + subcomponent_name) + sc_type = trim(ncpkg%subcomponent_type) + case default + end select + end if ! ! -- return return diff --git a/src/Utilities/Idm/SourceContext.f90 b/src/Utilities/Idm/SourceContext.f90 new file mode 100644 index 00000000000..b833f31e6ec --- /dev/null +++ b/src/Utilities/Idm/SourceContext.f90 @@ -0,0 +1,222 @@ +!> @brief This module contains the InputModelContextodule +!! +!! +!< +module InputModelContextModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: LINELENGTH, LENMODELNAME, LENTIMESERIESNAME + use SimModule, only: store_error, store_error_filename + use SimVariablesModule, only: errmsg + use ListModule, only: ListType + use InputDefinitionModule, only: InputParamDefinitionType + use NC4ModelInputsModule, only: NC4ModelInputsType + + implicit none + private + public :: ModelContextType + public :: AddModelNC4Context + public :: GetModelNC4Context + public :: ModelContextDestroy + + type(ListType) :: model_context_list + + !> @brief type for storing model context + !! + !! This type is used to store a list of context objects + !! associated with a model. Add additional context types as + !! appropriate. + !! + !< + type :: ModelContextType + character(len=LENMODELNAME) :: modelname !< name of model + character(len=LINELENGTH) :: modelfname !< name of model input file + type(NC4ModelInputsType), pointer :: nc4_context + contains + procedure :: init => modelctx_init + procedure :: destroy => modelctx_destroy + end type ModelContextType + +contains + + !> @brief model context init + !! + !< + subroutine modelctx_init(this, modelname, modelfname) + class(ModelContextType), intent(inout) :: this + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + ! + this%modelname = modelname + this%modelfname = modelfname + ! + nullify (this%nc4_context) + ! + return + end subroutine modelctx_init + + !> @brief destroy model context object + !! + !< + subroutine modelctx_destroy(this) + class(ModelContextType), intent(inout) :: this + ! + if (associated(this%nc4_context)) then + call this%nc4_context%destroy() + end if + ! + return + end subroutine modelctx_destroy + + !> @brief add model context object to list + !! + !< + subroutine AddModelContext(model_context) + ! -- dummy variables + class(ModelContextType), pointer, intent(inout) :: model_context + ! -- local variables + class(*), pointer :: obj + ! + obj => model_context + call model_context_list%Add(obj) + ! + ! -- return + return + end subroutine AddModelContext + + !> @brief get model context object from list + !! + !< + function GetModelContext(modelname) result(res) + ! -- dummy variables + character(len=*), intent(in) :: modelname + class(ModelContextType), pointer :: res + ! -- local variables + class(*), pointer :: obj + class(ModelContextType), pointer :: ctx + integer(I4B) :: n + ! + ! -- initialize res + res => null() + ! + ! -- get the object from the list + do n = 1, model_context_list%Count() + obj => model_context_list%GetItem(n) + if (associated(obj)) then + select type (obj) + class is (ModelContextType) + ctx => obj + if (ctx%modelname == modelname) then + res => obj + exit + end if + end select + end if + end do + ! + ! -- return + return + end function GetModelContext + + !> @brief cleanup model context objects + !! + !< + subroutine ModelContextDestroy() + ! -- dummy variables + ! -- local variables + class(*), pointer :: obj + class(ModelContextType), pointer :: ctx + integer(I4B) :: n + ! + ! -- get the object from the list + do n = 1, model_context_list%Count() + obj => model_context_list%GetItem(n) + if (associated(obj)) then + select type (obj) + class is (ModelContextType) + ctx => obj + call ctx%destroy() + deallocate (ctx) + nullify (ctx) + end select + end if + end do + ! + call model_context_list%clear() + ! + ! -- return + return + end subroutine ModelContextDestroy + + !> @brief get model context object from list + !! + !< + subroutine AddModelNC4Context(modelname, modelfname, nc4_context) + ! -- dummy variables + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + type(NC4ModelInputsType), pointer, intent(in) :: nc4_context + ! -- local variables + class(ModelContextType), pointer :: ctx + ! + ! -- initialize + nullify (ctx) + ! + ctx => GetModelContext(modelname) + ! + if (associated(ctx)) then + ctx%nc4_context => nc4_context + else + allocate (ctx) + call ctx%init(modelname, modelfname) + ctx%nc4_context => nc4_context + call AddModelContext(ctx) + end if + ! + ! -- return + return + end subroutine AddModelNC4Context + + !> @brief get model context object from list + !! + !< + function GetModelNC4Context(modelname) result(nc4_context) + ! -- dummy variables + character(len=*), intent(in) :: modelname + ! -- result + type(NC4ModelInputsType), pointer :: nc4_context + ! -- local variables + class(*), pointer :: obj + class(ModelContextType), pointer :: ctx + integer(I4B) :: n + ! + ! -- initialize res + nc4_context => null() + ! + ! -- get the object from the list + do n = 1, model_context_list%Count() + obj => model_context_list%GetItem(n) + if (associated(obj)) then + select type (obj) + class is (ModelContextType) + ctx => obj + if (ctx%modelname == modelname) then + nc4_context => ctx%nc4_context + exit + end if + end select + end if + end do + ! + ! -- set error if not found + if (.not. associated(nc4_context)) then + errmsg = 'Programming error. NC4 Model context not found. & + &Model='//trim(modelname)//'.' + call store_error(errmsg, .true.) + end if + ! + ! -- return + return + end function GetModelNC4Context + +end module InputModelContextModule diff --git a/src/Utilities/Idm/SourceLoad.F90 b/src/Utilities/Idm/SourceLoad.F90 index 236607806c7..e6c60898d58 100644 --- a/src/Utilities/Idm/SourceLoad.F90 +++ b/src/Utilities/Idm/SourceLoad.F90 @@ -19,6 +19,7 @@ module SourceLoadModule public :: create_pkg_loader public :: open_source_file public :: load_modelnam, load_simnam + public :: create_context contains @@ -52,7 +53,7 @@ function create_pkg_loader(component_type, subcomponent_type, pkgname, & source_type = package_source_type(filename) ! ! -- set source loader for model package - loader => package_loader(source_type) + loader => package_loader(source_type, modelfname) ! ! -- initialize loader call loader%init(mf6_input, modelname, modelfname, filename) @@ -63,11 +64,18 @@ end function create_pkg_loader !> @brief allocate source model package static loader !< - function package_loader(source_type) result(loader) + function package_loader(source_type, modelfname) result(loader) use InputLoadTypeModule, only: StaticPkgLoadBaseType use IdmMf6FileModule, only: Mf6FileStaticPkgLoadType +#if defined(__WITH_NETCDF__) + use IdmNetCDFFileModule, only: NC4StaticPkgLoadType +#endif character(len=*), intent(inout) :: source_type + character(len=*), intent(in) :: modelfname class(Mf6FileStaticPkgLoadType), pointer :: mf6file_loader +#if defined(__WITH_NETCDF__) + class(NC4StaticPkgLoadType), pointer :: nc4_loader +#endif class(StaticPkgLoadBaseType), pointer :: loader ! ! -- initialize @@ -78,11 +86,23 @@ function package_loader(source_type) result(loader) case ('MF6FILE') allocate (mf6file_loader) loader => mf6file_loader + case ('NETCDF4') +#if defined(__WITH_NETCDF__) + allocate (nc4_loader) + loader => nc4_loader +#else + write (errmsg, '(a)') & + 'Cannot load package inputs. NetCDF4 input file provided & + &but NetCDF4 libraries not available.' + call store_error(errmsg) + call store_error_filename(modelfname) +#endif case default write (errmsg, '(a)') & 'Simulation package input source type "'//trim(source_type)// & '" not currently supported.' - call store_error(errmsg, .true.) + call store_error(errmsg) + call store_error_filename(modelfname) end select ! ! -- return @@ -108,6 +128,16 @@ function open_source_file(pkgtype, filename, modelfname, iout) result(fd) select case (source_type) case ('MF6FILE') fd = open_mf6file(pkgtype, filename, modelfname, iout) + case ('NETCDF4') +#if defined(__WITH_NETCDF__) + ! -- no-op, don't provide ncid fd to a package +#else + write (errmsg, '(a)') & + 'Cannot open source input file. NetCDF4 input file provided & + &but NetCDF4 libraries not available.' + call store_error(errmsg) + call store_error_filename(modelfname) +#endif case default end select ! @@ -130,12 +160,25 @@ subroutine load_modelnam(mtype, mfname, mname, iout) source_type = package_source_type(mfname) ! ! -- create description of input - mf6_input = getModflowInput(mtype, idm_component_type(mtype), 'NAM', & - mname, 'NAM', mfname) + mf6_input = getModflowInput(mtype, idm_component_type(mtype), & + 'NAM', mname, 'NAM') ! select case (source_type) case ('MF6FILE') call input_load(mfname, mf6_input, simfile, iout) + case ('NETCDF4') +#if defined(__WITH_NETCDF__) + write (errmsg, '(a)') & + 'NetCDF4 Model name files not currently supported.' + call store_error(errmsg) + call store_error_filename(simfile) +#else + write (errmsg, '(a)') & + 'Cannot load name file. NetCDF4 input file provided & + &but but NetCDF4 libraries not available.' + call store_error(errmsg) + call store_error_filename(simfile) +#endif case default end select ! @@ -162,7 +205,7 @@ subroutine load_simnam() call sim_message(line, skipafter=1) ! ! -- create description of input - mf6_input = getModflowInput('NAM6', 'SIM', 'NAM', 'SIM', 'NAM', simfile) + mf6_input = getModflowInput('NAM6', 'SIM', 'NAM', 'SIM', 'NAM') ! ! -- open namfile and load to input context call input_load(simfile, mf6_input, simfile, iout) @@ -172,4 +215,182 @@ subroutine load_simnam() return end subroutine load_simnam + !> @brief create model context + !< + subroutine create_context(modeltype, component_type, modelname, & + modelfname, pkglist, iout) + ! -- modules + use ModelPackageInputsModule, only: LoadablePackageType + ! -- drummy + character(len=*), intent(in) :: modeltype + character(len=*), intent(in) :: component_type + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + type(LoadablePackageType), dimension(:), & + allocatable, intent(in) :: pkglist + integer(I4B), intent(in) :: iout + ! -- local + ! + ! -- create netcdf4 context + call create_nc4_context(modeltype, component_type, modelname, & + modelfname, pkglist, iout) + ! + ! -- return + return + end subroutine create_context + + !> @brief create model netcdf4 context + !< + subroutine create_nc4_context(modeltype, component_type, modelname, & + modelfname, pkglist, iout) + use ModelPackageInputsModule, only: LoadablePackageType + use NC4ModelInputsModule, only: NC4ModelInputsType + use InputModelContextModule, only: AddModelNC4Context +#if defined(__WITH_NETCDF__) + use LoadNetCDFDataModule, only: nc4_pkg_context + use IdmNetCDFFileModule, only: open_ncfile +#endif + ! -- drummy + character(len=*), intent(in) :: modeltype + character(len=*), intent(in) :: component_type + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + type(LoadablePackageType), dimension(:), & + allocatable, intent(in) :: pkglist + integer(I4B), intent(in) :: iout + ! -- local + character(len=LINELENGTH) :: nc4_fname + type(NC4ModelInputsType), pointer :: nc4_context +#if defined(__WITH_NETCDF__) + integer(I4B) :: ncid +#endif + ! + ! -- allocate context object + allocate (nc4_context) + ! + ! -- set model nc filename if provided + nc4_fname = nc4_filename(pkglist, modelfname) + ! + if (nc4_fname /= '') then +#if defined(__WITH_NETCDF__) + ! + ! -- open nc4 input file + ncid = open_ncfile(nc4_fname, iout) + ! + ! -- init model context object + call nc4_context%init(modeltype, component_type, modelname, nc4_fname, ncid) + ! + ! -- add NETCDF4 packages to context + call set_nc4_pkglist(component_type, modelfname, pkglist, & + nc4_context, nc4_fname) + ! + ! -- read the file and build the context + call nc4_pkg_context(nc4_context, iout) + ! +#else + write (errmsg, '(a)') & + 'Cannot load model packages. NetCDF4 input file & + &specified in model namefile options block but & + &NetCDF4 libraries are not available.' + call store_error(errmsg) + call store_error_filename(modelfname) +#endif + else + ! + ! -- initialize object + call nc4_context%init(modeltype, component_type, modelname, '', 0) + ! + end if + ! + ! -- add context to model context list + call AddModelNC4Context(modelname, modelfname, nc4_context) + ! + ! -- return + return + end subroutine create_nc4_context + +#if defined(__WITH_NETCDF__) + subroutine set_nc4_pkglist(component_type, modelfname, pkglist, & + nc4_context, nc4_fname) + use SourceCommonModule, only: package_source_type + use IdmDfnSelectorModule, only: idm_integrated + use ModelPackageInputsModule, only: LoadablePackageType + use NC4ModelInputsModule, only: NC4ModelInputsType + ! -- drummy + character(len=*), intent(in) :: component_type + character(len=*), intent(in) :: modelfname + type(LoadablePackageType), dimension(:), & + allocatable, intent(in) :: pkglist + type(NC4ModelInputsType), pointer, & + intent(inout) :: nc4_context + character(len=*) :: nc4_fname + ! -- local + integer(I4B) :: m, n + ! + ! -- add NETCDF4 packages to context list + do n = 1, size(pkglist) + ! -- load package instances + do m = 1, pkglist(n)%pnum + ! + if (package_source_type(pkglist(n)%filenames(m)) == 'NETCDF4') then + ! + if (idm_integrated(component_type, pkglist(n)%subcomponent_type)) then + ! + ! -- add pkg instance to context object + call nc4_context%add(pkglist(n)%pkgtype, & + pkglist(n)%pkgnames(m), & + pkglist(n)%subcomponent_type) + ! + else + errmsg = 'NetCDF4 is unsupported for package type "'// & + trim(pkglist(n)%pkgtype)//'".' + call store_error(errmsg) + call store_error_filename(modelfname) + end if + end if + end do + end do + ! + ! -- return + return + end subroutine set_nc4_pkglist +#endif + + !> @brief set model netcdf4 input filename + !< + function nc4_filename(pkglist, modelfname) result(nc4_fname) + use ModelPackageInputsModule, only: LoadablePackageType + use SourceCommonModule, only: package_source_type + ! -- drummy + type(LoadablePackageType), dimension(:), & + allocatable, intent(in) :: pkglist + character(len=*), intent(in) :: modelfname + ! -- local + character(len=LINELENGTH) :: nc4_fname + integer(I4B) :: m, n + ! + ! -- initialize + nc4_fname = '' + ! + ! -- verify single model nc file + do n = 1, size(pkglist) + do m = 1, pkglist(n)%pnum + if (package_source_type(pkglist(n)%filenames(m)) == 'NETCDF4') then + if (nc4_fname == '') then + nc4_fname = pkglist(n)%filenames(m) + else if (nc4_fname /= pkglist(n)%filenames(m)) then + nc4_fname = '' + errmsg = 'Multiple *.nc model input files detected in packages & + &block. Only one model NetCDF4 input supported.' + call store_error(errmsg) + call store_error_filename(modelfname) + end if + end if + end do + end do + ! + ! -- return + return + end function nc4_filename + end module SourceLoadModule diff --git a/src/Utilities/Idm/netcdf/IdmNetCDFFile.f90 b/src/Utilities/Idm/netcdf/IdmNetCDFFile.f90 new file mode 100644 index 00000000000..f609d962a8b --- /dev/null +++ b/src/Utilities/Idm/netcdf/IdmNetCDFFile.f90 @@ -0,0 +1,588 @@ +!> @brief This module contains the IdmNetCDFFileModule +!! +!! This module contains the high-level routines for loading +!! a MODFLOW netcdf input file to the input context. +!! +!< +module IdmNetCDFFileModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: LINELENGTH, LENCOMPONENTNAME, LENVARNAME + use SimModule, only: store_error, store_error_filename + use SimVariablesModule, only: errmsg + use ModflowInputModule, only: ModflowInputType + use InputLoadTypeModule, only: StaticPkgLoadBaseType, DynamicPkgLoadBaseType + use NC4ModelInputsModule, only: NC4BlockVariablesType + use NC4ModelInputsModule, only: NC4ModelInputsType, NC4ModelPackageInputType + use BoundInputContextModule, only: BoundInputContextType + use ArrayHandlersModule, only: expandarray + use InputLoadTypeModule, only: DynamicPkgLoadType + use SourceCommonModule, only: ReadStateVarType + + implicit none + private + public :: open_ncfile + public :: NC4StaticPkgLoadType + + !> @brief base abstract type for netcdf source dynamic load + !! + !< + type, abstract, extends(DynamicPkgLoadType) :: NCDynamicPkgLoadBaseType + contains + procedure(nc_period_load_if), deferred :: rp + end type NCDynamicPkgLoadBaseType + + abstract interface + subroutine nc_period_load_if(this, ncid, ncpkg) + import NCDynamicPkgLoadBaseType, NC4ModelPackageInputType, I4B + class(NCDynamicPkgLoadBaseType), intent(inout) :: this + integer(I4B), intent(in) :: ncid + type(NC4ModelPackageInputType), pointer, intent(inout) :: ncpkg + end subroutine + end interface + + type, extends(StaticPkgLoadBaseType) :: NC4StaticPkgLoadType + integer(I4B) :: ncid + type(NC4ModelPackageInputType), pointer :: ncpkg => null() + contains + procedure :: init => static_init + procedure :: load => static_load + end type NC4StaticPkgLoadType + + type, extends(DynamicPkgLoadBaseType) :: NC4DynamicPkgLoadType + integer(I4B), pointer :: iper => null() + integer(I4B), pointer :: ionper => null() + integer(I4B) :: iiper + integer(I4B) :: ncid + type(NC4ModelPackageInputType), pointer :: ncpkg => null() + class(NCDynamicPkgLoadBaseType), pointer :: loader => null() + contains + procedure :: init => dynamic_init + procedure :: set => dynamic_set + procedure :: df => dynamic_df + procedure :: rp => dynamic_load + procedure :: create_loader => dynamic_create_loader + procedure :: destroy => dynamic_destroy + end type NC4DynamicPkgLoadType + + !> @brief NetCDF list based dynamic loader type + !< + type, extends(NCDynamicPkgLoadBaseType) :: NCListInputType + integer(I4B) :: ncid + integer(I4B) :: iiper + integer(I4B) :: nparam + character(len=LENVARNAME), dimension(:), & + allocatable :: param_names + type(BoundInputContextType) :: bndctx + contains + procedure :: init => inlist_init + procedure :: rp => inlist_rp + procedure :: validate => inlist_validate + procedure :: destroy => inlist_destroy + end type NCListInputType + + !> @brief NetCDF array based dynamic loader type + !< + type, extends(NCDynamicPkgLoadBaseType) :: NCGridInputType + integer(I4B) :: ncid + integer(I4B) :: iiper + integer(I4B) :: nparam !< number of in scope dynamic parameters + character(len=LENVARNAME), dimension(:), & + allocatable :: param_names !< dynamic param names + type(ReadStateVarType), dimension(:), allocatable :: param_reads !< read states for current load + type(BoundInputContextType) :: bndctx !< boundary package input context + contains + procedure :: init => ingrid_init + procedure :: alloc => ingrid_alloc + procedure :: rp => ingrid_rp + procedure :: validate => ingrid_validate + procedure :: destroy => ingrid_destroy + end type NCGridInputType + +contains + + subroutine static_init(this, mf6_input, modelname, modelfname, source) + use InputModelContextModule, only: GetModelNC4Context + class(NC4StaticPkgLoadType), intent(inout) :: this + type(ModflowInputType), intent(in) :: mf6_input + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + character(len=*), intent(in) :: source + type(NC4ModelInputsType), pointer :: nc4_context + ! + call this%StaticPkgLoadType%init(mf6_input, modelname, modelfname, source) + ! + ! -- set context from input file + nc4_context => GetModelNC4Context(modelname) + ! + ! -- cache ncid for model netcdf4 file + this%ncid = nc4_context%ncid + ! + ! -- set package description from context + this%ncpkg => & + nc4_context%get_package(mf6_input%component_name, & + mf6_input%subcomponent_name) + ! -- return + return + end subroutine static_init + + function static_load(this, iout) result(period_loader) + use LoadNetCDFDataModule, only: input_load + class(NC4StaticPkgLoadType), intent(inout) :: this + integer(I4B), intent(in) :: iout + class(DynamicPkgLoadBaseType), pointer :: period_loader + class(NC4DynamicPkgLoadType), pointer :: nc4_loader => null() + ! + ! -- initialize + nullify (period_loader) + ! + ! -- load model package to input context + call input_load(this%mf6_input, this%ncpkg, this%ncid, this%sourcename, iout) + ! + ! -- check if package is dynamic + if (this%iperblock > 0) then + ! + ! -- create dynamic loader + allocate (nc4_loader) + ! + ! -- initailze loader + call nc4_loader%init(this%mf6_input, this%modelname, & + this%modelfname, this%sourcename, & + this%iperblock, iout) + ! + ! -- set period data + call nc4_loader%set(this%ncpkg, this%ncid) + ! + ! -- set returned base pointer + period_loader => nc4_loader + ! + end if + ! + ! -- return + return + end function static_load + + subroutine dynamic_init(this, mf6_input, modelname, modelfname, & + source, iperblock, iout) + use MemoryManagerModule, only: mem_allocate + class(NC4DynamicPkgLoadType), intent(inout) :: this + type(ModflowInputType), intent(in) :: mf6_input + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + character(len=*), intent(in) :: source + integer(I4B), intent(in) :: iperblock + integer(I4B), intent(in) :: iout + ! + call this%DynamicPkgLoadType%init(mf6_input, modelname, modelfname, & + source, iperblock, iout) + ! + ! -- allocate iper and ionper + call mem_allocate(this%iper, 'IPER', this%mf6_input%mempath) + call mem_allocate(this%ionper, 'IONPER', this%mf6_input%mempath) + ! + ! -- initialize package + this%iper = 0 + this%ionper = 0 + this%iiper = 0 + ! + ! -- return + return + end subroutine dynamic_init + + subroutine dynamic_set(this, ncpkg, ncid) + use MemoryManagerModule, only: mem_allocate + class(NC4DynamicPkgLoadType), intent(inout) :: this + type(NC4ModelPackageInputType), pointer, intent(in) :: ncpkg + integer(I4B), intent(in) :: ncid + ! + ! -- set context + this%ncpkg => ncpkg + this%ncid = ncid + ! + ! -- verify block variable iper list + if (.not. allocated(ncpkg%ipers)) then + errmsg = 'Required NetCDF package variable IPERS not found for package "' & + //trim(this%mf6_input%subcomponent_name)//'".' + call store_error(errmsg) + call store_error_filename(this%sourcename) + end if + ! + ! -- create the loader + call this%create_loader() + ! + ! -- return + return + end subroutine dynamic_set + + subroutine dynamic_df(this) + use TdisModule, only: nper + class(NC4DynamicPkgLoadType), intent(inout) :: this + ! + ! -- set first ionper + if (size(this%ncpkg%ipers) > this%iiper) then + this%iiper = this%iiper + 1 + this%ionper = this%ncpkg%ipers(this%iiper) + else + this%ionper = nper + 1 + end if + ! + ! -- return + return + end subroutine dynamic_df + + subroutine dynamic_load(this) + use TdisModule, only: kper, nper + class(NC4DynamicPkgLoadType), intent(inout) :: this + ! + ! -- check if ready to load + if (this%ionper /= kper) return + ! + ! -- package dynamic load + call this%loader%rp(this%ncid, this%ncpkg) + ! + ! -- update loaded iper + this%iper = kper + ! + ! -- read next ionper + if (size(this%ncpkg%ipers) > this%iiper) then + this%iiper = this%iiper + 1 + this%ionper = this%ncpkg%ipers(this%iiper) + else + this%ionper = nper + 1 + end if + ! + ! -- return + return + end subroutine dynamic_load + + !> @brief allocate a dynamic loader based on load context + !< + subroutine dynamic_create_loader(this) + ! -- dummy + class(NC4DynamicPkgLoadType), intent(inout) :: this + class(NCListInputType), pointer :: list_loader + class(NCGridInputType), pointer :: grid_loader + ! + ! -- allocate and set loader + if (this%readasarrays) then + allocate (grid_loader) + this%loader => grid_loader + else + allocate (list_loader) + this%loader => list_loader + end if + ! + ! -- initialize loader + call this%loader%init(this%mf6_input, & + this%modelname, & + this%modelfname, & + this%sourcename, & + this%iperblock, & + this%iout) + ! + ! -- validate package netcdf inputs + if (this%readasarrays) then + call grid_loader%validate(this%ncpkg) + else + call list_loader%validate(this%ncpkg) + end if + ! + ! -- return + return + end subroutine dynamic_create_loader + + subroutine dynamic_destroy(this) + class(NC4DynamicPkgLoadType), intent(inout) :: this + ! + ! -- destroy and deallocate loader + call this%loader%destroy() + deallocate (this%loader) + nullify (this%loader) + ! + ! -- deallocate input context + call this%DynamicPkgLoadType%destroy() + ! + ! -- return + return + end subroutine dynamic_destroy + + !> @brief build the type list with all model package instances + !< + function open_ncfile(nc_fname, iout) result(ncid) + ! -- modules + use MemoryManagerExtModule, only: mem_set_value + use LoadNetCDFDataModule, only: nc_fopen + ! -- dummy + character(len=*) :: nc_fname + integer(I4B) :: iout + ! -- result + integer(I4B) :: ncid + ! -- local + logical(LGP) :: exists + ! + ! -- initialize + ncid = 0 + ! + ! -- check if NETCDF file exists + inquire (file=nc_fname, exist=exists) + if (.not. exists) then + write (errmsg, '(a,a,a)') 'Specified NetCDF4 input file does & + ¬ exist [file=', trim(nc_fname), '].' + call store_error(errmsg, .true.) + end if + ! + ! -- open + ncid = nc_fopen(nc_fname, iout) + ! + ! -- return + return + end function open_ncfile + + subroutine inlist_init(this, mf6_input, modelname, modelfname, & + source, iperblock, iout) + ! -- modules + ! -- dummy + class(NCListInputType), intent(inout) :: this + type(ModflowInputType), intent(in) :: mf6_input + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + character(len=*), intent(in) :: source + integer(I4B), intent(in) :: iperblock + integer(I4B), intent(in) :: iout + ! + ! --initialize + call this%DynamicPkgLoadType%init(mf6_input, modelname, modelfname, & + source, iperblock, iout) + this%iiper = 0 + ! + ! -- initialize context + call this%bndctx%init(this%mf6_input, .false.) + ! + ! -- allocate bound params + call this%bndctx%alloc(this%sourcename) + call this%bndctx%filtered_params(this%param_names, this%nparam) + ! + ! -- return + return + end subroutine inlist_init + + subroutine inlist_rp(this, ncid, ncpkg) + ! -- modules + use LoadNetCDFDataModule, only: nc4_rp_list + ! -- dummy + class(NCListInputType), intent(inout) :: this + integer(I4B), intent(in) :: ncid + type(NC4ModelPackageInputType), pointer, intent(inout) :: ncpkg + ! + ! -- increment iiper + this%iiper = this%iiper + 1 + ! + ! -- dynamic load + call nc4_rp_list(this%mf6_input, ncpkg, this%bndctx, ncid, & + this%iiper, this%sourcename, this%iout) + ! + ! -- return + return + end subroutine inlist_rp + + subroutine inlist_validate(this, ncpkg) + ! -- modules + use MemoryManagerModule, only: mem_setptr + use InputDefinitionModule, only: InputParamDefinitionType + use DefinitionSelectModule, only: get_param_definition_type + ! -- dummy + class(NCListInputType), intent(inout) :: this + type(NC4ModelPackageInputType), pointer, intent(inout) :: ncpkg + ! -- local + type(InputParamDefinitionType), pointer :: idt + integer(I4B) :: iparam, ivar + character(len=LINELENGTH) :: tagname, varname + logical(LGP) :: found + ! + ! -- verify required package variables are defined + do iparam = 1, this%nparam + tagname = this%param_names(iparam) + ! + idt => & + get_param_definition_type(this%mf6_input%param_dfns, & + this%mf6_input%component_type, & + this%mf6_input%subcomponent_type, & + 'PERIOD', tagname, '') + ! + if (idt%required) then + found = .false. + do ivar = 1, ncpkg%blocklist(this%iperblock)%varnum + varname = ncpkg%blocklist(this%iperblock)%varnames(ivar) + if (varname == tagname) then + found = .true. + exit + end if + end do + ! + if (.not. found) then + errmsg = 'Required PERIOD package variable not found => '//trim(tagname) + call store_error(errmsg) + call store_error_filename(this%sourcename) + end if + end if + end do + ! + ! -- return + return + end subroutine inlist_validate + + subroutine inlist_destroy(this) + ! -- modules + ! -- dummy + class(NCListInputType), intent(inout) :: this + if (allocated(this%param_names)) deallocate (this%param_names) + call this%bndctx%destroy() + end subroutine inlist_destroy + + subroutine ingrid_init(this, mf6_input, modelname, modelfname, & + source, iperblock, iout) + ! -- modules + ! -- dummy + class(NCGridInputType), intent(inout) :: this + type(ModflowInputType), intent(in) :: mf6_input + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: modelfname + character(len=*), intent(in) :: source + integer(I4B), intent(in) :: iperblock + integer(I4B), intent(in) :: iout + ! + ! -- initialize context + call this%DynamicPkgLoadType%init(mf6_input, modelname, modelfname, & + source, iperblock, iout) + ! + this%iiper = 0 + ! + call this%bndctx%init(this%mf6_input, .true.) + ! + ! -- allocate bound params + call this%alloc() + ! + ! -- return + return + end subroutine ingrid_init + + subroutine ingrid_alloc(this) + ! -- modules + use MemoryManagerModule, only: mem_setptr + ! -- dummy + class(NCGridInputType), intent(inout) :: this !< StressGridInputType + character(len=LENVARNAME) :: rs_varname + integer(I4B), pointer :: intvar + integer(I4B) :: iparam + ! + ! -- allocate period dfn params + call this%bndctx%alloc(this%sourcename) + ! + ! -- set in scope param names + call this%bndctx%filtered_params(this%param_names, this%nparam) + ! + ! -- allocate and set param_reads pointer array + allocate (this%param_reads(this%nparam)) + ! + ! store read state variable pointers + do iparam = 1, this%nparam + ! -- allocate and store name of read state variable + rs_varname = this%bndctx%rsv_alloc(this%param_names(iparam)) + call mem_setptr(intvar, rs_varname, this%mf6_input%mempath) + this%param_reads(iparam)%invar => intvar + this%param_reads(iparam)%invar = 0 + end do + ! + ! -- return + return + end subroutine ingrid_alloc + + subroutine ingrid_rp(this, ncid, ncpkg) + ! -- modules + use LoadNetCDFDataModule, only: nc4_rp_array + ! -- dummy + class(NCGridInputType), intent(inout) :: this + integer(I4B), intent(in) :: ncid + type(NC4ModelPackageInputType), pointer, intent(inout) :: ncpkg + integer(I4B) :: iparam + ! + ! -- increment iiper + this%iiper = this%iiper + 1 + ! + ! -- reset read state vars + do iparam = 1, this%nparam + this%param_reads(iparam)%invar = 0 + end do + ! + ! -- dynamic load + call nc4_rp_array(this%mf6_input, ncpkg, this%param_names, & + this%param_reads, this%bndctx, & + ncid, this%sourcename, this%iout) + ! + ! -- return + return + end subroutine ingrid_rp + + subroutine ingrid_validate(this, ncpkg) + ! -- modules + use InputDefinitionModule, only: InputParamDefinitionType + use DefinitionSelectModule, only: get_param_definition_type + ! -- dummy + class(NCGridInputType), intent(inout) :: this + type(NC4ModelPackageInputType), pointer, intent(inout) :: ncpkg + ! -- local + type(InputParamDefinitionType), pointer :: idt + character(len=LINELENGTH) :: tagname, varname + integer(I4B) :: iparam, ivar + logical(LGP) :: found + ! + do iparam = 1, this%nparam + tagname = this%param_names(iparam) + ! + idt => & + get_param_definition_type(this%mf6_input%param_dfns, & + this%mf6_input%component_type, & + this%mf6_input%subcomponent_type, & + 'PERIOD', tagname, '') + ! + if (idt%required) then + found = .false. + do ivar = 1, ncpkg%blocklist(this%iperblock)%varnum + varname = ncpkg%blocklist(this%iperblock)%varnames(ivar) + if (varname == tagname) then + found = .true. + exit + end if + end do + ! + if (.not. found) then + errmsg = 'Required PERIOD variable not found: "'//trim(tagname)// & + '" for package "'//trim(this%mf6_input%subcomponent_name)//'".' + call store_error(errmsg) + call store_error_filename(this%sourcename) + end if + end if + ! + end do + ! + ! -- return + return + end subroutine ingrid_validate + + subroutine ingrid_destroy(this) + ! -- modules + ! -- dummy + class(NCGridInputType), intent(inout) :: this + ! + ! -- deallocate allocatables + if (allocated(this%param_names)) deallocate (this%param_names) + if (allocated(this%param_reads)) deallocate (this%param_reads) + ! + ! -- destroy bnd context object + call this%bndctx%destroy() + ! + ! -- return + return + end subroutine ingrid_destroy + +end module IdmNetCDFFileModule diff --git a/src/Utilities/Idm/netcdf/LoadNetCDFData.f90 b/src/Utilities/Idm/netcdf/LoadNetCDFData.f90 new file mode 100644 index 00000000000..a1e39d45ab5 --- /dev/null +++ b/src/Utilities/Idm/netcdf/LoadNetCDFData.f90 @@ -0,0 +1,1933 @@ +module LoadNetCDFDataModule + + use netcdf + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: LINELENGTH, LENMEMPATH, LENVARNAME, & + LENAUXNAME, LENBOUNDNAME, LENMODELNAME, & + LENPACKAGENAME, LENCOMPONENTNAME, DNODATA, & + INODATA, LENPACKAGETYPE, LENMEMADDRESS + use SimVariablesModule, only: errmsg + use SimModule, only: store_error, store_error_filename + use ArrayHandlersModule, only: expandarray + use CharacterStringModule, only: CharacterStringType + use InputDefinitionModule, only: InputParamDefinitionType + use DefinitionSelectModule, only: get_param_definition_type, & + package_scoped_param_dfn, & + get_aggregate_definition_type + use ModflowInputModule, only: ModflowInputType + use MemoryManagerModule, only: mem_allocate, mem_reallocate, mem_setptr + use MemoryHelperModule, only: create_mem_path + use SimVariablesModule, only: idm_context + use IdmLoggerModule, only: idm_log_var, idm_log_header, idm_log_close + use NC4ModelInputsModule, only: NC4BlockVariablesType, & + NC4ModelInputsType, & + NC4ModelPackageInputType + + implicit none + private + public nc_fopen, nc_fclose + public input_load, nc4_rp_list, nc4_rp_array + public nc4_pkg_context + + integer(I4B), parameter :: IDM_NETCDF4_MAX_DIM = 6 + +contains + + function nc_fopen(nc_fname, iout) result(ncid) + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- return + integer(I4B) :: ncid + ! -- local + ! + ! -- initialize + ncid = -1 + ! + ! -- open netcdf file + call nf_verify(nf90_open(nc_fname, NF90_NOWRITE, ncid), ncid, iout) + write (iout, '(a,i0)') 'Model NetCDF4 Input "'//trim(nc_fname)// & + '" opened NF90_NOWRITE, ncid=', ncid + ! + ! -- return + return + end function nc_fopen + + subroutine nc_fclose(ncid, iout) + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: iout + ! -- local + ! + ! -- close netcdf file + call nf_verify(nf90_close(ncid), ncid, iout) + ! + ! -- return + return + end subroutine nc_fclose + + subroutine input_load(mf6_input, nc4pkg, ncid, nc_fname, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(NC4ModelPackageInputType), pointer, intent(in) :: nc4pkg + integer(I4B), intent(in) :: ncid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), dimension(:), contiguous, pointer :: mshape + integer(I4B) :: iblock + ! + ! -- set pointer to model shape + if (mf6_input%pkgtype(1:3) /= 'DIS') then + call mem_setptr(mshape, 'MODEL_SHAPE', mf6_input%component_mempath) + else + mshape => null() + end if + ! + ! -- log package load header + call idm_log_header(mf6_input%component_name, & + mf6_input%subcomponent_name, iout) + ! + ! -- load blocks until period + do iblock = 1, size(nc4pkg%blocklist) + if (nc4pkg%blocklist(iblock)%blockname == 'PERIOD') then + ! + exit + ! + else + ! + if (nc4pkg%blocklist(iblock)%aggregate) then + ! -- load list stype input + call load_aggregate_block(mf6_input, nc4pkg, mshape, & + iblock, ncid, nc_fname, iout) + else + ! -- load block variables + call load_block_var(mf6_input, nc4pkg, mshape, & + iblock, ncid, nc_fname, iout) + end if + ! + end if + end do + ! + ! -- close logging + call idm_log_close(mf6_input%component_name, & + mf6_input%subcomponent_name, iout) + ! + ! -- return + return + end subroutine input_load + + subroutine load_aggregate_block(mf6_input, nc4pkg, mshape, & + iblock, ncid, nc_fname, iout) + use InputOutputModule, only: parseline + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(NC4ModelPackageInputType), pointer, intent(in) :: nc4pkg + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), intent(in) :: iblock + integer(I4B), intent(in) :: ncid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + type(InputParamDefinitionType), pointer :: idt, idta + character(len=16), dimension(:), allocatable :: words + character(len=:), allocatable :: parse_str + integer(I4B) :: icol, nwords + logical(LGP) :: found + ! + ! -- set input definition for this block + idta => get_aggregate_definition_type(mf6_input%aggregate_dfns, & + mf6_input%component_type, & + mf6_input%subcomponent_type, & + mf6_input%block_dfns(iblock)%blockname) + ! + ! -- identify variable names + parse_str = trim(idta%datatype)//' ' + call parseline(parse_str, nwords, words) + ! + ! + ! -- Skip first col which is RECARRAY + do icol = 2, nwords + ! + found = load_aggregate_var(mf6_input, nc4pkg, mshape, iblock, idta, & + words(icol), ncid, nc_fname, iout) + ! + ! -- ensure required params found + if (.not. found) then + if (idt%required) then + errmsg = 'Required input parameter "'//trim(idt%tagname)// & + '" not found for package "'// & + trim(nc4pkg%subcomponent_name)//'".' + call store_error(errmsg) + call store_error_filename(nc_fname) + end if + end if + ! + end do + ! + ! -- cleanup + if (allocated(parse_str)) deallocate (parse_str) + if (allocated(words)) deallocate (words) + ! + ! -- return + return + end subroutine load_aggregate_block + + function load_aggregate_var(mf6_input, nc4pkg, mshape, iblock, idta, & + tagname, ncid, nc_fname, iout) result(found) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(NC4ModelPackageInputType), pointer, intent(in) :: nc4pkg + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), intent(in) :: iblock + type(InputParamDefinitionType), pointer, intent(in) :: idta + character(len=*), intent(in) :: tagname + integer(I4B), intent(in) :: ncid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- return + logical(LGP) :: found + ! -- local + type(InputParamDefinitionType), pointer :: idt + type(InputParamDefinitionType), target :: idtp + integer(I4B) :: iparam + ! + ! -- initialize + found = .false. + ! + do iparam = 1, nc4pkg%blocklist(iblock)%varnum + ! + ! -- set param definition + idt => get_param_definition_type(mf6_input%param_dfns, & + mf6_input%component_type, & + mf6_input%subcomponent_type, & + nc4pkg%blocklist(iblock)%blockname, & + tagname, nc_fname) + + if (nc4pkg%blocklist(iblock)%varnames(iparam) == trim(tagname)) then + ! + ! -- set found + found = .true. + ! + ! -- adjust definition to reflect list style input + idtp = aggregate_param_dfn(idta, idt, nc_fname) + idt => idtp + ! + ! -- load + if (idt%datatype == 'IRREGINT1D') then + call load_varint1d_type(mf6_input, mshape, idt, ncid, & + nc4pkg%blocklist(iblock)%varids(iparam), iout) + else + call load_static_var(mf6_input, idt, mshape, ncid, & + nc4pkg%blocklist(iblock)%varids(iparam), & + nc_fname, iout) + end if + ! + exit + end if + end do + ! + ! -- return + return + end function load_aggregate_var + + subroutine load_block_var(mf6_input, nc4pkg, mshape, & + iblock, ncid, nc_fname, iout) + use SourceCommonModule, only: set_model_shape, mem_allocate_naux + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(NC4ModelPackageInputType), pointer, intent(in) :: nc4pkg + integer(I4B), dimension(:), contiguous, pointer, intent(inout) :: mshape + integer(I4B), intent(in) :: iblock + integer(I4B), intent(in) :: ncid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + type(InputParamDefinitionType), pointer :: idt + integer(I4B) :: iparam + ! + ! -- load variables in block order + do iparam = 1, nc4pkg%blocklist(iblock)%varnum + ! + ! -- TODO: ensure nc file mf6_param is always using tagname + idt => & + get_param_definition_type(mf6_input%param_dfns, & + mf6_input%component_type, & + mf6_input%subcomponent_type, & + nc4pkg%blocklist(iblock)%blockname, & + nc4pkg%blocklist(iblock)%varnames(iparam), & + nc_fname) + ! -- load the variable data + call load_static_var(mf6_input, idt, mshape, ncid, & + nc4pkg%blocklist(iblock)%varids(iparam), & + nc_fname, iout) + ! + end do + ! + ! -- block post processing + select case (nc4pkg%blocklist(iblock)%blockname) + case ('OPTIONS') + ! -- allocate naux and set to 0 if not allocated + do iparam = 1, size(mf6_input%param_dfns) + idt => mf6_input%param_dfns(iparam) + ! + if (idt%blockname == 'OPTIONS' .and. & + idt%tagname == 'AUXILIARY') then + call mem_allocate_naux(mf6_input%mempath) + exit + end if + end do + case ('DIMENSIONS') + ! -- set model shape if discretization dimensions have been read + if (mf6_input%pkgtype(1:3) == 'DIS') then + call set_model_shape(mf6_input%pkgtype, nc_fname, & + mf6_input%component_mempath, & + mf6_input%mempath, mshape) + end if + case default + end select + ! + ! -- return + return + end subroutine load_block_var + + function aggregate_param_dfn(idta, idtp, nc_fname) result(idt) + type(InputParamDefinitionType), pointer, intent(in) :: idta + type(InputParamDefinitionType), pointer, intent(in) :: idtp + character(len=*), intent(in) :: nc_fname + type(InputParamDefinitionType) :: idt !< input data type object describing this record + ! + ! + ! -- copy from input dfn + idt%component_type = trim(idtp%component_type) + idt%subcomponent_type = trim(idtp%subcomponent_type) + idt%blockname = trim(idtp%blockname) + idt%tagname = trim(idtp%tagname) + idt%mf6varname = trim(idtp%mf6varname) + idt%shape = idta%shape + idt%required = idtp%required + idt%in_record = idtp%in_record + idt%preserve_case = idtp%preserve_case + idt%layered = idtp%layered + idt%timeseries = idtp%timeseries + ! + ! -- set datatype + select case (idtp%datatype) + case ('INTEGER') + idt%datatype = 'INTEGER1D' + case ('DOUBLE') + idt%datatype = 'DOUBLE1D' + case ('INTEGER1D') + if (idt%mf6varname == 'ICVERT') then + ! -- jagged array dependent on NCVERT 1d shape array + idt%datatype = 'IRREGINT1D' + idt%shape = idtp%shape + else + errmsg = 'IDM UNIMPLEMENTED LoadNetCDFData::aggregate_param_dfn & + &tagname='//trim(idt%tagname) + call store_error(errmsg) + call store_error_filename(nc_fname) + end if + case default + errmsg = 'IDM UNIMPLEMENTED LoadNetCDFData::aggregate_param_dfn & + &datatype='//trim(idtp%datatype) + call store_error(errmsg) + call store_error_filename(nc_fname) + end select + ! + ! -- return + return + end function aggregate_param_dfn + + function format_input_str(input_str, len) result(str) + use InputOutputModule, only: lowcase, upcase + character(len=*), intent(in) :: input_str + integer(I4B), intent(in) :: len + character(len=len) :: str + integer(I4B) :: idx + ! + ! -- initialize + idx = 0 + str = trim(input_str) + ! + ! -- identify '"' character if present + idx = index(str, '"') + ! + if (idx > 0) then + str = str(1:idx - 1) + end if + ! + call upcase(str) + ! + ! -- return + return + end function format_input_str + + subroutine load_static_var(mf6_input, idt, mshape, ncid, varid, nc_fname, iout) + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + ! + ! -- allocate and load data type + select case (idt%datatype) + case ('KEYWORD') + call load_keyword_type(mf6_input, idt, ncid, varid, iout) + case ('STRING') + if (idt%shape == 'NAUX') then + call load_auxvar_names(mf6_input, idt, ncid, varid, iout) + else + call load_string_type(mf6_input, idt, ncid, varid, iout) + end if + case ('INTEGER') + call load_integer_type(mf6_input, idt, ncid, varid, iout) + case ('INTEGER1D') + call load_integer1d_type(mf6_input, mshape, idt, ncid, & + varid, nc_fname, iout) + case ('INTEGER2D') + call load_integer2d_type(mf6_input, mshape, idt, ncid, varid, iout) + case ('INTEGER3D') + call load_integer3d_type(mf6_input, mshape, idt, ncid, varid, iout) + case ('DOUBLE') + call load_double_type(mf6_input, idt, ncid, varid, iout) + case ('DOUBLE1D') + call load_double1d_type(mf6_input, mshape, idt, ncid, varid, nc_fname, iout) + case ('DOUBLE2D') + call load_double2d_type(mf6_input, mshape, idt, ncid, varid, iout) + case ('DOUBLE3D') + call load_double3d_type(mf6_input, mshape, idt, ncid, varid, iout) + case default + if (idt%datatype(1:6) == 'RECORD') then + call load_record_type(mf6_input, idt, ncid, varid, nc_fname, iout) + else + errmsg = 'IDM UNIMPLEMENTED LoadNetCDFData::load_static_var & + &datatype='//trim(idt%datatype) + call store_error(errmsg) + call store_error_filename(nc_fname) + end if + end select + ! + ! -- return + return + end subroutine load_static_var + + !> @brief load type keyword + !< + subroutine load_keyword_type(mf6_input, idt, ncid, varid, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), pointer :: intvar + ! + call mem_allocate(intvar, idt%mf6varname, mf6_input%mempath) + intvar = 1 + call idm_log_var(intvar, idt%mf6varname, mf6_input%mempath, & + idt%datatype, iout) + ! + return + end subroutine load_keyword_type + + !> @brief load type record + !! + !! Currently only supports file record types + !< + subroutine load_record_type(mf6_input, idt, ncid, varid, nc_fname, iout) + use InputOutputModule, only: parseline + use MemoryManagerModule, only: get_isize + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + type(InputParamDefinitionType), pointer :: io_idt + character(len=40), dimension(:), allocatable :: words + character(len=:), allocatable :: parse_str + character(len=LINELENGTH) :: cstr + type(CharacterStringType), dimension(:), contiguous, & + pointer :: charstr1d + integer(I4B) :: nwords, isize, idx + logical(LGP) :: found + ! + ! -- initialization + found = .false. + ! + ! -- set split string + parse_str = trim(idt%datatype)//' ' + ! + ! -- split + call parseline(parse_str, nwords, words) + ! + ! -- a filein/fileout record tag definition has 4 tokens + if (nwords == 4) then + ! + ! -- verify third definition token is FILEIN/FILEOUT + if (words(3) == 'FILEIN' .or. words(3) == 'FILEOUT') then + ! + ! -- matches, read and load file name + io_idt => & + get_param_definition_type(mf6_input%param_dfns, & + mf6_input%component_type, & + mf6_input%subcomponent_type, & + 'OPTIONS', words(4), nc_fname) + ! + ! -- io tag loaded + found = .true. + end if + end if + ! + ! -- load the data if found + if (found) then + ! + if (words(3) == 'FILEIN') then + ! + ! -- FILEIN record types support 1+ file specifications + call get_isize(io_idt%mf6varname, mf6_input%mempath, isize) + ! + ! -- allocate / reallocate + if (isize < 0) then + ! -- does not exist, allocate + call mem_allocate(charstr1d, LINELENGTH, 1, io_idt%mf6varname, & + mf6_input%mempath) + idx = 1 + else + ! -- exists, reallocate + call mem_setptr(charstr1d, io_idt%mf6varname, mf6_input%mempath) + call mem_reallocate(charstr1d, LINELENGTH, isize + 1, & + io_idt%mf6varname, mf6_input%mempath) + idx = isize + 1 + end if + ! + ! -- read and set the data + call nf_verify(nf90_get_var(ncid, varid, cstr), ncid, iout) + charstr1d(idx) = trim(cstr) + ! + else if (words(3) == 'FILEOUT') then + ! -- FILEOUT record types support a single file specification + call load_string_type(mf6_input, io_idt, ncid, varid, iout) + end if + ! + else + errmsg = 'IDM UNIMPLEMENTED LoadNetCDFData::load_record_type & + &tagname='//trim(idt%tagname) + call store_error(errmsg) + call store_error_filename(nc_fname) + end if + ! + ! -- cleanup + if (allocated(parse_str)) deallocate (parse_str) + if (allocated(words)) deallocate (words) + ! + return + end subroutine load_record_type + + !> @brief load type string + !< + subroutine load_string_type(mf6_input, idt, ncid, varid, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + character(len=LINELENGTH), pointer :: cstr + ! + call mem_allocate(cstr, LINELENGTH, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, cstr), ncid, iout) + call idm_log_var(cstr, idt%mf6varname, mf6_input%mempath, iout) + ! + return + end subroutine load_string_type + + subroutine load_auxvar_names(mf6_input, idt, ncid, varid, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + character(len=LENAUXNAME), dimension(:), pointer :: cstr1d + type(CharacterStringType), dimension(:), & + pointer, contiguous :: charstr1d !< variable for allocation + integer(I4B), pointer :: naux + character(len=LINELENGTH) :: varname, dimname + integer(I4B) :: vartype, ndims, nattr, n + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM) :: dimids, dimlens + ! + ! -- inquire for variable info + call nf_verify(nf90_inquire_variable(ncid, varid, varname, vartype, ndims, & + dimids, nattr), ncid, iout) + ! + ! -- set variable dimensions + do n = 1, ndims + ! -- dim 1 is strlen, dim 2 is naux + call nf_verify(nf90_inquire_dimension(ncid, dimids(n), dimname, & + dimlens(n)), ncid, iout) + end do + ! + ! -- allocate and set naux + call mem_allocate(naux, idt%shape, mf6_input%mempath) + naux = dimlens(2) + ! + ! -- allocate local string array + allocate (cstr1d(naux)) + ! + ! -- load data + call nf_verify(nf90_get_var(ncid, varid, cstr1d), ncid, iout) + ! + ! -- allocate charstring array + call mem_allocate(charstr1d, LENAUXNAME, naux, & + idt%mf6varname, mf6_input%mempath) + ! + ! --copy local string input to charstr array + do n = 1, naux + charstr1d(n) = trim(cstr1d(n)) + end do + ! + ! -- cleanup + deallocate (cstr1d) + ! + return + end subroutine load_auxvar_names + + !> @brief load type integer + !< + subroutine load_integer_type(mf6_input, idt, ncid, varid, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), pointer :: intvar + ! + call mem_allocate(intvar, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, intvar), ncid, iout) + call idm_log_var(intvar, idt%mf6varname, mf6_input%mempath, & + idt%datatype, iout) + ! + return + end subroutine load_integer_type + + !> @brief load type 1d integer + !< + subroutine load_integer1d_type(mf6_input, mshape, idt, ncid, & + varid, nc_fname, iout) + use SourceCommonModule, only: get_shape_from_string + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), dimension(:), contiguous, pointer :: int1d + integer(I4B), dimension(:, :), contiguous, pointer :: int2d + integer(I4B), dimension(:, :, :), contiguous, pointer :: int3d + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B) :: nvals, i, j, k, n + ! + ! -- initialize + n = 0 + nvals = 0 + ! + if (idt%shape == 'NODES') then + ! -- set number of values + nvals = product(mshape) + ! + ! -- allocate managed memory + call mem_allocate(int1d, nvals, idt%mf6varname, mf6_input%mempath) + ! + if (size(mshape) == 3) then + ! -- allocate local array + allocate (int3d(mshape(3), mshape(2), mshape(1))) + ! -- read the data into local array + call nf_verify(nf90_get_var(ncid, varid, int3d), ncid, iout) + ! + ! -- copy to flat array + do k = 1, size(int3d, dim=3) + do j = 1, size(int3d, dim=2) + do i = 1, size(int3d, dim=1) + n = n + 1 + int1d(n) = int3d(i, j, k) + end do + end do + end do + ! + ! -- deallocate local array + deallocate (int3d) + ! + else if (size(mshape) == 2) then + ! -- allocate local array + allocate (int2d(mshape(2), mshape(1))) + ! -- read the data into local array + call nf_verify(nf90_get_var(ncid, varid, int2d), ncid, iout) + ! + ! -- copy to flat array + do j = 1, size(int2d, dim=2) + do i = 1, size(int2d, dim=1) + n = n + 1 + int1d(n) = int2d(i, j) + end do + end do + ! + ! -- deallocate local array + deallocate (int2d) + ! + else + write (errmsg, '(a,i0)') & + 'IDM UNIMPLEMENTED LoadNetCDFData::load_integer1d_type & + &size mshape=', size(mshape) + call store_error(errmsg) + call store_error_filename(nc_fname) + end if + else + ! -- interpret shape + call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath) + ! + ! -- set nvals + nvals = array_shape(1) + ! + ! -- allocate managed memory + call mem_allocate(int1d, nvals, idt%mf6varname, mf6_input%mempath) + ! + ! -- read and set data + call nf_verify(nf90_get_var(ncid, varid, int1d), ncid, iout) + ! + end if + ! + ! -- log variable + call idm_log_var(int1d, idt%mf6varname, mf6_input%mempath, iout) + ! + return + end subroutine load_integer1d_type + + !> @brief load intvector type + !< + subroutine load_varint1d_type(mf6_input, mshape, idt, ncid, varid, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), dimension(:), pointer, contiguous :: int1d, intvector_shape + character(len=LINELENGTH) :: varname, dimname + integer(I4B) :: vartype, ndims, nattr, n, numvals + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM) :: dimids, dimlens + ! + ! -- inquire for variable info + call nf_verify(nf90_inquire_variable(ncid, varid, varname, vartype, ndims, & + dimids, nattr), ncid, iout) + ! + ! -- set variable dimensions + do n = 1, ndims + ! + call nf_verify(nf90_inquire_dimension(ncid, dimids(n), dimname, & + dimlens(n)), ncid, iout) + end do + ! + ! -- set pointer to shape array + call mem_setptr(intvector_shape, idt%shape, mf6_input%mempath) + ! + ! -- set total number of values + numvals = sum(intvector_shape) + ! + ! -- TODO check numvals against dimlens(1)? + ! + ! -- allocate managed memory + call mem_allocate(int1d, dimlens(1), idt%mf6varname, mf6_input%mempath) + ! + ! -- read data into managed memory + call nf_verify(nf90_get_var(ncid, varid, int1d), ncid, iout) + ! + ! -- log variable + call idm_log_var(int1d, idt%mf6varname, mf6_input%mempath, iout) + ! + ! -- return + return + end subroutine load_varint1d_type + + !> @brief load type 2d integer + !< + subroutine load_integer2d_type(mf6_input, mshape, idt, ncid, varid, iout) + use SourceCommonModule, only: get_shape_from_string + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), dimension(:, :), pointer, contiguous :: int2d + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B) :: nsize1, nsize2 + ! + ! -- determine the array shape from the input data defintion + ! idt%shape, which looks like "NCOL, NROW, NLAY" + call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) + ! + ! -- allocate, read and log + call mem_allocate(int2d, nsize1, nsize2, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, int2d), ncid, iout) + call idm_log_var(int2d, idt%mf6varname, mf6_input%mempath, iout) + ! + ! -- return + return + end subroutine load_integer2d_type + + !> @brief load type 3d integer + !< + subroutine load_integer3d_type(mf6_input, mshape, idt, ncid, varid, iout) + use SourceCommonModule, only: get_shape_from_string + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B), dimension(:, :, :), pointer, contiguous :: int3d + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B) :: nsize1, nsize2, nsize3 + ! + ! -- determine the array shape from the input data defintion + ! idt%shape, which looks like "NCOL, NROW, NLAY" + call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) + nsize3 = array_shape(3) + ! + ! -- allocate, read and log + call mem_allocate(int3d, nsize1, nsize2, nsize3, idt%mf6varname, & + mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, int3d), ncid, iout) + call idm_log_var(int3d, idt%mf6varname, mf6_input%mempath, iout) + ! + ! -- return + return + end subroutine load_integer3d_type + + !> @brief load type double + !< + subroutine load_double_type(mf6_input, idt, ncid, varid, iout) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + real(DP), pointer :: dblvar + ! + call mem_allocate(dblvar, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, dblvar), ncid, iout) + call idm_log_var(dblvar, idt%mf6varname, mf6_input%mempath, iout) + ! + ! -- return + return + end subroutine load_double_type + + !> @brief load type 1d double + !< + subroutine load_double1d_type(mf6_input, mshape, idt, ncid, & + varid, nc_fname, iout) + use SourceCommonModule, only: get_shape_from_string + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + real(DP), dimension(:), contiguous, pointer :: dbl1d + real(DP), dimension(:, :), contiguous, pointer :: dbl2d + real(DP), dimension(:, :, :), contiguous, pointer :: dbl3d + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B) :: nvals, i, j, k, n + ! + ! -- initialize + n = 0 + nvals = 0 + ! + if (idt%shape == 'NODES') then + ! -- set number of values + nvals = product(mshape) + ! + ! -- allocate managed memory + call mem_allocate(dbl1d, nvals, idt%mf6varname, mf6_input%mempath) + ! + if (size(mshape) == 3) then + ! -- allocate local array + allocate (dbl3d(mshape(3), mshape(2), mshape(1))) + ! -- read the data into local array + call nf_verify(nf90_get_var(ncid, varid, dbl3d), ncid, iout) + ! + ! -- copy to flat array + do k = 1, size(dbl3d, dim=3) + do j = 1, size(dbl3d, dim=2) + do i = 1, size(dbl3d, dim=1) + n = n + 1 + dbl1d(n) = dbl3d(i, j, k) + end do + end do + end do + ! + ! -- deallocate local array + deallocate (dbl3d) + ! + else if (size(mshape) == 2) then + ! -- allocate local array + allocate (dbl2d(mshape(2), mshape(1))) + ! -- read the data into local array + call nf_verify(nf90_get_var(ncid, varid, dbl2d), ncid, iout) + ! + ! -- copy to flat array + do j = 1, size(dbl2d, dim=2) + do i = 1, size(dbl2d, dim=1) + n = n + 1 + dbl1d(n) = dbl2d(i, j) + end do + end do + ! + ! -- deallocate local array + deallocate (dbl2d) + ! + else + write (errmsg, '(a,i0)') & + 'IDM UNIMPLEMENTED LoadNetCDFData::load_double1d_type & + &size mshape=', size(mshape) + call store_error(errmsg) + call store_error_filename(nc_fname) + end if + else + ! -- interpret shape + call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath) + ! + ! -- set nvals + nvals = array_shape(1) + ! + ! -- allocate managed memory + call mem_allocate(dbl1d, nvals, idt%mf6varname, mf6_input%mempath) + ! + ! -- read and set data + call nf_verify(nf90_get_var(ncid, varid, dbl1d), ncid, iout) + end if + ! + ! -- log variable + call idm_log_var(dbl1d, idt%mf6varname, mf6_input%mempath, iout) + ! + ! -- return + return + end subroutine load_double1d_type + + !> @brief load type 2d double + !< + subroutine load_double2d_type(mf6_input, mshape, idt, ncid, varid, iout) + use SourceCommonModule, only: get_shape_from_string + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + real(DP), dimension(:, :), pointer, contiguous :: dbl2d + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B) :: nsize1, nsize2 + ! + ! -- determine the array shape from the input dfn shape + call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) + ! + call mem_allocate(dbl2d, nsize1, nsize2, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, dbl2d), ncid, iout) + call idm_log_var(dbl2d, idt%mf6varname, mf6_input%mempath, iout) + ! + ! -- return + return + end subroutine load_double2d_type + + !> @brief load type 3d double + !< + subroutine load_double3d_type(mf6_input, mshape, idt, ncid, varid, iout) + use SourceCommonModule, only: get_shape_from_string + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + type(InputParamDefinitionType), intent(in) :: idt + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + real(DP), dimension(:, :, :), pointer, contiguous :: dbl3d + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B) :: nsize1, nsize2, nsize3 + ! + ! -- determine the array shape from the input dfn shape + call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) + nsize3 = array_shape(3) + ! + call mem_allocate(dbl3d, nsize1, nsize2, nsize3, idt%mf6varname, & + mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, dbl3d), ncid, iout) + call idm_log_var(dbl3d, idt%mf6varname, mf6_input%mempath, iout) + ! + return + end subroutine load_double3d_type + + subroutine nc4_rp_list(mf6_input, ncpkg, bndctx, ncid, iiper, nc_fname, iout) + use TdisModule, only: kper + use BoundInputContextModule, only: BoundInputContextType + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(NC4ModelPackageInputType), intent(in) :: ncpkg + type(BoundInputContextType), intent(in) :: bndctx + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: iiper + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + type(InputParamDefinitionType), pointer :: idt + integer(I4B), dimension(:), contiguous, pointer :: mshape + integer(I4B) :: iparam, iblock, nbound, prev_nbound, varid + character(len=LINELENGTH) :: tagname + ! + ! -- initialize + nbound = 0 + prev_nbound = 0 + ! + ! -- set mshape input context pointer + call mem_setptr(mshape, 'MODEL_SHAPE', mf6_input%component_mempath) + ! + do iblock = 1, size(ncpkg%blocklist) + if (ncpkg%blocklist(iblock)%blockname == 'PERIOD') then + do iparam = 1, ncpkg%blocklist(iblock)%varnum + ! + tagname = ncpkg%blocklist(iblock)%varnames(iparam) + ! + idt => & + get_param_definition_type(mf6_input%param_dfns, & + mf6_input%component_type, & + mf6_input%subcomponent_type, & + 'PERIOD', tagname, nc_fname) + ! + varid = ncpkg%blocklist(iblock)%varids(iparam) + nbound = load_var_rp_list(mf6_input, ncid, idt, mshape, & + varid, iiper, nc_fname, iout) + ! + if (prev_nbound /= 0 .and. nbound /= prev_nbound) then + write (errmsg, '(a,i0,a)') & + 'Package input nbound inconsistency in period ', kper, & + ' for data in package "'//trim(ncpkg%subcomponent_name)//'". ' + call store_error(errmsg) + call store_error_filename(nc_fname) + else + prev_nbound = nbound + end if + ! + end do + end if + end do + ! + ! -- set nbound + bndctx%nbound = nbound + ! + ! -- return + return + end subroutine nc4_rp_list + + function load_var_rp_list(mf6_input, ncid, idt, mshape, & + varid, iiper, nc_fname, iout) result(nbound) + use BoundInputContextModule, only: BoundInputContextType + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iiper + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- result + integer(I4B) :: nbound + ! -- local + integer(I4B) :: n + character(len=LINELENGTH) :: varname, dimname + integer(I4B) :: vartype, ndims, nattr + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM) :: dimids, dimlens + ! + ! -- initialize + nbound = 0 + ! + ! -- inquire for variable info + call nf_verify(nf90_inquire_variable(ncid, varid, varname, vartype, ndims, & + dimids, nattr), ncid, iout) + ! + ! -- set variable dimensions + do n = 1, ndims + ! + call nf_verify(nf90_inquire_dimension(ncid, dimids(n), dimname, & + dimlens(n)), ncid, iout) + end do + ! + ! -- read and load variable data + select case (idt%datatype) + case ('STRING') + nbound = load_charstr1d_rp_list(mf6_input, ncid, idt, mshape, & + dimlens, varid, iiper, iout) + ! + case ('INTEGER') + nbound = load_int1d_rp_list(mf6_input, ncid, idt, mshape, & + dimlens, varid, iiper, iout) + ! + case ('INTEGER1D') + nbound = load_int2d_rp_list(mf6_input, ncid, idt, mshape, & + dimlens, varid, iiper, iout) + ! + case ('DOUBLE') + nbound = load_dbl1d_rp_list(mf6_input, ncid, idt, mshape, & + dimlens, varid, iiper, iout) + ! + case ('DOUBLE1D') + nbound = load_dbl2d_rp_list(mf6_input, ncid, idt, mshape, & + dimlens, varid, iiper, iout) + ! + case default + errmsg = 'IDM UNIMPLEMENTED LoadNetCDFData::load_var_rp_list & + &datatype='//trim(idt%datatype) + call store_error(errmsg) + call store_error_filename(nc_fname) + end select + ! + ! -- return + return + end function load_var_rp_list + + function load_charstr1d_rp_list(mf6_input, ncid, idt, mshape, dimlens, & + varid, iiper, iout) result(nbound) + use InputOutputModule, only: lowcase, upcase + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iiper + integer(I4B), intent(in) :: iout + ! -- result + integer(I4B) :: nbound + ! -- local + type(CharacterStringType), dimension(:), & + contiguous, pointer :: charstr1d + character(len=LENBOUNDNAME), dimension(:), & + contiguous, pointer :: cstr1d + integer(I4B) :: n + ! + ! -- initialize + nbound = 0 + ! + ! -- allocate local character string array + allocate (cstr1d(dimlens(2))) + ! + ! -- read string array + call nf_verify(nf90_get_var(ncid, varid, cstr1d, & + start=(/1, iiper/), & + count=(/dimlens(1), dimlens(2), 1/)), ncid, iout) + ! + ! -- copy local data to managed memory + call mem_setptr(charstr1d, idt%mf6varname, mf6_input%mempath) + ! + do n = 1, size(cstr1d) + charstr1d(n) = format_input_str(cstr1d(n), LENBOUNDNAME) + if (charstr1d(n) == '') then + nbound = n - 1 + exit + end if + end do + ! + if (nbound == 0) then + nbound = size(cstr1d) + end if + ! + ! -- deallocate local array + deallocate (cstr1d) + ! + ! -- return + return + end function load_charstr1d_rp_list + + function load_int1d_rp_list(mf6_input, ncid, idt, mshape, dimlens, & + varid, iiper, iout) result(nbound) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iiper + integer(I4B), intent(in) :: iout + ! -- result + integer(I4B) :: nbound + ! -- local + integer(I4B), dimension(:), pointer, contiguous :: int1d + integer(I4B) :: n + ! + nbound = 0 + ! + ! -- read and load data to managed memory + call mem_setptr(int1d, idt%mf6varname, mf6_input%mempath) + ! + call nf_verify(nf90_get_var(ncid, varid, int1d, & + start=(/1, iiper/), & + count=(/dimlens(1), 1/)), ncid, iout) + ! + ! -- set nbound + do n = 1, size(int1d) + if (int1d(n) == INODATA) then + nbound = n - 1 + exit + end if + end do + ! + if (nbound == 0) then + nbound = size(int1d) + end if + ! + ! -- return + return + end function load_int1d_rp_list + + function load_int2d_rp_list(mf6_input, ncid, idt, mshape, dimlens, & + varid, iiper, iout) result(nbound) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iiper + integer(I4B), intent(in) :: iout + ! -- result + integer(I4B) :: nbound + ! -- local + integer(I4B), dimension(:, :), pointer, contiguous :: int2d + integer(I4B) :: n, m + ! + nbound = 0 + ! + ! -- read and load data to managed memory + call mem_setptr(int2d, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, int2d, & + start=(/1, 1, iiper/), & + count=(/dimlens(1), dimlens(2), 1/)), ncid, iout) + ! + ! -- set nbound + outer: do n = 1, size(int2d, 2) + do m = 1, size(int2d, 1) + if (int2d(m, n) == INODATA) then + nbound = n - 1 + exit outer + end if + end do + end do outer + ! + if (nbound == 0) then + nbound = size(int2d, 2) + end if + ! + ! -- return + return + end function load_int2d_rp_list + + function load_dbl1d_rp_list(mf6_input, ncid, idt, mshape, dimlens, & + varid, iiper, iout) result(nbound) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iiper + integer(I4B), intent(in) :: iout + ! -- result + integer(I4B) :: nbound + ! -- local + real(DP), dimension(:), pointer, contiguous :: dbl1d + integer(I4B) :: n + ! + nbound = 0 + ! + ! -- read and load data to managed memory + call mem_setptr(dbl1d, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, dbl1d, & + start=(/1, iiper/), & + count=(/dimlens(1), 1/)), ncid, iout) + ! + ! -- set nbound + do n = 1, size(dbl1d) + if (dbl1d(n) == DNODATA) then + nbound = n - 1 + exit + end if + end do + ! + if (nbound == 0) then + nbound = size(dbl1d) + end if + ! + ! -- return + return + end function load_dbl1d_rp_list + + function load_dbl2d_rp_list(mf6_input, ncid, idt, mshape, dimlens, & + varid, iiper, iout) result(nbound) + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iiper + integer(I4B), intent(in) :: iout + ! -- result + integer(I4B) :: nbound + ! -- local + real(DP), dimension(:, :), pointer, contiguous :: dbl2d + integer(I4B) :: n, m + ! + nbound = 0 + ! + ! -- read and load data to managed memory + call mem_setptr(dbl2d, idt%mf6varname, mf6_input%mempath) + call nf_verify(nf90_get_var(ncid, varid, dbl2d, & + start=(/1, 1, iiper/), & + count=(/dimlens(1), dimlens(2), 1/)), ncid, iout) + ! + ! -- set nbound + outer: do n = 1, size(dbl2d, 2) + do m = 1, size(dbl2d, 1) + if (dbl2d(m, n) == DNODATA) then + nbound = n - 1 + exit outer + end if + end do + end do outer + ! + if (nbound == 0) then + nbound = size(dbl2d, 2) + end if + ! + ! -- return + return + end function load_dbl2d_rp_list + + subroutine nc4_rp_array(mf6_input, ncpkg, param_names, param_reads, & + bndctx, ncid, nc_fname, iout) + use BoundInputContextModule, only: BoundInputContextType + use SourceCommonModule, only: ReadStateVarType + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + type(NC4ModelPackageInputType), intent(in) :: ncpkg + character(len=LENVARNAME), dimension(:), & + allocatable :: param_names !< dynamic param names + type(ReadStateVarType), dimension(:), allocatable :: param_reads + type(BoundInputContextType), intent(in) :: bndctx + integer(I4B), intent(in) :: ncid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + type(InputParamDefinitionType), pointer :: idt + integer(I4B) :: iparam, iblock, varid + character(len=LINELENGTH) :: tagname + ! + do iblock = 1, size(ncpkg%blocklist) + if (ncpkg%blocklist(iblock)%blockname == 'PERIOD') then + do iparam = 1, ncpkg%blocklist(iblock)%varnum + ! + ! -- set varid + varid = ncpkg%blocklist(iblock)%varids(iparam) + ! + ! -- set tagname + tagname = ncpkg%blocklist(iblock)%varnames(iparam) + ! + idt => & + get_param_definition_type(mf6_input%param_dfns, & + mf6_input%component_type, & + mf6_input%subcomponent_type, & + 'PERIOD', tagname, nc_fname) + ! + call load_var_rp_array(mf6_input, param_names, & + param_reads, ncid, idt, varid, & + nc_fname, iout) + end do + end if + end do + ! + ! -- return + return + end subroutine nc4_rp_array + + subroutine load_var_rp_array(mf6_input, param_names, param_reads, & + ncid, idt, varid, nc_fname, iout) + use TdisModule, only: kper + use BoundInputContextModule, only: BoundInputContextType + use SourceCommonModule, only: ReadStateVarType + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + character(len=LENVARNAME), dimension(:), & + allocatable :: param_names !< dynamic param names + type(ReadStateVarType), dimension(:), allocatable :: param_reads + integer(I4B), intent(in) :: ncid + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), intent(in) :: varid + character(len=*), intent(in) :: nc_fname + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B) :: n, load_idx, len + character(len=LINELENGTH) :: varname, dimname + integer(I4B) :: vartype, ndims, nattr + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM) :: dimids, dimlens + integer(I4B), dimension(:), allocatable :: ipers + ! + ! -- index for loading kper data + load_idx = 0 + ! + ! -- inquire for variable info + call nf_verify(nf90_inquire_variable(ncid, varid, varname, vartype, ndims, & + dimids, nattr), ncid, iout) + ! + ! -- set variable dimensions + do n = 1, ndims + ! + ! -- (1+: data, last-1: maxbound, last:ipers) + call nf_verify(nf90_inquire_dimension(ncid, dimids(n), dimname, & + dimlens(n)), ncid, iout) + end do + ! + ! -- set current load index + if (nf90_inquire_attribute(ncid, varid, "mf6_iper", & + len=len) == NF90_NOERR) then + ! -- allocate local array + allocate (ipers(len)) + ! -- read variable ipers array + if (nf90_get_att(ncid, varid, "mf6_iper", & + ipers) == NF90_NOERR) then + ! -- check load periods + do n = 1, size(ipers) + ! + if (ipers(n) == kper) then + ! -- load, save index + load_idx = n + exit + end if + end do + end if + end if + ! + ! -- read and load variable data if load_idx defined + if (load_idx > 0) then + select case (idt%datatype) + !case ('INTEGER1D') + ! call load_int3d_rp_array(mf6_input, ncid, idt, & + ! dimlens, varid, load_idx, iout) + ! + case ('DOUBLE1D') + call load_dbl3d_rp_array(mf6_input, ncid, param_names, param_reads, & + idt, dimlens, varid, load_idx, iout) + ! + case default + errmsg = 'IDM UNIMPLEMENTED LoadNetCDFData::load_var_rp_array & + &datatype='//trim(idt%datatype) + call store_error(errmsg) + call store_error_filename(nc_fname) + end select + end if + ! + ! -- return + return + end subroutine load_var_rp_array + +! subroutine load_int3d_rp_array(mf6_input, ncid, idt, & +! dimlens, varid, load_idx, iout) +! ! -- dummy +! type(ModflowInputType), intent(in) :: mf6_input +! integer(I4B), intent(in) :: ncid +! type(InputParamDefinitionType), pointer, intent(in) :: idt +! integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens +! integer(I4B), intent(in) :: varid +! integer(I4B), intent(in) :: load_idx +! integer(I4B), intent(in) :: iout +! ! -- local +! integer(I4B), dimension(:), pointer, contiguous :: int1d +! ! +! ! -- return +! return +! end subroutine load_int3d_rp_array + + subroutine load_dbl3d_rp_array(mf6_input, ncid, param_names, param_reads, & + idt, dimlens, varid, load_idx, iout) + ! -- modules + use SourceCommonModule, only: ReadStateVarType + use ArrayHandlersModule, only: ifind + use SourceCommonModule, only: ifind_charstr + ! -- dummy + type(ModflowInputType), intent(in) :: mf6_input + integer(I4B), intent(in) :: ncid + character(len=LENVARNAME), dimension(:), & + allocatable :: param_names !< dynamic param names + type(ReadStateVarType), dimension(:), allocatable :: param_reads + type(InputParamDefinitionType), pointer, intent(in) :: idt + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM), intent(in) :: dimlens + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: load_idx + integer(I4B), intent(in) :: iout + ! -- local + real(DP), dimension(:), pointer, contiguous :: dbl1d + real(DP), dimension(:, :, :), pointer, contiguous :: dbl3d + integer(I4B), dimension(:), pointer, contiguous :: ialayer + integer(I4B) :: n, i, j, k, layer, iparam, ilayer + ! + ! -- initialize + n = 0 + layer = 1 + ilayer = 0 + nullify (ialayer) + ! + ! -- TEMP set layer index array + ilayer = ifind(param_names, 'IRCH') + if (ilayer <= 0) iparam = ifind(param_names, 'IEVT') + if (ilayer > 0) then + call mem_setptr(ialayer, param_names(ilayer), mf6_input%mempath) + end if + ! + ! -- set iparam to name index + iparam = ifind(param_names, idt%tagname) + ! + ! -- set pointer to managed memory variable + call mem_setptr(dbl1d, idt%mf6varname, mf6_input%mempath) + ! + ! -- allocate local array + allocate (dbl3d(dimlens(3), dimlens(2), dimlens(1))) + ! + ! -- read and load data to local array + ! -- dims: (1: ncol, 2: nrow, 3: nlay, 4: iper) + call nf_verify(nf90_get_var(ncid, varid, dbl3d, & + start=(/1, 1, 1, load_idx/), & + count=(/dimlens(1), dimlens(2), & + dimlens(3), 1/)), ncid, iout) + ! + ! -- copy data from local array to managed memory + do k = 1, size(dbl3d, dim=3) + do j = 1, size(dbl3d, dim=2) + do i = 1, size(dbl3d, dim=1) + if (n < size(dbl1d)) then + n = n + 1 + else + n = 1 + layer = layer + 1 + end if + if (dbl3d(i, j, k) /= DNODATA) then + dbl1d(n) = dbl3d(i, j, k) + if (ilayer > 0) ialayer(n) = layer + end if + end do + end do + end do + ! + ! -- set read state vars + if (iparam > 0) then + param_reads(iparam)%invar = 1 + end if + if (ilayer > 0) then + param_reads(ilayer)%invar = 1 + end if + ! + ! -- deallocate local array + deallocate (dbl3d) + ! + ! -- return + return + end subroutine load_dbl3d_rp_array + + subroutine nf_verify(res, ncid, iout) + integer(I4B), intent(in) :: res + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: iout + ! -- local variables + character(len=LINELENGTH) :: errstr + ! + ! -- strings are set for a subset of errors + ! but the exit status will always be reported + if (res /= NF90_NOERR) then + ! + select case (res) + case (NF90_EINVAL) ! (-36) + errstr = 'Invalid Argument' + case (NF90_EPERM) ! (-37) + errstr = 'Write to read only' + case (NF90_EINVALCOORDS) ! (-40) + errstr = 'Index exceeds dimension bound' + case (NF90_ENAMEINUSE) ! (-42) + errstr = 'String match to name in use' + case (NF90_ENOTATT) ! (-43) + errstr = 'Attribute not found' + case (NF90_EBADDIM) ! (-46) + errstr = 'Invalid dimension id or name' + case (NF90_ENOTVAR) ! (-49) + errstr = 'Variable not found' + case (NF90_ENOTNC) ! (-51) + errstr = 'Not a netcdf file' + case (NF90_ECHAR) ! (-56) + errstr = 'Attempt to convert between text & numbers' + case (NF90_EEDGE) ! (-57) + errstr = 'Edge+start exceeds dimension bound' + case (NF90_ESTRIDE) ! (-58) + errstr = 'Illegal stride' + case (NF90_EBADNAME) ! (-59) + errstr = 'Attribute or variable name contains illegal characters' + case default + errstr = '' + end select + ! + if (errstr /= '') then + write (errmsg, '(a,a,a,i0,a,i0,a)') 'NETCDF4 library error [error="', & + trim(errstr), '", exit code=', res, ', ncid=', ncid, '].' + else + write (errmsg, '(a,i0,a,i0,a)') 'NETCDF4 library error [exit code=', & + res, ', ncid=', ncid, '].' + end if + ! + call store_error(errmsg, .true.) + end if + ! + ! -- return + return + end subroutine nf_verify + + !> @brief get package list + !< + function get_nc4_package(nc4_context, modelname, scname, & + ptype, sctype, iout) result(ncpkg) + ! -- dummy variables + type(NC4ModelInputsType), intent(in) :: nc4_context + character(len=*), intent(in) :: modelname + character(len=*), intent(in) :: scname + character(len=*), intent(in) :: ptype + character(len=*), intent(in) :: sctype + integer(I4B), intent(in) :: iout + ! -- return + class(NC4ModelPackageInputType), pointer :: ncpkg + ! -- local variables + class(NC4ModelPackageInputType), pointer :: pkg + integer(I4B) :: n + ! + ! -- initialize + ncpkg => null() + ! + do n = 1, nc4_context%pkglist%count() + pkg => nc4_context%get(n) + if (pkg%subcomponent_name == scname .and. & + pkg%component_name == modelname) then + ncpkg => pkg + exit + end if + end do + ! + if (.not. associated(ncpkg)) then + errmsg = 'NC package context not found:'//trim(modelname)//'/'//trim(scname) + call store_error(errmsg) + call store_error_filename(nc4_context%modelfname) + else + if (ncpkg%subcomponent_type /= sctype) then + ! -- pkgtype (EVT6 or RCH6) used in namefile packages block- + ! dfns were initialized to list but need to be array. + call ncpkg%reset(ptype, sctype) + end if + end if + ! + ! -- return + return + end function get_nc4_package + + !> @brief get package list + !< + function read_as_arrays(param_dfns) + ! -- dummy variables + type(InputParamDefinitionType), dimension(:), & + pointer, intent(in) :: param_dfns + ! -- return + logical(LGP) :: read_as_arrays + ! -- local variables + type(InputParamDefinitionType), pointer :: idt + integer(I4B) :: i + ! + ! -- initialize + read_as_arrays = .false. + ! + do i = 1, size(param_dfns) + idt => param_dfns(i) + if (param_dfns(i)%tagname == 'READASARRAYS') then + read_as_arrays = .true. + exit + end if + end do + ! + ! -- return + return + end function read_as_arrays + + subroutine add_block_var(ncpkg, ncid, varid, varname, ctype, sctype, & + mempath, modelfname, iout) + ! -- modules + use IdmDfnSelectorModule, only: param_definitions + ! -- dummy + type(NC4ModelPackageInputType), pointer, intent(in) :: ncpkg + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + character(len=*), intent(inout) :: varname + character(len=*), intent(in) :: ctype + character(len=*), intent(in) :: sctype + character(len=*), intent(in) :: mempath + character(len=*), intent(in) :: modelfname + integer(I4B), intent(in) :: iout + ! -- local + type(InputParamDefinitionType), dimension(:), & + pointer :: param_dfns + type(InputParamDefinitionType), pointer :: idt + integer(I4B) :: griddata + logical(LGP) :: array_input + ! + ! -- set dfns and input type + param_dfns => param_definitions(ctype, sctype) + array_input = read_as_arrays(param_dfns) + ! + ! -- assume tagname is unique to input dfn set + idt => package_scoped_param_dfn(param_dfns, ncpkg%component_type, & + sctype, varname, '') + ! + if (array_input) then + ! -- verify variable attribute mf6_iper is defined for dynamic array inputs + if (idt%blockname == 'PERIOD') then + if (nf90_inquire_attribute(ncid, varid, "mf6_iper", & + len=griddata) /= NF90_NOERR) then + errmsg = 'Required attribute mf6_iper not found for array & + &input param "'//trim(varname)//'".' + call store_error(errmsg) + call store_error_filename(modelfname) + end if + end if + end if + ! + ! -- add package variable + call ncpkg%add(idt%blockname, varname, varid) + write (iout, '(a)') 'NC4 input parameter "'//trim(varname)// & + '" discovered for package "'//trim(mempath)//'".' + ! + ! -- return + return + end subroutine add_block_var + + subroutine create_pkg_ipers(ncpkg, ncid, varid, varname, iout) + ! -- modules + ! -- dummy + type(NC4ModelPackageInputType), pointer, intent(in) :: ncpkg + integer(I4B), intent(in) :: ncid + integer(I4B), intent(in) :: varid + character(len=*), intent(inout) :: varname + integer(I4B), intent(in) :: iout + ! -- local + integer(I4B) :: vartype, ndims, nattr + integer(I4B), dimension(IDM_NETCDF4_MAX_DIM) :: dimids, dimlens + character(len=LINELENGTH) :: dimname + ! + ! -- inquire for IPER variable info + call nf_verify(nf90_inquire_variable(ncid, varid, varname, vartype, ndims, & + dimids, nattr), ncid, iout) + ! + ! -- allocate ipers array + if (ndims == 1) then + call nf_verify(nf90_inquire_dimension(ncid, dimids(1), dimname, & + dimlens(1)), ncid, iout) + ! -- allocate and set ipers array + allocate (ncpkg%ipers(dimlens(1))) + call nf_verify(nf90_get_var(ncid, varid, ncpkg%ipers), ncid, iout) + end if + ! + ! -- return + return + end subroutine create_pkg_ipers + + subroutine add_package_var(nc4_context, modelname, varid, iout) + ! -- modules + use MemoryHelperModule, only: split_mem_path + use SourceCommonModule, only: idm_subcomponent_type + ! -- dummy + type(NC4ModelInputsType), intent(inout) :: nc4_context + character(len=*), intent(in) :: modelname + integer(I4B), intent(in) :: varid + integer(I4B), intent(in) :: iout + ! -- local + type(NC4ModelPackageInputType), pointer :: ncpkg + character(len=LENCOMPONENTNAME) :: component, scname + character(len=LINELENGTH) :: input_str + character(len=LENPACKAGETYPE) :: sctype, pkgtype + character(len=LINELENGTH) :: mem_address + character(len=LENMEMPATH) :: mempath + character(len=LINELENGTH) :: varname + integer(I4B) :: ncid, idx + ! + ! -- initialize + ncid = nc4_context%ncid + varname = '' + ! + if (nf90_get_att(ncid, varid, "mf6_input", & + input_str) == NF90_NOERR) then + ! + idx = index(input_str, ':') + ! + if (idx > 0) then + ! -- split input string into package type and path tokens + pkgtype = input_str(1:idx - 1) + ! -- subcomponent type + sctype = idm_subcomponent_type(nc4_context%modeltype, pkgtype) + ! -- not a true address- tagname instead of mf6varname + mem_address = input_str(idx + 1:len_trim(input_str)) + else + ! TODO: invalid mf6_input string error + end if + ! + idx = index(mem_address, '/', back=.true.) + ! + if (idx > 0) then + varname = mem_address(idx + 1:len_trim(mem_address)) + ! -- split address into mempath and input tag + mempath = mem_address(1:idx - 1) + call split_mem_path(mempath, component, scname) + ! + ncpkg => get_nc4_package(nc4_context, component, scname, & + pkgtype, sctype, iout) + ! + if (varname == 'IPER') then + ! -- allocate and set package ipers + call create_pkg_ipers(ncpkg, ncid, varid, varname, iout) + ! + else + ! -- track package parameter + call add_block_var(ncpkg, ncid, varid, varname, & + nc4_context%component_type, sctype, & + mempath, nc4_context%modelfname, iout) + end if + else + ! TODO: invalid mf6_input string + end if + end if + ! + ! -- return + return + end subroutine add_package_var + + subroutine nc4_pkg_context(nc4_context, iout) + use InputOutputModule, only: lowcase, upcase + type(NC4ModelInputsType), intent(inout) :: nc4_context + integer(I4B), intent(in) :: iout + integer(I4B) :: ndim, nvar, nattr, unlimDimID, format + integer(I4B), dimension(:), allocatable :: varids + ! -- local + character(len=LENPACKAGETYPE) :: modeltype + character(len=LENCOMPONENTNAME) :: modelname + integer(I4B) :: iparam, ncid + ! + ncid = nc4_context%ncid + ! + ! -- verify expected mf6_modeltype file attribute + if (nf90_get_att(ncid, NF90_GLOBAL, "mf6_modeltype", & + modeltype) == NF90_NOERR) then + ! + call upcase(modeltype) + ! -- verify modeltype matches this model + if (nc4_context%modeltype /= modeltype) then + errmsg = 'NC4 input file model type does not match namefile model type' + call store_error(errmsg) + call store_error_filename(nc4_context%modelfname) + end if + ! + else + errmsg = 'NC4 input file global attribute "mf6_modeltype" not found.' + call store_error(errmsg) + call store_error_filename(nc4_context%modelfname) + ! + end if + ! + ! -- verify expected mf6_modelname file attribute + if (nf90_get_att(ncid, NF90_GLOBAL, "mf6_modelname", & + modelname) == NF90_NOERR) then + ! + call upcase(modelname) + ! -- verify modelname matches this model + if (nc4_context%modelname /= modelname) then + errmsg = 'NC4 input file model name does not match namefile model name' + call store_error(errmsg) + call store_error_filename(nc4_context%modelfname) + end if + ! + else + errmsg = 'NC4 input file global attribute "mf6_modelname" not found.' + call store_error(errmsg) + call store_error_filename(nc4_context%modelfname) + ! + end if + ! + ! -- inquire for root dataset info + call nf_verify(nf90_inquire(ncid, ndim, nvar, nattr, unlimdimid, format), & + ncid, iout) + ! + ! -- allocate and set varids + allocate (varids(nvar)) + call nf_verify(nf90_inq_varids(ncid, nvar, varids), ncid, iout) + ! + ! -- identify package variables; build block variable lists + do iparam = 1, nvar + ! + call add_package_var(nc4_context, modelname, varids(iparam), iout) + ! + end do + ! + ! -- cleanup + deallocate (varids) + ! + ! -- return + return + end subroutine nc4_pkg_context + +end module LoadNetCDFDataModule diff --git a/src/Utilities/Idm/netcdf/NC4ModelInputs.f90 b/src/Utilities/Idm/netcdf/NC4ModelInputs.f90 new file mode 100644 index 00000000000..4a7472407c8 --- /dev/null +++ b/src/Utilities/Idm/netcdf/NC4ModelInputs.f90 @@ -0,0 +1,426 @@ +!> @brief This module contains the NC4ModelInputsModule +!! +!! +!< +module NC4ModelInputsModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: LINELENGTH, LENCOMPONENTNAME, & + LENVARNAME, LENPACKAGETYPE, & + LENPACKAGENAME + use SimModule, only: store_error, store_error_filename + use SimVariablesModule, only: errmsg + use ListModule, only: ListType + use ArrayHandlersModule, only: expandarray + use InputDefinitionModule, only: InputParamDefinitionType, & + InputBlockDefinitionType + use IdmDfnSelectorModule, only: param_definitions, & + block_definitions + + implicit none + private + public :: NC4ModelPackageInputType + public :: NC4ModelInputsType + public :: NC4BlockVariablesType + + type :: NC4BlockVariablesType + ! -- block + character(len=LENCOMPONENTNAME) :: blockname + integer(I4B) :: idfn + logical(LGP) :: aggregate + ! -- input tagname; variable tags can exceed LENVARNAME + character(len=LINELENGTH), dimension(:), allocatable :: varnames + ! -- model nc4 dataset variable id + integer(I4B), dimension(:), allocatable :: varids + ! -- number of variables + integer(I4B) :: varnum + contains + procedure :: create => blockvars_create + procedure :: add => blockvars_add + procedure :: destroy => blockvars_destroy + end type NC4BlockVariablesType + + type :: NC4ModelPackageInputType + type(NC4BlockVariablesType), dimension(:), allocatable :: blocklist + integer(I4B), dimension(:), allocatable :: ipers + character(len=LENCOMPONENTNAME) :: component_type + character(len=LENCOMPONENTNAME) :: subcomponent_type + character(len=LENCOMPONENTNAME) :: component_name + character(len=LENCOMPONENTNAME) :: subcomponent_name + contains + procedure :: init => nc4pkg_init + procedure :: reset => nc4pkg_reset + procedure :: destroy => nc4pkg_destroy + procedure :: add => nc4pkg_add + end type NC4ModelPackageInputType + + type :: NC4ModelInputsType + character(len=LENPACKAGETYPE) :: modeltype + character(len=LENCOMPONENTNAME) :: component_type + character(len=LENCOMPONENTNAME) :: modelname + character(len=LINELENGTH) :: modelfname + integer(I4B) :: ncid + type(ListType) :: pkglist + contains + procedure :: init => nc4modelinputs_init + procedure :: add => nc4modelinputs_add + procedure :: get => nc4modelinputs_get + procedure :: get_package => nc4modelinputs_get_package + procedure :: destroy => nc4modelinputs_destroy + end type NC4ModelInputsType + +contains + + !> @brief create a new package type + !< + subroutine blockvars_create(this, blockname, idfn, aggregate) + ! -- modules + ! -- dummy + class(NC4BlockVariablesType) :: this + character(len=*), intent(in) :: blockname + integer(I4B), intent(in) :: idfn + logical(LGP), intent(in) :: aggregate + ! -- local + ! + ! -- initialize + this%blockname = blockname + this%idfn = idfn + this%aggregate = aggregate + this%varnum = 0 + ! + ! -- allocate arrays + allocate (this%varnames(0)) + allocate (this%varids(0)) + ! + ! -- return + return + end subroutine blockvars_create + + !> @brief add a new package instance to this package type + !< + subroutine blockvars_add(this, varname, varid) + ! -- modules + ! -- dummy + class(NC4BlockVariablesType) :: this + character(len=*), intent(in) :: varname + integer(I4B), intent(in) :: varid + ! -- local + ! + ! -- reallocate + call expandarray(this%varnames) + call expandarray(this%varids) + ! + ! -- add new package instance + this%varnum = this%varnum + 1 + this%varnames(this%varnum) = varname + this%varids(this%varnum) = varid + ! + ! -- return + return + end subroutine blockvars_add + + !> @brief deallocate object + !< + subroutine blockvars_destroy(this) + ! -- modules + ! -- dummy + class(NC4BlockVariablesType) :: this + ! -- local + ! + ! -- deallocate dynamic arrays + if (allocated(this%varnames)) deallocate (this%varnames) + if (allocated(this%varids)) deallocate (this%varids) + ! + ! -- return + return + end subroutine blockvars_destroy + + !> @brief initialize model package inputs object + !< + subroutine nc4pkg_init(this, pkgtype, component_type, subcomponent_type, & + component_name, subcomponent_name, filename) + ! -- modules + ! -- dummy + class(NC4ModelPackageInputType) :: this + character(len=*) :: pkgtype + character(len=*) :: component_type + character(len=*) :: subcomponent_type + character(len=*) :: component_name + character(len=*) :: subcomponent_name + character(len=*) :: filename + ! -- local + character(len=100) :: blockname + type(InputBlockDefinitionType), dimension(:), & + pointer :: block_dfns + integer(I4B) :: iblock + ! + ! -- initialize object + this%component_type = component_type + this%subcomponent_type = subcomponent_type + this%component_name = component_name + this%subcomponent_name = subcomponent_name + ! + ! -- set pointer to block definition list + block_dfns => block_definitions(component_type, subcomponent_type) + ! + ! allocate blocklist + allocate (this%blocklist(size(block_dfns))) + ! + ! -- create blocklist types + do iblock = 1, size(block_dfns) + blockname = block_dfns(iblock)%blockname + call this%blocklist(iblock)%create(blockname, iblock, & + block_dfns(iblock)%aggregate) + end do + ! + ! -- return + return + end subroutine nc4pkg_init + + !> @brief reset model nc4 package blocklist + !! + !! This routine resets the blocklist. Needed for array based inputs + !! because the model namefile packages block will specify RCH6, for + !! example, in reference to either the array or list based definition + !! set. This routine is called when the definition set should be array + !! but was initialized as list. It should only be called once per + !! package when the first input variable is read. + !! + !< + subroutine nc4pkg_reset(this, pkgtype, subcomponent_type) + ! -- modules + ! -- dummy + class(NC4ModelPackageInputType) :: this + character(len=*) :: pkgtype + character(len=*) :: subcomponent_type + ! -- local + character(len=100) :: blockname + type(InputBlockDefinitionType), dimension(:), & + pointer :: block_dfns + integer(I4B) :: n, iblock + ! + ! -- destroy blocklist + do n = 1, size(this%blocklist) + call this%blocklist(n)%destroy() + end do + ! + if (allocated(this%blocklist)) deallocate (this%blocklist) + ! + ! -- reset subcomponent_type + this%subcomponent_type = subcomponent_type + ! + ! -- set pointer to block definition list + block_dfns => block_definitions(this%component_type, this%subcomponent_type) + ! + ! allocate blocklist + allocate (this%blocklist(size(block_dfns))) + ! + ! -- create blocklist types + do iblock = 1, size(block_dfns) + blockname = block_dfns(iblock)%blockname + call this%blocklist(iblock)%create(blockname, iblock, & + block_dfns(iblock)%aggregate) + end do + ! + ! -- return + return + end subroutine nc4pkg_reset + + !> @brief add package variable info to package type list + !< + subroutine nc4pkg_add(this, blockname, varname, varid) + ! -- modules + use InputOutputModule, only: lowcase, upcase + ! -- dummy + class(NC4ModelPackageInputType) :: this + character(len=*), intent(inout) :: blockname + character(len=*), intent(inout) :: varname + integer(I4B) :: varid + ! -- local + integer(I4B) :: n + ! + ! -- set to upper case + call upcase(blockname) + call upcase(varname) + ! + ! -- locate index of pkgtype in pkglist + do n = 1, size(this%blocklist) + if (this%blocklist(n)%blockname == blockname) then + call this%blocklist(n)%add(varname, varid) + exit + end if + end do + ! + ! -- return + return + end subroutine nc4pkg_add + + !> @brief deallocate object + !< + subroutine nc4pkg_destroy(this) + ! -- modules + ! -- dummy + class(NC4ModelPackageInputType) :: this + ! -- local + integer(I4B) :: n + ! + ! -- + do n = 1, size(this%blocklist) + call this%blocklist(n)%destroy() + end do + ! + if (allocated(this%ipers)) deallocate (this%ipers) + if (allocated(this%blocklist)) deallocate (this%blocklist) + ! + ! -- return + return + end subroutine nc4pkg_destroy + + !> @brief init object + !< + subroutine nc4modelinputs_init(this, modeltype, component_type, & + modelname, modelfname, ncid) + ! -- modules + ! -- dummy + class(NC4ModelInputsType) :: this + character(len=*) :: modeltype + character(len=*) :: component_type + character(len=*) :: modelname + character(len=*) :: modelfname + integer(I4B) :: ncid + ! -- local + ! + ! -- initialize object + this%modeltype = modeltype + this%component_type = component_type + this%modelname = modelname + this%modelfname = modelfname + this%ncid = ncid + ! + ! -- return + return + end subroutine nc4modelinputs_init + + !> @brief add to package list + !< + subroutine nc4modelinputs_add(this, pkgtype, pkgname, sc_type) + ! -- dummy variables + class(NC4ModelInputsType) :: this + character(len=*) :: pkgtype + character(len=*) :: pkgname + character(len=*) :: sc_type + class(NC4ModelPackageInputType), pointer :: nc4pkg + ! -- local variables + class(*), pointer :: obj + ! + ! -- + allocate (nc4pkg) + call nc4pkg%init(pkgtype, this%component_type, & + sc_type, this%modelname, & + pkgname, this%modelfname) + ! + obj => nc4pkg + call this%pkglist%Add(obj) + ! + ! -- return + return + end subroutine nc4modelinputs_add + + !> @brief get package list + !< + function nc4modelinputs_get(this, idx) result(nc4pkg) + ! -- dummy variables + class(NC4ModelInputsType) :: this + integer(I4B), intent(in) :: idx + ! -- return + class(NC4ModelPackageInputType), pointer :: nc4pkg + ! -- local variables + class(*), pointer :: obj + ! + ! -- initialize + nc4pkg => null() + ! + obj => this%pkglist%GetItem(idx) + if (associated(obj)) then + select type (obj) + class is (NC4ModelPackageInputType) + nc4pkg => obj + end select + end if + ! + ! -- return + return + end function nc4modelinputs_get + + !> @brief get package list + !< + subroutine nc4modelinputs_destroy(this) + ! -- dummy variables + class(NC4ModelInputsType) :: this + ! -- local variables + class(NC4ModelPackageInputType), pointer :: nc4pkg + integer(I4B) :: n + ! + ! -- initialize + nc4pkg => null() + ! + do n = 1, this%pkglist%count() + nc4pkg => this%get(n) + call nc4pkg%destroy() + deallocate (nc4pkg) + nullify (nc4pkg) + end do + ! + call this%pkglist%clear() + ! + ! -- return + return + end subroutine nc4modelinputs_destroy + + !> @brief get package list + !< + function nc4modelinputs_get_package(this, component_name, subcomponent_name) & + result(nc4_pkg) + ! -- dummy variables + class(NC4ModelInputsType) :: this + character(len=*), intent(in) :: component_name + character(len=*), intent(in) :: subcomponent_name + ! -- return + class(NC4ModelPackageInputType), pointer :: nc4_pkg + ! -- local variables + class(*), pointer :: obj + class(NC4ModelPackageInputType), pointer :: pkg + integer(I4B) :: n + ! + ! -- initialize + nc4_pkg => null() + ! + do n = 1, this%pkglist%count() + obj => this%pkglist%GetItem(n) + if (associated(obj)) then + select type (obj) + class is (NC4ModelPackageInputType) + pkg => obj + if (component_name == pkg%component_name .and. & + subcomponent_name == pkg%subcomponent_name) then + nc4_pkg => obj + exit + end if + end select + end if + end do + ! + ! -- set error if not found + if (.not. associated(nc4_pkg)) then + errmsg = 'Programming error. NC4 Model package not found. & + &Model='//trim(component_name)//', package='// & + trim(subcomponent_name)//'.' + call store_error(errmsg) + call store_error_filename(this%modelfname) + end if + ! + ! -- return + return + end function nc4modelinputs_get_package + +end module NC4ModelInputsModule diff --git a/src/meson.build b/src/meson.build index c5ee53ab891..9e780534dec 100644 --- a/src/meson.build +++ b/src/meson.build @@ -170,6 +170,7 @@ modflow_sources = files( 'Utilities' / 'Idm' / 'ModelPackageInputs.f90', 'Utilities' / 'Idm' / 'ModflowInput.f90', 'Utilities' / 'Idm' / 'SourceCommon.f90', + 'Utilities' / 'Idm' / 'SourceContext.f90', 'Utilities' / 'Idm' / 'SourceLoad.F90', 'Utilities' / 'Idm' / 'mf6blockfile' / 'AsciiInputLoadType.f90', 'Utilities' / 'Idm' / 'mf6blockfile' / 'IdmMf6File.f90', @@ -178,6 +179,7 @@ modflow_sources = files( 'Utilities' / 'Idm' / 'mf6blockfile' / 'StressListInput.f90', 'Utilities' / 'Idm' / 'mf6blockfile' / 'StructArray.f90', 'Utilities' / 'Idm' / 'mf6blockfile' / 'StructVector.f90', + 'Utilities' / 'Idm' / 'netcdf' / 'NC4ModelInputs.f90', 'Utilities' / 'Idm' / 'selector' / 'IdmDfnSelector.f90', 'Utilities' / 'Idm' / 'selector' / 'IdmGwfDfnSelector.f90', 'Utilities' / 'Idm' / 'selector' / 'IdmGwtDfnSelector.f90', @@ -272,6 +274,11 @@ modflow_mpi_sources = files( 'Solution' / 'ParallelSolution.f90' ) +modflow_netcdf_sources = files( + 'Utilities' / 'Idm' / 'netcdf' / 'IdmNetCDFFile.f90', + 'Utilities' / 'Idm' / 'netcdf' / 'LoadNetCDFData.f90', +) + if with_petsc modflow_sources += modflow_petsc_sources endif @@ -280,6 +287,10 @@ if with_mpi modflow_sources += modflow_mpi_sources endif +if with_netcdf + modflow_sources += modflow_netcdf_sources +endif + mf6_external = static_library('mf6_external', external_libraries) message('MODFLOW 6 executable name: ' + buildname)