From c50ec80f9a7c402b770990941b520cf27559ea0c Mon Sep 17 00:00:00 2001 From: Eric Morway Date: Sat, 7 Oct 2023 04:09:38 -0700 Subject: [PATCH] refactor(tsp): Elevate FMI to generalized transport class (#1386) * refactor(tsp): Elevate FMI to generalized transport class * fprettify * Some missing ADV code as a result of not elevating it with this PR * more spruce up in ADV * Missing a deallocation in GwtInterfaceModel * Thought I ran build_makefiles.py, but maybe not --- make/makefile | 52 +- msvs/mf6core.vfproj | 7 +- src/Model/Connection/GwtInterfaceModel.f90 | 14 +- src/Model/GroundWaterTransport/gwt1.f90 | 268 +++------ src/Model/GroundWaterTransport/gwt1adv1.f90 | 24 +- src/Model/GroundWaterTransport/gwt1apt1.f90 | 4 +- src/Model/GroundWaterTransport/gwt1dsp1.f90 | 6 +- src/Model/GroundWaterTransport/gwt1ist1.f90 | 6 +- src/Model/GroundWaterTransport/gwt1lkt1.f90 | 4 +- src/Model/GroundWaterTransport/gwt1mst1.f90 | 6 +- src/Model/GroundWaterTransport/gwt1mvt1.f90 | 22 +- src/Model/GroundWaterTransport/gwt1mwt1.f90 | 4 +- src/Model/GroundWaterTransport/gwt1sft1.f90 | 4 +- src/Model/GroundWaterTransport/gwt1ssm1.f90 | 6 +- src/Model/GroundWaterTransport/gwt1uzt1.f90 | 4 +- .../ModelUtilities/FlowModelInterface.f90 | 20 +- src/Model/TransportModel/tsp1.f90 | 529 +++++++++++++++++- .../tsp1fmi1.f90} | 147 +++-- src/meson.build | 2 +- 19 files changed, 806 insertions(+), 323 deletions(-) rename src/Model/{GroundWaterTransport/gwt1fmi1.f90 => TransportModel/tsp1fmi1.f90} (90%) diff --git a/make/makefile b/make/makefile index 2743fb98deb..adffc5fe84e 100644 --- a/make/makefile +++ b/make/makefile @@ -1,40 +1,40 @@ -# makefile created by pymake (version 1.2.9.dev0) for the 'mf6' executable. +# makefile created by pymake (version 1.2.7) for the 'mf6' executable. include ./makedefaults # Define the source file directories SOURCEDIR1=../src -SOURCEDIR2=../src/Exchange -SOURCEDIR3=../src/Model -SOURCEDIR4=../src/Model/Geometry -SOURCEDIR5=../src/Model/TransportModel -SOURCEDIR6=../src/Model/ModelUtilities -SOURCEDIR7=../src/Model/Connection +SOURCEDIR2=../src/Distributed +SOURCEDIR3=../src/Exchange +SOURCEDIR4=../src/Model +SOURCEDIR5=../src/Model/Connection +SOURCEDIR6=../src/Model/Geometry +SOURCEDIR7=../src/Model/GroundWaterFlow SOURCEDIR8=../src/Model/GroundWaterTransport -SOURCEDIR9=../src/Model/GroundWaterFlow -SOURCEDIR10=../src/Distributed +SOURCEDIR9=../src/Model/ModelUtilities +SOURCEDIR10=../src/Model/TransportModel SOURCEDIR11=../src/Solution -SOURCEDIR12=../src/Solution/PETSc -SOURCEDIR13=../src/Solution/LinearMethods +SOURCEDIR12=../src/Solution/LinearMethods +SOURCEDIR13=../src/Solution/PETSc SOURCEDIR14=../src/Timing SOURCEDIR15=../src/Utilities -SOURCEDIR16=../src/Utilities/TimeSeries -SOURCEDIR17=../src/Utilities/Libraries -SOURCEDIR18=../src/Utilities/Libraries/rcm -SOURCEDIR19=../src/Utilities/Libraries/sparsekit -SOURCEDIR20=../src/Utilities/Libraries/sparskit2 +SOURCEDIR16=../src/Utilities/ArrayRead +SOURCEDIR17=../src/Utilities/Idm +SOURCEDIR18=../src/Utilities/Idm/mf6blockfile +SOURCEDIR19=../src/Utilities/Idm/selector +SOURCEDIR20=../src/Utilities/Libraries SOURCEDIR21=../src/Utilities/Libraries/blas SOURCEDIR22=../src/Utilities/Libraries/daglib -SOURCEDIR23=../src/Utilities/Idm -SOURCEDIR24=../src/Utilities/Idm/selector -SOURCEDIR25=../src/Utilities/Idm/mf6blockfile +SOURCEDIR23=../src/Utilities/Libraries/rcm +SOURCEDIR24=../src/Utilities/Libraries/sparsekit +SOURCEDIR25=../src/Utilities/Libraries/sparskit2 SOURCEDIR26=../src/Utilities/Matrix -SOURCEDIR27=../src/Utilities/Vector +SOURCEDIR27=../src/Utilities/Memory SOURCEDIR28=../src/Utilities/Observation SOURCEDIR29=../src/Utilities/OutputControl -SOURCEDIR30=../src/Utilities/Memory -SOURCEDIR31=../src/Utilities/ArrayRead +SOURCEDIR30=../src/Utilities/TimeSeries +SOURCEDIR31=../src/Utilities/Vector VPATH = \ ${SOURCEDIR1} \ @@ -207,7 +207,7 @@ $(OBJDIR)/IndexMap.o \ $(OBJDIR)/VirtualModel.o \ $(OBJDIR)/BaseExchange.o \ $(OBJDIR)/UzfCellGroup.o \ -$(OBJDIR)/gwt1fmi1.o \ +$(OBJDIR)/tsp1fmi1.o \ $(OBJDIR)/OutputControlData.o \ $(OBJDIR)/gwf3ic8.o \ $(OBJDIR)/Xt3dInterface.o \ @@ -220,6 +220,9 @@ $(OBJDIR)/ImsLinearSettings.o \ $(OBJDIR)/ConvergenceSummary.o \ $(OBJDIR)/CellWithNbrs.o \ $(OBJDIR)/NumericalExchange.o \ +$(OBJDIR)/gwf3disv8.o \ +$(OBJDIR)/gwf3disu8.o \ +$(OBJDIR)/gwf3dis8.o \ $(OBJDIR)/gwf3uzf8.o \ $(OBJDIR)/gwt1apt1.o \ $(OBJDIR)/GwtSpc.o \ @@ -240,7 +243,6 @@ $(OBJDIR)/SparseMatrix.o \ $(OBJDIR)/LinearSolverBase.o \ $(OBJDIR)/ims8reordering.o \ $(OBJDIR)/VirtualExchange.o \ -$(OBJDIR)/gwf3disu8.o \ $(OBJDIR)/GridSorting.o \ $(OBJDIR)/DisConnExchange.o \ $(OBJDIR)/CsrUtils.o \ @@ -258,8 +260,6 @@ $(OBJDIR)/gwt1ist1.o \ $(OBJDIR)/gwt1dsp1.o \ $(OBJDIR)/gwt1cnc1.o \ $(OBJDIR)/gwt1adv1.o \ -$(OBJDIR)/gwf3disv8.o \ -$(OBJDIR)/gwf3dis8.o \ $(OBJDIR)/gwf3api8.o \ $(OBJDIR)/gwf3wel8.o \ $(OBJDIR)/gwf3rch8.o \ diff --git a/msvs/mf6core.vfproj b/msvs/mf6core.vfproj index 0949f90e7d7..1a428652844 100644 --- a/msvs/mf6core.vfproj +++ b/msvs/mf6core.vfproj @@ -171,7 +171,6 @@ - @@ -184,8 +183,7 @@ - - + @@ -209,7 +207,8 @@ - + + diff --git a/src/Model/Connection/GwtInterfaceModel.f90 b/src/Model/Connection/GwtInterfaceModel.f90 index 81d08b1a064..04fab6facfe 100644 --- a/src/Model/Connection/GwtInterfaceModel.f90 +++ b/src/Model/Connection/GwtInterfaceModel.f90 @@ -1,11 +1,12 @@ module GwtInterfaceModelModule use KindModule, only: I4B, DP + use ConstantsModule, only: DONE use MemoryManagerModule, only: mem_allocate, mem_deallocate, mem_reallocate use MemoryHelperModule, only: create_mem_path use NumericalModelModule, only: NumericalModelType use GwtModule, only: GwtModelType, CastAsGwtModel use GwfDisuModule, only: disu_cr, CastAsDisuType - use GwtFmiModule, only: fmi_cr, GwtFmiType + use TspFmiModule, only: fmi_cr, TspFmiType use GwtAdvModule, only: adv_cr, GwtAdvType use GwtAdvOptionsModule, only: GwtAdvOptionsType use GwtDspModule, only: dsp_cr, GwtDspType @@ -25,6 +26,7 @@ module GwtInterfaceModelModule integer(i4B), pointer :: iAdvScheme => null() !< the advection scheme: 0 = up, 1 = central, 2 = tvd integer(i4B), pointer :: ixt3d => null() !< xt3d setting: 0 = off, 1 = lhs, 2 = rhs + real(DP), pointer :: ieqnsclfac => null() !< governing eqn scaling factor: 1: GWT, >1: GWE class(GridConnectionType), pointer :: gridConnection => null() !< The grid connection class will provide the interface grid class(GwtModelType), private, pointer :: owner => null() !< the real GWT model for which the exchange coefficients @@ -59,6 +61,7 @@ subroutine gwtifmod_cr(this, name, iout, gridConn) ! defaults this%iAdvScheme = 0 this%ixt3d = 0 + this%ieqnsclfac = DONE this%iout = iout this%gridConnection => gridConn @@ -79,8 +82,10 @@ subroutine gwtifmod_cr(this, name, iout, gridConn) ! create dis and packages call disu_cr(this%dis, this%name, '', -1, this%iout) - call fmi_cr(this%fmi, this%name, 0, this%iout) - call adv_cr(this%adv, this%name, adv_unit, this%iout, this%fmi) + call fmi_cr(this%fmi, this%name, 0, this%iout, this%ieqnsclfac, & + this%depvartype) + call adv_cr(this%adv, this%name, adv_unit, this%iout, this%fmi, & + this%ieqnsclfac) call dsp_cr(this%dsp, this%name, '', -dsp_unit, this%iout, this%fmi) call gwt_obs_cr(this%obs, inobs) @@ -94,6 +99,7 @@ subroutine allocate_scalars(this, modelname) call mem_allocate(this%iAdvScheme, 'ADVSCHEME', this%memoryPath) call mem_allocate(this%ixt3d, 'IXT3D', this%memoryPath) + call mem_allocate(this%ieqnsclfac, 'IEQNSCLFAC', this%memoryPath) end subroutine allocate_scalars @@ -192,6 +198,7 @@ subroutine gwtifmod_da(this) ! this call mem_deallocate(this%iAdvScheme) call mem_deallocate(this%ixt3d) + call mem_deallocate(this%ieqnsclfac) ! gwt packages call this%dis%dis_da() @@ -219,6 +226,7 @@ subroutine gwtifmod_da(this) call mem_deallocate(this%inmvt) call mem_deallocate(this%inoc) call mem_deallocate(this%inobs) + call mem_deallocate(this%eqnsclfac) ! base call this%NumericalModelType%model_da() diff --git a/src/Model/GroundWaterTransport/gwt1.f90 b/src/Model/GroundWaterTransport/gwt1.f90 index 29bc7bbadce..6c34835ad4a 100644 --- a/src/Model/GroundWaterTransport/gwt1.f90 +++ b/src/Model/GroundWaterTransport/gwt1.f90 @@ -8,14 +8,15 @@ module GwtModule use KindModule, only: DP, I4B - use ConstantsModule, only: LENFTYPE, LENMEMPATH, DZERO, LENPAKLOC + use ConstantsModule, only: LENFTYPE, LENMEMPATH, DZERO, DONE, & + LENPAKLOC, LENVARNAME use VersionModule, only: write_listfile_header use NumericalModelModule, only: NumericalModelType use TransportModelModule, only: TransportModelType use BaseModelModule, only: BaseModelType use BndModule, only: BndType, AddBndToList, GetBndFromList use GwtIcModule, only: GwtIcType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use GwtAdvModule, only: GwtAdvType use GwtDspModule, only: GwtDspType use GwtSsmModule, only: GwtSsmType @@ -32,11 +33,13 @@ module GwtModule public :: gwt_cr public :: GwtModelType public :: CastAsGwtModel + character(len=LENVARNAME), parameter :: dvt = 'CONCENTRATION ' !< dependent variable type, varies based on model type + character(len=LENVARNAME), parameter :: dvu = 'MASS ' !< dependent variable unit of measure, either "mass" or "energy" + character(len=LENVARNAME), parameter :: dvua = 'M ' !< abbreviation of the dependent variable unit of measure, either "M" or "E" type, extends(TransportModelType) :: GwtModelType type(GwtIcType), pointer :: ic => null() ! initial conditions package - type(GwtFmiType), pointer :: fmi => null() ! flow model interface type(GwtMstType), pointer :: mst => null() ! mass storage and transfer package type(GwtAdvType), pointer :: adv => null() ! advection package type(GwtDspType), pointer :: dsp => null() ! dispersion package @@ -44,9 +47,7 @@ module GwtModule type(GwtMvtType), pointer :: mvt => null() ! mover transport package type(GwtOcType), pointer :: oc => null() ! output control package type(GwtObsType), pointer :: obs => null() ! observation package - type(BudgetType), pointer :: budget => null() ! budget object integer(I4B), pointer :: inic => null() ! unit number IC - integer(I4B), pointer :: infmi => null() ! unit number FMI integer(I4B), pointer :: inmvt => null() ! unit number MVT integer(I4B), pointer :: inmst => null() ! unit number MST integer(I4B), pointer :: inadv => null() ! unit number ADV @@ -72,18 +73,16 @@ module GwtModule procedure :: model_da => gwt_da procedure :: model_bdentry => gwt_bdentry procedure :: allocate_scalars - procedure, private :: package_create - procedure, private :: ftype_check procedure :: get_iasym => gwt_get_iasym procedure, private :: gwt_ot_flow procedure, private :: gwt_ot_flowja procedure, private :: gwt_ot_dv procedure, private :: gwt_ot_bdsummary procedure, private :: gwt_ot_obs - procedure, private :: create_packages + procedure :: create_packages => create_gwt_packages procedure, private :: create_bndpkgs - procedure, private :: create_lstfile - procedure, private :: log_namfile_options + procedure, private :: package_create + end type GwtModelType contains @@ -104,6 +103,7 @@ subroutine gwt_cr(filename, id, modelname) integer(I4B), intent(in) :: id character(len=*), intent(in) :: modelname ! -- local + integer(I4B) :: indis type(GwtModelType), pointer :: this class(BaseModelType), pointer :: model character(len=LENMEMPATH) :: input_mempath @@ -118,6 +118,10 @@ subroutine gwt_cr(filename, id, modelname) ! ! -- Allocate scalars and add model to basemodellist call this%allocate_scalars(modelname) + ! + ! -- set labels for transport model - needed by create_packages() below + call this%set_tsp_labels(this%macronym, dvt, dvu, dvua) + ! model => this call AddBaseModelToList(basemodellist, model) ! @@ -138,24 +142,19 @@ subroutine gwt_cr(filename, id, modelname) found%print_flows) call mem_set_value(this%ipakcb, 'SAVE_FLOWS', input_mempath, found%save_flows) ! - ! -- create the list file - call this%create_lstfile(lst_fname, filename, found%list) - ! ! -- activate save_flows if found if (found%save_flows) then this%ipakcb = -1 end if ! - ! -- log set options - if (this%iout > 0) then - call this%log_namfile_options(found) - end if - ! ! -- Create utility objects call budget_cr(this%budget, this%name) ! + ! -- Call parent class routine + call this%tsp_cr(filename, id, modelname, indis) + ! ! -- create model packages - call this%create_packages() + call this%create_packages(indis) ! ! -- return return @@ -165,7 +164,6 @@ end subroutine gwt_cr ! ! (1) call df routines for each package ! (2) set variables and pointers - ! !< subroutine gwt_df(this) ! -- modules @@ -296,6 +294,15 @@ subroutine gwt_ar(this) if (this%inssm > 0) call this%ssm%ssm_ar(this%dis, this%ibound, this%x) if (this%inobs > 0) call this%obs%gwt_obs_ar(this%ic, this%x, this%flowja) ! + ! -- Set governing equation scale factor. Note that this scale factor + ! -- cannot be set arbitrarily. For solute transport, it must be set + ! -- to 1. Setting it to a different value will NOT automatically + ! -- scale all the terms of the governing equation correctly by that + ! -- value. This is because much of the coding in the associated + ! -- packages implicitly assumes the governing equation for solute + ! -- transport is scaled by 1. (effectively unscaled). + this%eqnsclfac = DONE + ! ! -- Call dis_ar to write binary grid file !call this%dis%dis_ar(this%npf%icelltype) ! @@ -658,7 +665,9 @@ subroutine gwt_ot_obs(this) call packobj%bnd_bd_obs() call packobj%bnd_ot_obs() end do - + ! + ! -- Return + return end subroutine gwt_ot_obs !> @brief Save flows @@ -691,7 +700,7 @@ subroutine gwt_ot_flow(this, icbcfl, ibudfl, icbcun) if (this%inmvt > 0) then call this%mvt%mvt_ot_saveflow(icbcfl, ibudfl) end if - + ! ! -- Print GWF flows ! no need to print flowja ! no need to print mst @@ -703,7 +712,7 @@ subroutine gwt_ot_flow(this, icbcfl, ibudfl, icbcun) packobj => GetBndFromList(this%bndlist, ip) call packobj%bnd_ot_model_flows(icbcfl=icbcfl, ibudfl=ibudfl, icbcun=0) end do - + ! ! -- Print advanced package flows do ip = 1, this%bndlist%Count() packobj => GetBndFromList(this%bndlist, ip) @@ -712,7 +721,9 @@ subroutine gwt_ot_flow(this, icbcfl, ibudfl, icbcun) if (this%inmvt > 0) then call this%mvt%mvt_ot_printflow(icbcfl, ibudfl) end if - + ! + ! -- Return + return end subroutine gwt_ot_flow !> @brief Write intercell flows @@ -756,16 +767,18 @@ subroutine gwt_ot_dv(this, idvsave, idvprint, ipflag) integer(I4B), intent(inout) :: ipflag class(BndType), pointer :: packobj integer(I4B) :: ip - + ! ! -- Print advanced package dependent variables do ip = 1, this%bndlist%Count() packobj => GetBndFromList(this%bndlist, ip) call packobj%bnd_ot_dv(idvsave, idvprint) end do - + ! ! -- save head and print head call this%oc%oc_ot(ipflag) - + ! + ! -- Return + return end subroutine gwt_ot_dv !> @brief Print budget summary @@ -777,28 +790,29 @@ subroutine gwt_ot_bdsummary(this, ibudfl, ipflag) integer(I4B), intent(inout) :: ipflag class(BndType), pointer :: packobj integer(I4B) :: ip - ! ! -- Package budget summary do ip = 1, this%bndlist%Count() packobj => GetBndFromList(this%bndlist, ip) call packobj%bnd_ot_bdsummary(kstp, kper, this%iout, ibudfl) end do - + ! ! -- mover budget summary if (this%inmvt > 0) then call this%mvt%mvt_ot_bdsummary(ibudfl) end if - + ! ! -- model budget summary if (ibudfl /= 0) then ipflag = 1 call this%budget%budget_ot(kstp, kper, this%iout) end if - + ! ! -- Write to budget csv call this%budget%writecsv(totim) - + ! + ! -- Return + return end subroutine gwt_ot_bdsummary !> @brief Deallocate @@ -821,11 +835,10 @@ subroutine gwt_da(this) ! -- Internal flow packages deallocate call this%dis%dis_da() call this%ic%ic_da() - call this%fmi%fmi_da() - call this%adv%adv_da() call this%dsp%dsp_da() call this%ssm%ssm_da() call this%mst%mst_da() + call this%adv%adv_da() call this%mvt%mvt_da() call this%budget%budget_da() call this%oc%oc_da() @@ -834,11 +847,10 @@ subroutine gwt_da(this) ! -- Internal package objects deallocate (this%dis) deallocate (this%ic) - deallocate (this%fmi) - deallocate (this%adv) deallocate (this%dsp) deallocate (this%ssm) deallocate (this%mst) + deallocate (this%adv) deallocate (this%mvt) deallocate (this%budget) deallocate (this%oc) @@ -853,7 +865,6 @@ subroutine gwt_da(this) ! ! -- Scalars call mem_deallocate(this%inic) - call mem_deallocate(this%infmi) call mem_deallocate(this%inadv) call mem_deallocate(this%indsp) call mem_deallocate(this%inssm) @@ -862,10 +873,13 @@ subroutine gwt_da(this) call mem_deallocate(this%inoc) call mem_deallocate(this%inobs) ! + ! -- Parent class members + call this%TransportModelType%tsp_da() + ! ! -- NumericalModelType call this%NumericalModelType%model_da() ! - ! -- return + ! -- Return return end subroutine gwt_da @@ -889,7 +903,7 @@ subroutine gwt_bdentry(this, budterm, budtxt, rowlabel) ! call this%budget%addentry(budterm, delt, budtxt, rowlabel=rowlabel) ! - ! -- return + ! -- Return return end subroutine gwt_bdentry @@ -922,7 +936,7 @@ function gwt_get_iasym(this) result(iasym) if (packobj%iasym /= 0) iasym = 1 end do ! - ! -- return + ! -- Return return end function gwt_get_iasym @@ -935,31 +949,29 @@ subroutine allocate_scalars(this, modelname) class(GwtModelType) :: this character(len=*), intent(in) :: modelname ! - ! -- allocate members from parent class - call this%NumericalModelType%allocate_scalars(modelname) + ! -- allocate parent class scalars + call this%allocate_tsp_scalars(modelname) ! ! -- allocate members that are part of model class call mem_allocate(this%inic, 'INIC', this%memoryPath) - call mem_allocate(this%infmi, 'INFMI', this%memoryPath) + call mem_allocate(this%inadv, 'INADV', this%memoryPath) call mem_allocate(this%inmvt, 'INMVT', this%memoryPath) call mem_allocate(this%inmst, 'INMST', this%memoryPath) - call mem_allocate(this%inadv, 'INADV', this%memoryPath) call mem_allocate(this%indsp, 'INDSP', this%memoryPath) call mem_allocate(this%inssm, 'INSSM', this%memoryPath) call mem_allocate(this%inoc, 'INOC ', this%memoryPath) call mem_allocate(this%inobs, 'INOBS', this%memoryPath) ! this%inic = 0 - this%infmi = 0 + this%inadv = 0 this%inmvt = 0 this%inmst = 0 - this%inadv = 0 this%indsp = 0 this%inssm = 0 this%inoc = 0 this%inobs = 0 ! - ! -- return + ! -- Return return end subroutine allocate_scalars @@ -1033,61 +1045,23 @@ subroutine package_create(this, filtyp, ipakid, ipaknum, pakname, inunit, & end do call AddBndToList(this%bndlist, packobj) ! - ! -- return + ! -- Return return end subroutine package_create - !> @brief Make sure required input files have been specified - !< - subroutine ftype_check(this, indis) - ! -- modules - use ConstantsModule, only: LINELENGTH - use SimModule, only: store_error, count_errors, store_error_filename - ! -- dummy - class(GwtModelType) :: this - integer(I4B), intent(in) :: indis - ! -- local - character(len=LINELENGTH) :: errmsg - ! - ! -- Check for IC6, DIS(u), and MST. Stop if not present. - if (this%inic == 0) then - write (errmsg, '(a)') & - 'Initial conditions (IC6) package not specified.' - call store_error(errmsg) - end if - if (indis == 0) then - write (errmsg, '(a)') & - 'Discretization (DIS6 or DISU6) package not specified.' - call store_error(errmsg) - end if - if (this%inmst == 0) then - write (errmsg, '(a)') 'Mass storage and transfer (MST6) & - &package not specified.' - call store_error(errmsg) - end if - ! - if (count_errors() > 0) then - write (errmsg, '(a)') 'Required package(s) not specified.' - call store_error(errmsg) - call store_error_filename(this%filename) - end if - ! - ! -- return - return - end subroutine ftype_check - !> @brief Cast to GwtModelType function CastAsGwtModel(model) result(gwtmodel) class(*), pointer :: model !< The object to be cast class(GwtModelType), pointer :: gwtmodel !< The GWT model - + ! gwtmodel => null() if (.not. associated(model)) return select type (model) type is (GwtModelType) gwtmodel => model end select - + ! -- Return + return end function CastAsGwtModel !> @brief Source package info and begin to process @@ -1115,7 +1089,7 @@ subroutine create_bndpkgs(this, bndpkgs, pkgtypes, pkgnames, & character(len=LENMEMPATH) :: mempath integer(I4B), pointer :: inunit integer(I4B) :: n - + ! if (allocated(bndpkgs)) then ! ! -- create stress packages @@ -1143,13 +1117,13 @@ subroutine create_bndpkgs(this, bndpkgs, pkgtypes, pkgnames, & deallocate (bndpkgs) end if ! - ! -- return + ! -- Return return end subroutine create_bndpkgs !> @brief Source package info and begin to process !< - subroutine create_packages(this) + subroutine create_gwt_packages(this, indis) ! -- modules use ConstantsModule, only: LINELENGTH, LENPACKAGENAME use CharacterStringModule, only: CharacterStringType @@ -1157,11 +1131,7 @@ subroutine create_packages(this) use MemoryManagerModule, only: mem_setptr use MemoryHelperModule, only: create_mem_path use SimVariablesModule, only: idm_context - use GwfDisModule, only: dis_cr - use GwfDisvModule, only: disv_cr - use GwfDisuModule, only: disu_cr use GwtIcModule, only: ic_cr - use GwtFmiModule, only: fmi_cr use GwtMstModule, only: mst_cr use GwtAdvModule, only: adv_cr use GwtDspModule, only: dsp_cr @@ -1171,6 +1141,7 @@ subroutine create_packages(this) use GwtObsModule, only: gwt_obs_cr ! -- dummy class(GwtModelType) :: this + integer(I4B), intent(in) :: indis ! -- local type(CharacterStringType), dimension(:), contiguous, & pointer :: pkgtypes => null() @@ -1187,7 +1158,6 @@ subroutine create_packages(this) integer(I4B), pointer :: inunit integer(I4B), dimension(:), allocatable :: bndpkgs integer(I4B) :: n - integer(I4B) :: indis = 0 ! DIS enabled flag character(len=LENMEMPATH) :: mempathdsp = '' ! ! -- set input memory paths, input/model and input/model/namfile @@ -1209,19 +1179,8 @@ subroutine create_packages(this) ! ! -- create dis package first as it is a prerequisite for other packages select case (pkgtype) - case ('DIS6') - indis = 1 - call dis_cr(this%dis, this%name, mempath, indis, this%iout) - case ('DISV6') - indis = 1 - call disv_cr(this%dis, this%name, mempath, indis, this%iout) - case ('DISU6') - indis = 1 - call disu_cr(this%dis, this%name, mempath, indis, this%iout) case ('IC6') this%inic = inunit - case ('FMI6') - this%infmi = inunit case ('MVT6') this%inmvt = inunit case ('MST6') @@ -1248,9 +1207,9 @@ subroutine create_packages(this) ! ! -- Create packages that are tied directly to model call ic_cr(this%ic, this%name, this%inic, this%iout, this%dis) - call fmi_cr(this%fmi, this%name, this%infmi, this%iout) call mst_cr(this%mst, this%name, this%inmst, this%iout, this%fmi) - call adv_cr(this%adv, this%name, this%inadv, this%iout, this%fmi) + call adv_cr(this%adv, this%name, this%inadv, this%iout, this%fmi, & + this%eqnsclfac) call dsp_cr(this%dsp, this%name, mempathdsp, this%indsp, this%iout, & this%fmi) call ssm_cr(this%ssm, this%name, this%inssm, this%iout, this%fmi) @@ -1259,95 +1218,12 @@ subroutine create_packages(this) call gwt_obs_cr(this%obs, this%inobs) ! ! -- Check to make sure that required ftype's have been specified - call this%ftype_check(indis) + call this%ftype_check(indis, this%inmst, this%inic) ! call this%create_bndpkgs(bndpkgs, pkgtypes, pkgnames, mempaths, inunits) - - end subroutine create_packages - - subroutine create_lstfile(this, lst_fname, model_fname, defined) - ! -- modules - use KindModule, only: LGP - use InputOutputModule, only: openfile, getunit - ! -- dummy - class(GwtModelType) :: this - character(len=*), intent(inout) :: lst_fname - character(len=*), intent(in) :: model_fname - logical(LGP), intent(in) :: defined - ! -- local - integer(I4B) :: i, istart, istop ! - ! -- set list file name if not provided - if (.not. defined) then - ! - ! -- initialize - lst_fname = ' ' - istart = 0 - istop = len_trim(model_fname) - ! - ! -- identify '.' character position from back of string - do i = istop, 1, -1 - if (model_fname(i:i) == '.') then - istart = i - exit - end if - end do - ! - ! -- if not found start from string end - if (istart == 0) istart = istop + 1 - ! - ! -- set list file name - lst_fname = model_fname(1:istart) - istop = istart + 3 - lst_fname(istart:istop) = '.lst' - end if - ! - ! -- create the list file - this%iout = getunit() - call openfile(this%iout, 0, lst_fname, 'LIST', filstat_opt='REPLACE') - ! - ! -- write list file header - call write_listfile_header(this%iout, 'GROUNDWATER TRANSPORT MODEL (GWT)') - ! - ! -- return + ! -- Return return - end subroutine create_lstfile - - !> @brief Write model namfile options to list file - !< - subroutine log_namfile_options(this, found) - use GwfNamInputModule, only: GwfNamParamFoundType - class(GwtModelType) :: this - type(GwfNamParamFoundType), intent(in) :: found - - write (this%iout, '(1x,a)') 'NAMEFILE OPTIONS:' - - if (found%newton) then - write (this%iout, '(4x,a)') & - 'NEWTON-RAPHSON method enabled for the model.' - if (found%under_relaxation) then - write (this%iout, '(4x,a,a)') & - 'NEWTON-RAPHSON UNDER-RELAXATION based on the bottom ', & - 'elevation of the model will be applied to the model.' - end if - end if - - if (found%print_input) then - write (this%iout, '(4x,a)') 'STRESS PACKAGE INPUT WILL BE PRINTED '// & - 'FOR ALL MODEL STRESS PACKAGES' - end if - - if (found%print_flows) then - write (this%iout, '(4x,a)') 'PACKAGE FLOWS WILL BE PRINTED '// & - 'FOR ALL MODEL PACKAGES' - end if - - if (found%save_flows) then - write (this%iout, '(4x,a)') & - 'FLOWS WILL BE SAVED TO BUDGET FILE SPECIFIED IN OUTPUT CONTROL' - end if - - write (this%iout, '(1x,a)') 'END NAMEFILE OPTIONS:' - end subroutine log_namfile_options + end subroutine create_gwt_packages end module GwtModule diff --git a/src/Model/GroundWaterTransport/gwt1adv1.f90 b/src/Model/GroundWaterTransport/gwt1adv1.f90 index 0e9f4bdb487..ba54de606d3 100644 --- a/src/Model/GroundWaterTransport/gwt1adv1.f90 +++ b/src/Model/GroundWaterTransport/gwt1adv1.f90 @@ -4,7 +4,7 @@ module GwtAdvModule use ConstantsModule, only: DONE, DZERO, DHALF, DTWO use NumericalPackageModule, only: NumericalPackageType use BaseDisModule, only: DisBaseType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use GwtAdvOptionsModule, only: GwtAdvOptionsType use MatrixBaseModule @@ -17,7 +17,8 @@ module GwtAdvModule integer(I4B), pointer :: iadvwt => null() !< advection scheme (0 up, 1 central, 2 tvd) integer(I4B), dimension(:), pointer, contiguous :: ibound => null() !< pointer to model ibound - type(GwtFmiType), pointer :: fmi => null() !< pointer to fmi object + type(TspFmiType), pointer :: fmi => null() !< pointer to fmi object + real(DP), pointer :: eqnsclfac => null() !< governing equation scale factor; =1. for solute; =rhow*cpw for energy contains @@ -38,7 +39,7 @@ module GwtAdvModule contains - subroutine adv_cr(advobj, name_model, inunit, iout, fmi) + subroutine adv_cr(advobj, name_model, inunit, iout, fmi, eqnsclfac) ! ****************************************************************************** ! adv_cr -- Create a new ADV object ! ****************************************************************************** @@ -50,7 +51,8 @@ subroutine adv_cr(advobj, name_model, inunit, iout, fmi) character(len=*), intent(in) :: name_model integer(I4B), intent(in) :: inunit integer(I4B), intent(in) :: iout - type(GwtFmiType), intent(in), target :: fmi + type(TspFmiType), intent(in), target :: fmi + real(DP), intent(in), pointer :: eqnsclfac !< governing equation scale factor ! ------------------------------------------------------------------------------ ! ! -- Create the object @@ -66,6 +68,7 @@ subroutine adv_cr(advobj, name_model, inunit, iout, fmi) advobj%inunit = inunit advobj%iout = iout advobj%fmi => fmi + advobj%eqnsclfac => eqnsclfac ! ! -- Return return @@ -74,7 +77,7 @@ end subroutine adv_cr subroutine adv_df(this, adv_options) class(GwtAdvType) :: this type(GwtAdvOptionsType), optional, intent(in) :: adv_options !< the optional options, for when not constructing from file - ! local + ! -- local character(len=*), parameter :: fmtadv = & "(1x,/1x,'ADV-- ADVECTION PACKAGE, VERSION 1, 8/25/2017', & &' INPUT READ FROM UNIT ', i0, //)" @@ -96,7 +99,9 @@ subroutine adv_df(this, adv_options) ! --set options from input arg this%iadvwt = adv_options%iAdvScheme end if - + ! + ! -- Return + return end subroutine adv_df subroutine adv_ar(this, dis, ibound) @@ -110,7 +115,7 @@ subroutine adv_ar(this, dis, ibound) ! -- dummy class(GwtAdvType) :: this class(DisBaseType), pointer, intent(in) :: dis - integer(I4B), dimension(:), pointer, contiguous :: ibound + integer(I4B), dimension(:), pointer, contiguous, intent(in) :: ibound ! -- local ! -- formats ! ------------------------------------------------------------------------------ @@ -152,7 +157,7 @@ subroutine adv_fc(this, nodes, matrix_sln, idxglo, cnew, rhs) if (this%dis%con%mask(ipos) == 0) cycle m = this%dis%con%ja(ipos) if (this%ibound(m) == 0) cycle - qnm = this%fmi%gwfflowja(ipos) + qnm = this%fmi%gwfflowja(ipos) * this%eqnsclfac omega = this%adv_weight(this%iadvwt, ipos, n, m, qnm) call matrix_sln%add_value_pos(idxglo(ipos), qnm * (DONE - omega)) call matrix_sln%add_value_pos(idxglo(idiag), qnm * omega) @@ -269,6 +274,7 @@ function advqtvd(this, n, m, iposnm, cnew) result(qtvd) if (smooth > DZERO) then alimiter = DTWO * smooth / (DONE + smooth) qtvd = DHALF * alimiter * qnm * (cnew(idn) - cnew(iup)) + qtvd = qtvd * this%eqnsclfac end if end if ! @@ -303,7 +309,7 @@ subroutine adv_cq(this, cnew, flowja) do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1 m = this%dis%con%ja(ipos) if (this%ibound(m) == 0) cycle - qnm = this%fmi%gwfflowja(ipos) + qnm = this%fmi%gwfflowja(ipos) * this%eqnsclfac omega = this%adv_weight(this%iadvwt, ipos, n, m, qnm) flowja(ipos) = flowja(ipos) + qnm * omega * cnew(n) + & qnm * (DONE - omega) * cnew(m) diff --git a/src/Model/GroundWaterTransport/gwt1apt1.f90 b/src/Model/GroundWaterTransport/gwt1apt1.f90 index 6f50995ac4c..aa4c28778b5 100644 --- a/src/Model/GroundWaterTransport/gwt1apt1.f90 +++ b/src/Model/GroundWaterTransport/gwt1apt1.f90 @@ -44,7 +44,7 @@ module GwtAptModule use SimModule, only: store_error, store_error_unit, count_errors use SimVariablesModule, only: errmsg use BndModule, only: BndType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use BudgetObjectModule, only: BudgetObjectType, budgetobject_cr use BudgetTermModule, only: BudgetTermType use TableModule, only: TableType, table_cr @@ -93,7 +93,7 @@ module GwtAptModule dimension(:), pointer, contiguous :: featname => null() real(DP), dimension(:), pointer, contiguous :: concfeat => null() !< concentration of the feature real(DP), dimension(:, :), pointer, contiguous :: lauxvar => null() !< auxiliary variable - type(GwtFmiType), pointer :: fmi => null() !< pointer to fmi object + type(TspFmiType), pointer :: fmi => null() !< pointer to fmi object real(DP), dimension(:), pointer, contiguous :: qsto => null() !< mass flux due to storage change real(DP), dimension(:), pointer, contiguous :: ccterm => null() !< mass flux required to maintain constant concentration integer(I4B), pointer :: idxbudfjf => null() !< index of flow ja face in flowbudptr diff --git a/src/Model/GroundWaterTransport/gwt1dsp1.f90 b/src/Model/GroundWaterTransport/gwt1dsp1.f90 index acf48a88052..97f45dba3cf 100644 --- a/src/Model/GroundWaterTransport/gwt1dsp1.f90 +++ b/src/Model/GroundWaterTransport/gwt1dsp1.f90 @@ -4,7 +4,7 @@ module GwtDspModule use ConstantsModule, only: DONE, DZERO, DHALF, DPI use NumericalPackageModule, only: NumericalPackageType use BaseDisModule, only: DisBaseType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use Xt3dModule, only: Xt3dType, xt3d_cr use GwtDspOptionsModule, only: GwtDspOptionsType use MatrixBaseModule @@ -17,7 +17,7 @@ module GwtDspModule type, extends(NumericalPackageType) :: GwtDspType integer(I4B), dimension(:), pointer, contiguous :: ibound => null() ! pointer to GWT model ibound - type(GwtFmiType), pointer :: fmi => null() ! pointer to GWT fmi object + type(TspFmiType), pointer :: fmi => null() ! pointer to GWT fmi object real(DP), dimension(:), pointer, contiguous :: thetam => null() ! pointer to GWT storage porosity (voids per aquifer volume) real(DP), dimension(:), pointer, contiguous :: diffc => null() ! molecular diffusion coefficient for each cell real(DP), dimension(:), pointer, contiguous :: alh => null() ! longitudinal horizontal dispersivity @@ -88,7 +88,7 @@ subroutine dsp_cr(dspobj, name_model, input_mempath, inunit, iout, fmi) character(len=*), intent(in) :: input_mempath integer(I4B), intent(in) :: inunit integer(I4B), intent(in) :: iout - type(GwtFmiType), intent(in), target :: fmi + type(TspFmiType), intent(in), target :: fmi ! -- locals ! -- formats character(len=*), parameter :: fmtdsp = & diff --git a/src/Model/GroundWaterTransport/gwt1ist1.f90 b/src/Model/GroundWaterTransport/gwt1ist1.f90 index 99ed8128e4d..ad0a91ee4a0 100644 --- a/src/Model/GroundWaterTransport/gwt1ist1.f90 +++ b/src/Model/GroundWaterTransport/gwt1ist1.f90 @@ -19,7 +19,7 @@ module GwtIstModule LENBUDTXT, DHNOFLO use BndModule, only: BndType use BudgetModule, only: BudgetType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use GwtMstModule, only: GwtMstType, get_zero_order_decay use OutputControlDataModule, only: OutputControlDataType use MatrixBaseModule @@ -49,7 +49,7 @@ module GwtIstModule !< type, extends(BndType) :: GwtIstType - type(GwtFmiType), pointer :: fmi => null() !< pointer to fmi object + type(TspFmiType), pointer :: fmi => null() !< pointer to fmi object type(GwtMstType), pointer :: mst => null() !< pointer to mst object integer(I4B), pointer :: icimout => null() !< unit number for binary cim output @@ -116,7 +116,7 @@ subroutine ist_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, & character(len=*), intent(in) :: pakname !< name of the package ! -- local type(GwtIstType), pointer :: istobj - type(GwtFmiType), pointer :: fmi + type(TspFmiType), pointer :: fmi type(GwtMstType), pointer :: mst ! ! -- allocate the object and assign values to object variables diff --git a/src/Model/GroundWaterTransport/gwt1lkt1.f90 b/src/Model/GroundWaterTransport/gwt1lkt1.f90 index 98ef40abcd0..8221859bf81 100644 --- a/src/Model/GroundWaterTransport/gwt1lkt1.f90 +++ b/src/Model/GroundWaterTransport/gwt1lkt1.f90 @@ -37,7 +37,7 @@ module GwtLktModule use ConstantsModule, only: DZERO, DONE, LINELENGTH use SimModule, only: store_error use BndModule, only: BndType, GetBndFromList - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use LakModule, only: LakType use ObserveModule, only: ObserveType use GwtAptModule, only: GwtAptType, apt_process_obsID, & @@ -108,7 +108,7 @@ subroutine lkt_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, & integer(I4B), intent(in) :: iout character(len=*), intent(in) :: namemodel character(len=*), intent(in) :: pakname - type(GwtFmiType), pointer :: fmi + type(TspFmiType), pointer :: fmi ! -- local type(GwtLktType), pointer :: lktobj ! ------------------------------------------------------------------------------ diff --git a/src/Model/GroundWaterTransport/gwt1mst1.f90 b/src/Model/GroundWaterTransport/gwt1mst1.f90 index f962a5f2bef..48842e211b5 100644 --- a/src/Model/GroundWaterTransport/gwt1mst1.f90 +++ b/src/Model/GroundWaterTransport/gwt1mst1.f90 @@ -17,7 +17,7 @@ module GwtMstModule use MatrixBaseModule use NumericalPackageModule, only: NumericalPackageType use BaseDisModule, only: DisBaseType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType implicit none public :: GwtMstType @@ -60,7 +60,7 @@ module GwtMstModule ! ! -- misc integer(I4B), dimension(:), pointer, contiguous :: ibound => null() !< pointer to model ibound - type(GwtFmiType), pointer :: fmi => null() !< pointer to fmi object + type(TspFmiType), pointer :: fmi => null() !< pointer to fmi object contains @@ -100,7 +100,7 @@ subroutine mst_cr(mstobj, name_model, inunit, iout, fmi) character(len=*), intent(in) :: name_model !< name of the model integer(I4B), intent(in) :: inunit !< unit number of WEL package input file integer(I4B), intent(in) :: iout !< unit number of model listing file - type(GwtFmiType), intent(in), target :: fmi !< fmi package for this GWT model + type(TspFmiType), intent(in), target :: fmi !< fmi package for this GWT model ! ! -- Create the object allocate (mstobj) diff --git a/src/Model/GroundWaterTransport/gwt1mvt1.f90 b/src/Model/GroundWaterTransport/gwt1mvt1.f90 index 732b2e59ac3..1f045772acc 100644 --- a/src/Model/GroundWaterTransport/gwt1mvt1.f90 +++ b/src/Model/GroundWaterTransport/gwt1mvt1.f90 @@ -11,7 +11,7 @@ module GwtMvtModule use SimModule, only: store_error use BaseDisModule, only: DisBaseType use NumericalPackageModule, only: NumericalPackageType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use BudgetModule, only: BudgetType, budget_cr use BudgetObjectModule, only: BudgetObjectType, budgetobject_cr use TableModule, only: TableType, table_cr @@ -28,8 +28,8 @@ module GwtMvtModule integer(I4B), pointer :: maxpackages !< max number of packages integer(I4B), pointer :: ibudgetout => null() !< unit number for budget output file integer(I4B), pointer :: ibudcsv => null() !< unit number for csv budget output file - type(GwtFmiType), pointer :: fmi1 => null() !< pointer to fmi object for model 1 - type(GwtFmiType), pointer :: fmi2 => null() !< pointer to fmi object for model 2 (set to fmi1 for single model) + type(TspFmiType), pointer :: fmi1 => null() !< pointer to fmi object for model 1 + type(TspFmiType), pointer :: fmi2 => null() !< pointer to fmi object for model 2 (set to fmi1 for single model) type(BudgetType), pointer :: budget => null() !< mover transport budget object (used to write balance table) type(BudgetObjectType), pointer :: budobj => null() !< budget container (used to write binary file) type(BudgetObjectType), pointer :: mvrbudobj => null() !< pointer to the water mover budget object @@ -75,10 +75,10 @@ subroutine mvt_cr(mvt, name_model, inunit, iout, fmi1, gwfmodelname1, & character(len=*), intent(in) :: name_model integer(I4B), intent(in) :: inunit integer(I4B), intent(in) :: iout - type(GwtFmiType), intent(in), target :: fmi1 + type(TspFmiType), intent(in), target :: fmi1 character(len=*), intent(in), optional :: gwfmodelname1 character(len=*), intent(in), optional :: gwfmodelname2 - type(GwtFmiType), intent(in), target, optional :: fmi2 + type(TspFmiType), intent(in), target, optional :: fmi2 ! ------------------------------------------------------------------------------ ! ! -- Create the object @@ -251,8 +251,8 @@ subroutine mvt_fc(this, cnew1, cnew2) real(DP) :: q, cp real(DP), dimension(:), pointer :: concpak real(DP), dimension(:), contiguous, pointer :: cnew - type(GwtFmiType), pointer :: fmi_pr !< pointer to provider model fmi package - type(GwtFmiType), pointer :: fmi_rc !< pointer to receiver model fmi package + type(TspFmiType), pointer :: fmi_pr !< pointer to provider model fmi package + type(TspFmiType), pointer :: fmi_rc !< pointer to receiver model fmi package ! ------------------------------------------------------------------------------ ! ! -- Add mover QC terms to the receiver packages @@ -337,8 +337,8 @@ subroutine set_fmi_pr_rc(this, ibudterm, fmi_pr, fmi_rc) ! -- dummy class(GwtMvtType) :: this integer(I4B), intent(in) :: ibudterm - type(GwtFmiType), pointer :: fmi_pr - type(GwtFmiType), pointer :: fmi_rc + type(TspFmiType), pointer :: fmi_pr + type(TspFmiType), pointer :: fmi_rc fmi_pr => null() fmi_rc => null() @@ -817,8 +817,8 @@ subroutine mvt_fill_budobj(this, cnew1, cnew2) real(DP), intent(in), dimension(:), contiguous, target :: cnew1 real(DP), intent(in), dimension(:), contiguous, target :: cnew2 ! -- local - type(GwtFmiType), pointer :: fmi_pr - type(GwtFmiType), pointer :: fmi_rc + type(TspFmiType), pointer :: fmi_pr + type(TspFmiType), pointer :: fmi_rc real(DP), dimension(:), contiguous, pointer :: cnew integer(I4B) :: nbudterm integer(I4B) :: nlist diff --git a/src/Model/GroundWaterTransport/gwt1mwt1.f90 b/src/Model/GroundWaterTransport/gwt1mwt1.f90 index 15137d3a5c6..5b286a6b1fb 100644 --- a/src/Model/GroundWaterTransport/gwt1mwt1.f90 +++ b/src/Model/GroundWaterTransport/gwt1mwt1.f90 @@ -38,7 +38,7 @@ module GwtMwtModule use ConstantsModule, only: DZERO, LINELENGTH use SimModule, only: store_error use BndModule, only: BndType, GetBndFromList - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use MawModule, only: MawType use ObserveModule, only: ObserveType use GwtAptModule, only: GwtAptType, apt_process_obsID, & @@ -101,7 +101,7 @@ subroutine mwt_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, & integer(I4B), intent(in) :: iout character(len=*), intent(in) :: namemodel character(len=*), intent(in) :: pakname - type(GwtFmiType), pointer :: fmi + type(TspFmiType), pointer :: fmi ! -- local type(GwtMwtType), pointer :: mwtobj ! ------------------------------------------------------------------------------ diff --git a/src/Model/GroundWaterTransport/gwt1sft1.f90 b/src/Model/GroundWaterTransport/gwt1sft1.f90 index fe310f5eb42..362f615363a 100644 --- a/src/Model/GroundWaterTransport/gwt1sft1.f90 +++ b/src/Model/GroundWaterTransport/gwt1sft1.f90 @@ -36,7 +36,7 @@ module GwtSftModule use ConstantsModule, only: DZERO, DONE, LINELENGTH use SimModule, only: store_error use BndModule, only: BndType, GetBndFromList - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use SfrModule, only: SfrType use ObserveModule, only: ObserveType use GwtAptModule, only: GwtAptType, apt_process_obsID, & @@ -105,7 +105,7 @@ subroutine sft_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, & integer(I4B), intent(in) :: iout character(len=*), intent(in) :: namemodel character(len=*), intent(in) :: pakname - type(GwtFmiType), pointer :: fmi + type(TspFmiType), pointer :: fmi ! -- local type(GwtSftType), pointer :: sftobj ! ------------------------------------------------------------------------------ diff --git a/src/Model/GroundWaterTransport/gwt1ssm1.f90 b/src/Model/GroundWaterTransport/gwt1ssm1.f90 index e8684820918..b51164f3c74 100644 --- a/src/Model/GroundWaterTransport/gwt1ssm1.f90 +++ b/src/Model/GroundWaterTransport/gwt1ssm1.f90 @@ -15,7 +15,7 @@ module GwtSsmModule use SimVariablesModule, only: errmsg use NumericalPackageModule, only: NumericalPackageType use BaseDisModule, only: DisBaseType - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use TableModule, only: TableType, table_cr use GwtSpcModule, only: GwtSpcType use MatrixBaseModule @@ -41,7 +41,7 @@ module GwtSsmModule integer(I4B), dimension(:), pointer, contiguous :: iauxpak => null() !< aux col for concentration integer(I4B), dimension(:), pointer, contiguous :: ibound => null() !< pointer to model ibound real(DP), dimension(:), pointer, contiguous :: cnew => null() !< pointer to gwt%x - type(GwtFmiType), pointer :: fmi => null() !< pointer to fmi object + type(TspFmiType), pointer :: fmi => null() !< pointer to fmi object type(TableType), pointer :: outputtab => null() !< output table object type(GwtSpcType), dimension(:), pointer :: ssmivec => null() !< array of stress package concentration objects @@ -84,7 +84,7 @@ subroutine ssm_cr(ssmobj, name_model, inunit, iout, fmi) character(len=*), intent(in) :: name_model !< name of the model integer(I4B), intent(in) :: inunit !< fortran unit for input integer(I4B), intent(in) :: iout !< fortran unit for output - type(GwtFmiType), intent(in), target :: fmi !< GWT FMI package + type(TspFmiType), intent(in), target :: fmi !< GWT FMI package ! ! -- Create the object allocate (ssmobj) diff --git a/src/Model/GroundWaterTransport/gwt1uzt1.f90 b/src/Model/GroundWaterTransport/gwt1uzt1.f90 index c6be55aec38..83fbf2efa2b 100644 --- a/src/Model/GroundWaterTransport/gwt1uzt1.f90 +++ b/src/Model/GroundWaterTransport/gwt1uzt1.f90 @@ -30,7 +30,7 @@ module GwtUztModule use ConstantsModule, only: DZERO, DONE, LINELENGTH use SimModule, only: store_error use BndModule, only: BndType, GetBndFromList - use GwtFmiModule, only: GwtFmiType + use TspFmiModule, only: TspFmiType use UzfModule, only: UzfType use ObserveModule, only: ObserveType use GwtAptModule, only: GwtAptType, apt_process_obsID, & @@ -93,7 +93,7 @@ subroutine uzt_create(packobj, id, ibcnum, inunit, iout, namemodel, pakname, & integer(I4B), intent(in) :: iout character(len=*), intent(in) :: namemodel character(len=*), intent(in) :: pakname - type(GwtFmiType), pointer :: fmi + type(TspFmiType), pointer :: fmi ! -- local type(GwtUztType), pointer :: uztobj ! ------------------------------------------------------------------------------ diff --git a/src/Model/ModelUtilities/FlowModelInterface.f90 b/src/Model/ModelUtilities/FlowModelInterface.f90 index 68d1ce7f1c2..c55787568d1 100644 --- a/src/Model/ModelUtilities/FlowModelInterface.f90 +++ b/src/Model/ModelUtilities/FlowModelInterface.f90 @@ -2,7 +2,7 @@ module FlowModelInterfaceModule use KindModule, only: DP, I4B, LGP use ConstantsModule, only: DONE, DZERO, DHALF, LINELENGTH, LENBUDTXT, & - LENPACKAGENAME + LENPACKAGENAME, LENVARNAME use SimModule, only: store_error, store_error_unit use SimVariablesModule, only: errmsg use NumericalPackageModule, only: NumericalPackageType @@ -43,6 +43,8 @@ module FlowModelInterfaceModule type(PackageBudgetType), dimension(:), allocatable :: gwfpackages !< used to get flows between a package and gwf type(BudgetObjectType), pointer :: mvrbudobj => null() !< pointer to the mover budget budget object character(len=16), dimension(:), allocatable :: flowpacknamearray !< array of boundary package names (e.g. LAK-1, SFR-3, etc.) + character(len=LENVARNAME) :: depvartype = '' + contains procedure :: advance_bfr @@ -69,6 +71,7 @@ module FlowModelInterfaceModule contains !> @brief Define the flow model interface + !< subroutine fmi_df(this, dis) ! -- modules use SimModule, only: store_error @@ -120,6 +123,7 @@ subroutine fmi_df(this, dis) end subroutine fmi_df !> @brief Allocate the package + !< subroutine fmi_ar(this, ibound) ! -- modules use SimModule, only: store_error @@ -138,6 +142,7 @@ subroutine fmi_ar(this, ibound) end subroutine fmi_ar !> @brief Deallocate variables + !< subroutine fmi_da(this) ! -- modules use MemoryManagerModule, only: mem_deallocate @@ -183,6 +188,7 @@ subroutine fmi_da(this) end subroutine fmi_da !> @brief Allocate scalars + !< subroutine allocate_scalars(this) ! -- modules use MemoryManagerModule, only: mem_allocate, mem_setptr @@ -219,6 +225,7 @@ subroutine allocate_scalars(this) end subroutine allocate_scalars !> @brief Allocate arrays + !< subroutine allocate_arrays(this, nodes) use MemoryManagerModule, only: mem_allocate !modules @@ -281,6 +288,7 @@ subroutine allocate_arrays(this, nodes) end subroutine allocate_arrays !> @brief Read options from input file + !< subroutine read_options(this) ! -- modules use ConstantsModule, only: LINELENGTH, DEM6 @@ -322,6 +330,7 @@ subroutine read_options(this) end subroutine read_options !> @brief Read packagedata block from input file + !< subroutine read_packagedata(this) ! -- modules use OpenSpecModule, only: ACCESS, FORM @@ -421,6 +430,7 @@ subroutine read_packagedata(this) end subroutine read_packagedata !> @brief Initialize the budget file reader + !< subroutine initialize_bfr(this) ! -- modules class(FlowModelInterfaceType) :: this @@ -585,6 +595,7 @@ subroutine advance_bfr(this) end subroutine advance_bfr !> @brief Finalize the budget file reader + !< subroutine finalize_bfr(this) ! -- modules class(FlowModelInterfaceType) :: this @@ -596,6 +607,7 @@ subroutine finalize_bfr(this) end subroutine finalize_bfr !> @brief Initialize the head file reader + !< subroutine initialize_hfr(this) ! -- modules class(FlowModelInterfaceType) :: this @@ -609,6 +621,7 @@ subroutine initialize_hfr(this) end subroutine initialize_hfr !> @brief Advance the head file reader + !< subroutine advance_hfr(this) ! -- modules use TdisModule, only: kstp, kper @@ -701,6 +714,7 @@ subroutine advance_hfr(this) end subroutine advance_hfr !> @brief Finalize the head file reader + !< subroutine finalize_hfr(this) ! -- modules class(FlowModelInterfaceType) :: this @@ -816,6 +830,7 @@ subroutine initialize_gwfterms_from_bfr(this) end subroutine initialize_gwfterms_from_bfr !> @brief Initialize gwf terms from a GWF exchange + !< subroutine initialize_gwfterms_from_gwfbndlist(this) ! -- modules use BndModule, only: BndType, GetBndFromList @@ -883,7 +898,6 @@ end subroutine initialize_gwfterms_from_gwfbndlist !! gwfpackages is an array of PackageBudget objects. !! This routine allocates gwfpackages to the proper size and initializes some !! member variables. - !! !< subroutine allocate_gwfpackages(this, ngwfterms) ! -- modules @@ -920,6 +934,7 @@ subroutine allocate_gwfpackages(this, ngwfterms) end subroutine allocate_gwfpackages !> @brief Deallocate memory in the gwfpackages array + !< subroutine deallocate_gwfpackages(this) ! -- modules ! -- dummy @@ -937,6 +952,7 @@ subroutine deallocate_gwfpackages(this) end subroutine deallocate_gwfpackages !> @brief Find the package index for the package with the given name + !< subroutine get_package_index(this, name, idx) use BndModule, only: BndType, GetBndFromList class(FlowModelInterfaceType) :: this diff --git a/src/Model/TransportModel/tsp1.f90 b/src/Model/TransportModel/tsp1.f90 index f37082cc3cb..cef33d6bb0b 100644 --- a/src/Model/TransportModel/tsp1.f90 +++ b/src/Model/TransportModel/tsp1.f90 @@ -6,9 +6,14 @@ module TransportModelModule use KindModule, only: DP, I4B - use ConstantsModule, only: LENFTYPE + use VersionModule, only: write_listfile_header + use ConstantsModule, only: LENFTYPE, LINELENGTH, DZERO, LENPAKLOC, & + LENMEMPATH, LENVARNAME use SimVariablesModule, only: errmsg use NumericalModelModule, only: NumericalModelType + use TspFmiModule, only: TspFmiType + use BudgetModule, only: BudgetType + use MatrixBaseModule implicit none @@ -18,8 +23,530 @@ module TransportModelModule type, extends(NumericalModelType) :: TransportModelType + type(TspFmiType), pointer :: fmi => null() ! flow model interface + type(BudgetType), pointer :: budget => null() !< budget object + integer(I4B), pointer :: infmi => null() ! unit number FMI + real(DP), pointer :: eqnsclfac => null() !< constant factor by which all terms in the model's governing equation are scaled (divided) for formulation and solution + ! Labels that will be defined + character(len=LENVARNAME) :: tsptype = '' !< "solute" or "heat" + character(len=LENVARNAME) :: depvartype = '' !< "concentration" or "temperature" + character(len=LENVARNAME) :: depvarunit = '' !< "mass" or "energy" + character(len=LENVARNAME) :: depvarunitabbrev = '' !< "M" or "E" + contains + ! -- public + procedure, public :: tsp_cr + procedure, public :: tsp_df + procedure, public :: tsp_da + procedure, public :: tsp_ac + procedure, public :: tsp_mc + procedure, public :: tsp_ar + procedure, public :: tsp_rp + procedure, public :: tsp_ad + procedure, public :: tsp_fc + procedure, public :: tsp_cc + procedure, public :: tsp_cq + procedure, public :: tsp_bd + procedure, public :: allocate_tsp_scalars + procedure, public :: set_tsp_labels + procedure, public :: ftype_check + ! -- private + procedure, private :: create_lstfile + procedure, private :: create_tsp_packages + procedure, private :: log_namfile_options + end type TransportModelType +contains + + !> @brief Create a new generalized transport model object + !! + !! Create a new transport model that will be further refined into GWT or GWE + !< + subroutine tsp_cr(this, filename, id, modelname, indis) + ! -- modules + use GwfNamInputModule, only: GwfNamParamFoundType + use BudgetModule, only: budget_cr + ! -- dummy + class(TransportModelType) :: this + character(len=*), intent(in) :: filename + integer(I4B), intent(in) :: id + integer(I4B), intent(inout) :: indis + character(len=*), intent(in) :: modelname + ! -- local + character(len=LINELENGTH) :: lst_fname + type(GwfNamParamFoundType) :: found + ! + ! -- create the list file + call this%create_lstfile(lst_fname, filename, found%list) + ! + ! -- log set options + if (this%iout > 0) then + call this%log_namfile_options(found) + end if + ! + ! -- Create utility objects + call budget_cr(this%budget, this%name) + ! + ! -- create model packages + call this%create_tsp_packages(indis) + ! + ! -- Return + return + end subroutine tsp_cr + + !> @brief Generalized transport model define model + !! + !! This subroutine extended by either GWT or GWE. This routine calls the + !! define (df) routines for each attached package and sets variables and + !! pointers. + !< + subroutine tsp_df(this) + ! -- dummy variables + class(TransportModelType) :: this + ! + ! -- Return + return + end subroutine tsp_df + + !> @brief Generalized transport model add connections + !! + !! This subroutine extended by either GWT or GWE. This routine adds the + !! internal connections of this model to the sparse matrix + !< + subroutine tsp_ac(this, sparse) + ! -- modules + use SparseModule, only: sparsematrix + ! -- dummy variables + class(TransportModelType) :: this + type(sparsematrix), intent(inout) :: sparse + ! -- local +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_ac + + !> @brief Generalized transport model map coefficients + !! + !! This subroutine extended by either GWT or GWE. This routine maps the + !! positions of this models connections in the numerical solution coefficient + !! matrix. + !< + subroutine tsp_mc(this, matrix_sln) + ! -- dummy + class(TransportModelType) :: this + class(MatrixBaseType), pointer :: matrix_sln !< global system matrix + ! -- local +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_mc + + !> @brief Generalized transport model allocate and read + !! + !! This subroutine extended by either GWT or GWE. This routine calls + !! the allocate and reads (ar) routines of attached packages and allocates + !! memory for arrays required by the model object. + !< + subroutine tsp_ar(this) + ! -- dummy variables + class(TransportModelType) :: this +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_ar + + !> @brief Generalized transport model read and prepare + !! + !! This subroutine extended by either GWT or GWE. This routine calls + !! the read and prepare (rp) routines of attached packages. + !< + subroutine tsp_rp(this) + ! -- dummy variables + class(TransportModelType) :: this +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_rp + + !> @brief Generalized transport model time step advance + !! + !! This subroutine extended by either GWT or GWE. This routine calls + !! the advance time step (ad) routines of attached packages. + !< + subroutine tsp_ad(this) + ! -- dummy variables + class(TransportModelType) :: this +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_ad + + !> @brief Generalized transport model fill coefficients + !! + !! This subroutine extended by either GWT or GWE. This routine calls + !! the fill coefficients (fc) routines of attached packages. + !< + subroutine tsp_fc(this, kiter, matrix_sln, inwtflag) + ! -- dummy variables + class(TransportModelType) :: this + integer(I4B), intent(in) :: kiter + class(MatrixBaseType), pointer :: matrix_sln + integer(I4B), intent(in) :: inwtflag +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_fc + + !> @brief Generalized transport model final convergence check + !! + !! This subroutine extended by either GWT or GWE. This routine calls + !! the convergence check (cc) routines of attached packages. + !< + subroutine tsp_cc(this, innertot, kiter, iend, icnvgmod, cpak, ipak, dpak) + ! -- dummy + class(TransportModelType) :: this + integer(I4B), intent(in) :: innertot + integer(I4B), intent(in) :: kiter + integer(I4B), intent(in) :: iend + integer(I4B), intent(in) :: icnvgmod + character(len=LENPAKLOC), intent(inout) :: cpak + integer(I4B), intent(inout) :: ipak + real(DP), intent(inout) :: dpak + ! -- local +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_cc + + !> @brief Generalized transport model calculate flows + !! + !! This subroutine extended by either GWT or GWE. This routine calculates + !! intercell flows (flowja) + !< + subroutine tsp_cq(this, icnvg, isuppress_output) + ! -- dummy variables + class(TransportModelType) :: this + integer(I4B), intent(in) :: icnvg + integer(I4B), intent(in) :: isuppress_output + ! -- local +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_cq + + !> @brief Generalized transport model budget + !! + !! This subroutine extended by either GWT or GWE. This routine calculates + !! package contributions to model budget + !< + subroutine tsp_bd(this, icnvg, isuppress_output) + ! -- dummy + class(TransportModelType) :: this + integer(I4B), intent(in) :: icnvg + integer(I4B), intent(in) :: isuppress_output +! ------------------------------------------------------------------------------ + ! + ! -- Return + return + end subroutine tsp_bd + + !> @brief Allocate scalar variables for transport model + !! + !! Method to allocate memory for non-allocatable members. + !< + subroutine allocate_tsp_scalars(this, modelname) + ! -- modules + use MemoryManagerModule, only: mem_allocate + ! -- dummy + class(TransportModelType) :: this + character(len=*), intent(in) :: modelname +! ------------------------------------------------------------------------------ + ! + ! -- allocate members from (grand)parent class + call this%NumericalModelType%allocate_scalars(modelname) + ! + ! -- allocate members that are part of model class + call mem_allocate(this%infmi, 'INFMI', this%memoryPath) + call mem_allocate(this%eqnsclfac, 'EQNSCLFAC', this%memoryPath) + ! + this%infmi = 0 + this%eqnsclfac = DZERO + ! + ! -- Return + return + end subroutine allocate_tsp_scalars + + !> @brief Define the labels corresponding to the flavor of + !! transport model + !! + !! Set variable names according to type of transport model + !< + subroutine set_tsp_labels(this, tsptype, depvartype, depvarunit, & + depvarunitabbrev) + class(TransportModelType) :: this + character(len=*), intent(in), pointer :: tsptype !< type of model, default is GWT (alternative is GWE) + character(len=*), intent(in) :: depvartype !< dependent variable type, default is "CONCENTRATION" + character(len=*), intent(in) :: depvarunit !< units of dependent variable for writing to list file + character(len=*), intent(in) :: depvarunitabbrev !< abbreviation of associated units + ! + ! -- Set the model type + this%tsptype = tsptype + ! + ! -- Set the type of dependent variable being solved for + this%depvartype = depvartype + ! + ! -- Set the units associated with the dependent variable + this%depvarunit = depvarunit + ! + ! -- Set the units abbreviation + this%depvarunitabbrev = depvarunitabbrev + ! + ! -- Return + return + end subroutine set_tsp_labels + + !> @brief Deallocate memory + !! + !! Deallocate memmory at conclusion of model run + !< + subroutine tsp_da(this) + ! -- modules + use MemoryManagerModule, only: mem_deallocate + ! -- dummy + class(TransportModelType) :: this + ! -- local + ! + ! -- Internal flow packages deallocate + call this%fmi%fmi_da() + ! + ! -- Internal package objects + deallocate (this%fmi) + ! + ! -- Scalars + call mem_deallocate(this%infmi) + call mem_deallocate(this%eqnsclfac) + ! + ! -- Return + return + end subroutine tsp_da + + !> @brief Generalized tranpsort model routine + !! + !! Check to make sure required input files have been specified + !< + subroutine ftype_check(this, indis, inmst, inic) + ! -- modules + use ConstantsModule, only: LINELENGTH + use SimModule, only: store_error, count_errors, store_error_filename + ! -- dummy + class(TransportModelType) :: this + integer(I4B), intent(in) :: indis + integer(I4B), intent(in) :: inmst + integer(I4B), intent(in) :: inic + ! -- local + character(len=LINELENGTH) :: errmsg +! ------------------------------------------------------------------------------ + ! + ! -- Check for IC6, DIS(u), and MST. Stop if not present. + if (inic == 0) then + write (errmsg, '(a)') & + 'Initial conditions (IC6) package not specified.' + call store_error(errmsg) + end if + if (indis == 0) then + write (errmsg, '(a)') & + 'Discretization (DIS6 or DISU6) package not specified.' + call store_error(errmsg) + end if + if (inmst == 0) then + write (errmsg, '(a)') 'Mass storage and transfer (MST6) & + &package not specified.' + call store_error(errmsg) + end if + ! + if (count_errors() > 0) then + write (errmsg, '(a)') 'Required package(s) not specified.' + call store_error(errmsg) + call store_error_filename(this%filename) + end if + ! + ! -- Return + return + end subroutine ftype_check + + subroutine create_lstfile(this, lst_fname, model_fname, defined) + ! -- modules + use KindModule, only: LGP + use InputOutputModule, only: openfile, getunit + ! -- dummy + class(TransportModelType) :: this + character(len=*), intent(inout) :: lst_fname + character(len=*), intent(in) :: model_fname + logical(LGP), intent(in) :: defined + ! -- local + integer(I4B) :: i, istart, istop + ! + ! -- set list file name if not provided + if (.not. defined) then + ! + ! -- initialize + lst_fname = ' ' + istart = 0 + istop = len_trim(model_fname) + ! + ! -- identify '.' character position from back of string + do i = istop, 1, -1 + if (model_fname(i:i) == '.') then + istart = i + exit + end if + end do + ! + ! -- if not found start from string end + if (istart == 0) istart = istop + 1 + ! + ! -- set list file name + lst_fname = model_fname(1:istart) + istop = istart + 3 + lst_fname(istart:istop) = '.lst' + end if + ! + ! -- create the list file + this%iout = getunit() + call openfile(this%iout, 0, lst_fname, 'LIST', filstat_opt='REPLACE') + ! + ! -- write list file header + call write_listfile_header(this%iout, 'GROUNDWATER TRANSPORT MODEL (GWT)') + ! + ! -- Return + return + end subroutine create_lstfile + + !> @brief Write model namfile options to list file + !< + subroutine log_namfile_options(this, found) + use GwfNamInputModule, only: GwfNamParamFoundType + class(TransportModelType) :: this + type(GwfNamParamFoundType), intent(in) :: found + ! + write (this%iout, '(1x,a)') 'NAMEFILE OPTIONS:' + ! + if (found%newton) then + write (this%iout, '(4x,a)') & + 'NEWTON-RAPHSON method enabled for the model.' + if (found%under_relaxation) then + write (this%iout, '(4x,a,a)') & + 'NEWTON-RAPHSON UNDER-RELAXATION based on the bottom ', & + 'elevation of the model will be applied to the model.' + end if + end if + ! + if (found%print_input) then + write (this%iout, '(4x,a)') 'STRESS PACKAGE INPUT WILL BE PRINTED '// & + 'FOR ALL MODEL STRESS PACKAGES' + end if + ! + if (found%print_flows) then + write (this%iout, '(4x,a)') 'PACKAGE FLOWS WILL BE PRINTED '// & + 'FOR ALL MODEL PACKAGES' + end if + ! + if (found%save_flows) then + write (this%iout, '(4x,a)') & + 'FLOWS WILL BE SAVED TO BUDGET FILE SPECIFIED IN OUTPUT CONTROL' + end if + ! + write (this%iout, '(1x,a)') 'END NAMEFILE OPTIONS:' + ! + ! -- Return + return + end subroutine log_namfile_options + + !> @brief Source package info and begin to process + !< + subroutine create_tsp_packages(this, indis) + ! -- modules + use ConstantsModule, only: LINELENGTH, LENPACKAGENAME + use CharacterStringModule, only: CharacterStringType + use ArrayHandlersModule, only: expandarray + use MemoryManagerModule, only: mem_setptr + use MemoryHelperModule, only: create_mem_path + use SimVariablesModule, only: idm_context + use GwfDisModule, only: dis_cr + use GwfDisvModule, only: disv_cr + use GwfDisuModule, only: disu_cr + use TspFmiModule, only: fmi_cr + ! -- dummy + class(TransportModelType) :: this + integer(I4B), intent(inout) :: indis ! DIS enabled flag + ! -- local + type(CharacterStringType), dimension(:), contiguous, & + pointer :: pkgtypes => null() + type(CharacterStringType), dimension(:), contiguous, & + pointer :: pkgnames => null() + type(CharacterStringType), dimension(:), contiguous, & + pointer :: mempaths => null() + integer(I4B), dimension(:), contiguous, & + pointer :: inunits => null() + character(len=LENMEMPATH) :: model_mempath + character(len=LENFTYPE) :: pkgtype + character(len=LENPACKAGENAME) :: pkgname + character(len=LENMEMPATH) :: mempath + integer(I4B), pointer :: inunit + integer(I4B) :: n + ! + ! -- Initialize + indis = 0 + ! + ! -- set input memory paths, input/model and input/model/namfile + model_mempath = create_mem_path(component=this%name, context=idm_context) + ! + ! -- set pointers to model path package info + call mem_setptr(pkgtypes, 'PKGTYPES', model_mempath) + call mem_setptr(pkgnames, 'PKGNAMES', model_mempath) + call mem_setptr(mempaths, 'MEMPATHS', model_mempath) + call mem_setptr(inunits, 'INUNITS', model_mempath) + ! + do n = 1, size(pkgtypes) + ! + ! attributes for this input package + pkgtype = pkgtypes(n) + pkgname = pkgnames(n) + mempath = mempaths(n) + inunit => inunits(n) + ! + ! -- create dis package as it is a prerequisite for other packages + select case (pkgtype) + case ('DIS6') + indis = 1 + call dis_cr(this%dis, this%name, mempath, indis, this%iout) + case ('DISV6') + indis = 1 + call disv_cr(this%dis, this%name, mempath, indis, this%iout) + case ('DISU6') + indis = 1 + call disu_cr(this%dis, this%name, mempath, indis, this%iout) + case ('FMI6') + this%infmi = inunit + end select + end do + ! + ! -- Create packages that are tied directly to model + call fmi_cr(this%fmi, this%name, this%infmi, this%iout, this%eqnsclfac, & + this%depvartype) + ! + ! -- Return + return + end subroutine create_tsp_packages + end module TransportModelModule diff --git a/src/Model/GroundWaterTransport/gwt1fmi1.f90 b/src/Model/TransportModel/tsp1fmi1.f90 similarity index 90% rename from src/Model/GroundWaterTransport/gwt1fmi1.f90 rename to src/Model/TransportModel/tsp1fmi1.f90 index 4e33a8af884..fe92e15fdcf 100644 --- a/src/Model/GroundWaterTransport/gwt1fmi1.f90 +++ b/src/Model/TransportModel/tsp1fmi1.f90 @@ -1,8 +1,8 @@ -module GwtFmiModule +module TspFmiModule use KindModule, only: DP, I4B use ConstantsModule, only: DONE, DZERO, DHALF, LINELENGTH, LENBUDTXT, & - LENPACKAGENAME + LENPACKAGENAME, LENVARNAME use SimModule, only: store_error, store_error_unit use SimVariablesModule, only: errmsg use FlowModelInterfaceModule, only: FlowModelInterfaceType @@ -16,7 +16,7 @@ module GwtFmiModule implicit none private - public :: GwtFmiType + public :: TspFmiType public :: fmi_cr character(len=LENPACKAGENAME) :: text = ' GWTFMI' @@ -34,14 +34,16 @@ module GwtFmiModule type(BudgetObjectType), pointer :: ptr end type BudObjPtrArray - type, extends(FlowModelInterfaceType) :: GwtFmiType + type, extends(FlowModelInterfaceType) :: TspFmiType integer(I4B), dimension(:), pointer, contiguous :: iatp => null() !< advanced transport package applied to gwfpackages integer(I4B), pointer :: iflowerr => null() !< add the flow error correction real(DP), dimension(:), pointer, contiguous :: flowcorrect => null() !< mass flow correction + real(DP), pointer :: eqnsclfac => null() !< governing equation scale factor; =1. for solute; =rhow*cpw for energy type(DataAdvancedPackageType), & dimension(:), pointer, contiguous :: datp => null() type(BudObjPtrArray), dimension(:), allocatable :: aptbudobj !< flow budget objects for the advanced packages + contains procedure :: allocate_arrays => gwtfmi_allocate_arrays @@ -62,17 +64,20 @@ module GwtFmiModule procedure :: set_aptbudobj_pointer procedure :: read_packagedata => gwtfmi_read_packagedata - end type GwtFmiType + end type TspFmiType contains - !> @brief Create a new FMI object - subroutine fmi_cr(fmiobj, name_model, inunit, iout) + !> @breif Create a new FMI object + !< + subroutine fmi_cr(fmiobj, name_model, inunit, iout, eqnsclfac, depvartype) ! -- dummy - type(GwtFmiType), pointer :: fmiobj + type(TspFmiType), pointer :: fmiobj character(len=*), intent(in) :: name_model integer(I4B), intent(in) :: inunit integer(I4B), intent(in) :: iout + real(DP), intent(in), pointer :: eqnsclfac !< governing equation scale factor + character(len=LENVARNAME), intent(in) :: depvartype ! ! -- Create the object allocate (fmiobj) @@ -91,16 +96,23 @@ subroutine fmi_cr(fmiobj, name_model, inunit, iout) ! -- Initialize block parser call fmiobj%parser%Initialize(fmiobj%inunit, fmiobj%iout) ! + ! -- Assign label based on dependent variable + fmiobj%depvartype = depvartype + ! + ! -- Store pointer to governing equation scale factor + fmiobj%eqnsclfac => eqnsclfac + ! ! -- Return return end subroutine fmi_cr !> @brief Read and prepare + !< subroutine fmi_rp(this, inmvr) ! -- modules use TdisModule, only: kper, kstp ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer(I4B), intent(in) :: inmvr ! -- local ! -- formats @@ -126,12 +138,13 @@ subroutine fmi_rp(this, inmvr) return end subroutine fmi_rp - !> @brief Advance + !> @brief Advance routine for FMI object + !< subroutine fmi_ad(this, cnew) ! -- modules use ConstantsModule, only: DHDRY ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this real(DP), intent(inout), dimension(:) :: cnew ! -- local integer(I4B) :: n @@ -231,10 +244,12 @@ subroutine fmi_ad(this, cnew) return end subroutine fmi_ad - !> @brief Calculate coefficients and fill matrix and rhs + !> @brief Calculate coefficients and fill matrix and rhs terms associated + !! with FMI object + !< subroutine fmi_fc(this, nodes, cold, nja, matrix_sln, idxglo, rhs) ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer, intent(in) :: nodes real(DP), intent(in), dimension(nodes) :: cold integer(I4B), intent(in) :: nja @@ -261,9 +276,14 @@ subroutine fmi_fc(this, nodes, cold, nja, matrix_sln, idxglo, rhs) end subroutine fmi_fc !> @brief Calculate flow correction + !! + !! Where there is a flow imbalance for a given cell, a correction may be + !! applied if selected + !< subroutine fmi_cq(this, cnew, flowja) + ! -- modules ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this real(DP), intent(in), dimension(:) :: cnew real(DP), dimension(:), contiguous, intent(inout) :: flowja ! -- local @@ -279,7 +299,7 @@ subroutine fmi_cq(this, cnew, flowja) rate = DZERO idiag = this%dis%con%ia(n) if (this%ibound(n) > 0) then - rate = -this%gwfflowja(idiag) * cnew(n) + rate = -this%gwfflowja(idiag) * cnew(n) * this%eqnsclfac end if this%flowcorrect(n) = rate flowja(idiag) = flowja(idiag) + rate @@ -290,13 +310,14 @@ subroutine fmi_cq(this, cnew, flowja) return end subroutine fmi_cq - !> @brief Calculate budget terms + !> @brief Calculate budget terms associated with FMI object + !< subroutine fmi_bd(this, isuppress_output, model_budget) ! -- modules use TdisModule, only: delt use BudgetModule, only: BudgetType, rate_accumulator ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer(I4B), intent(in) :: isuppress_output type(BudgetType), intent(inout) :: model_budget ! -- local @@ -313,10 +334,11 @@ subroutine fmi_bd(this, isuppress_output, model_budget) return end subroutine fmi_bd - !> @brief Save budget terms + !> @brief Save budget terms associated with FMI object + !< subroutine fmi_ot_flow(this, icbcfl, icbcun) ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer(I4B), intent(in) :: icbcfl integer(I4B), intent(in) :: icbcun ! -- local @@ -354,11 +376,14 @@ subroutine fmi_ot_flow(this, icbcfl, icbcun) end subroutine fmi_ot_flow !> @brief Deallocate variables + !! + !! Deallocate memory associated with FMI object + !< subroutine gwtfmi_da(this) ! -- modules use MemoryManagerModule, only: mem_deallocate ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- todo: finalize hfr and bfr either here or in a finalize routine ! ! -- deallocate any memory stored with gwfpackages @@ -405,12 +430,15 @@ subroutine gwtfmi_da(this) return end subroutine gwtfmi_da - !> @brief Allocate scalars + !> @ brief Allocate scalars + !! + !! Allocate scalar variables for an FMI object + !< subroutine gwtfmi_allocate_scalars(this) ! -- modules use MemoryManagerModule, only: mem_allocate, mem_setptr ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- local ! ! -- allocate scalars in parent @@ -430,13 +458,16 @@ subroutine gwtfmi_allocate_scalars(this) return end subroutine gwtfmi_allocate_scalars - !> @brief Allocate arrays + !> @ brief Allocate arrays for FMI object + !! + !! Method to allocate arrays for the FMI package. + !< subroutine gwtfmi_allocate_arrays(this, nodes) use MemoryManagerModule, only: mem_allocate - !modules + ! -- modules use ConstantsModule, only: DZERO ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer(I4B), intent(in) :: nodes ! -- local integer(I4B) :: n @@ -458,10 +489,15 @@ subroutine gwtfmi_allocate_arrays(this, nodes) return end subroutine gwtfmi_allocate_arrays - !> @brief Calculate groundwater cell head saturation for end of last time step + !> @brief Calculate the previous saturation level + !! + !! Calculate the groundwater cell head saturation for the end of + !! the last time step + !< function gwfsatold(this, n, delt) result(satold) + ! -- modules ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer(I4B), intent(in) :: n real(DP), intent(in) :: delt ! -- result @@ -484,13 +520,14 @@ function gwfsatold(this, n, delt) result(satold) end function gwfsatold !> @brief Read options from input file + !< subroutine gwtfmi_read_options(this) ! -- modules use ConstantsModule, only: LINELENGTH, DEM6 use InputOutputModule, only: getunit, openfile, urdaux use SimModule, only: store_error, store_error_unit ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- local character(len=LINELENGTH) :: keyword integer(I4B) :: ierr @@ -529,11 +566,14 @@ subroutine gwtfmi_read_options(this) write (this%iout, '(1x,a)') 'END OF FMI OPTIONS' end if ! - ! -- return + ! -- Return return end subroutine gwtfmi_read_options - !> @brief Read packagedata block from input file + !> @brief Read PACKAGEDATA block + !! + !! Read packagedata block from input file + !< subroutine gwtfmi_read_packagedata(this) ! -- modules use OpenSpecModule, only: ACCESS, FORM @@ -541,7 +581,7 @@ subroutine gwtfmi_read_packagedata(this) use InputOutputModule, only: getunit, openfile, urdaux use SimModule, only: store_error, store_error_unit ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- local type(BudgetObjectType), pointer :: budobjptr character(len=LINELENGTH) :: keyword, fname @@ -659,7 +699,7 @@ subroutine gwtfmi_read_packagedata(this) write (this%iout, '(1x,a)') 'END OF FMI PACKAGEDATA' end if ! - ! -- return + ! -- Return return end subroutine gwtfmi_read_packagedata @@ -669,11 +709,10 @@ end subroutine gwtfmi_read_packagedata !! pointer budget object, and this routine will look through the budget !! objects managed by FMI and point to the one with the same name, such as !! LAK-1, SFR-1, etc. - !! !< subroutine set_aptbudobj_pointer(this, name, budobjptr) ! -- modules - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- dumm character(len=*), intent(in) :: name type(BudgetObjectType), pointer :: budobjptr @@ -688,17 +727,22 @@ subroutine set_aptbudobj_pointer(this, name, budobjptr) end if end do ! - ! -- return + ! -- Return return end subroutine set_aptbudobj_pointer - !> @brief Initialize terms and count unique terms/packages in file + !> @brief Initialize the groundwater flow terms based on the budget file + !! reader + !! + !! Initalize terms and figure out how many different terms and packages + !! are contained within the file + !< subroutine initialize_gwfterms_from_bfr(this) ! -- modules use MemoryManagerModule, only: mem_allocate use SimModule, only: store_error, store_error_unit, count_errors ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- local integer(I4B) :: nflowpack integer(I4B) :: i, ip @@ -788,16 +832,19 @@ subroutine initialize_gwfterms_from_bfr(this) call this%parser%StoreErrorUnit() end if ! - ! -- return + ! -- Return return end subroutine initialize_gwfterms_from_bfr - !> @brief Initialize flow terms from a gwf-gwt exchange + !> @brief Initialize groundwater flow terms from the groundwater budget + !! + !! Flows are coming from a gwf-gwt exchange object + !< subroutine initialize_gwfterms_from_gwfbndlist(this) ! -- modules use BndModule, only: BndType, GetBndFromList ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- local integer(I4B) :: ngwfpack integer(I4B) :: ngwfterms @@ -856,21 +903,22 @@ subroutine initialize_gwfterms_from_gwfbndlist(this) iterm = iterm + 1 end if end do + ! + ! -- Return return end subroutine initialize_gwfterms_from_gwfbndlist - !> @brief Allocate GWF packages + !> @brief Initialize an array for storing PackageBudget objects. !! !! This routine allocates gwfpackages (an array of PackageBudget !! objects) to the proper size and initializes member variables. - !! !< subroutine gwtfmi_allocate_gwfpackages(this, ngwfterms) ! -- modules use ConstantsModule, only: LENMEMPATH use MemoryManagerModule, only: mem_allocate ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this integer(I4B), intent(in) :: ngwfterms ! -- local integer(I4B) :: n @@ -898,15 +946,18 @@ subroutine gwtfmi_allocate_gwfpackages(this, ngwfterms) call this%gwfpackages(n)%initialize(memPath) end do ! - ! -- return + ! -- Return return end subroutine gwtfmi_allocate_gwfpackages - !> @brief Deallocate memory in the gwfpackages array + !> @brief Deallocate memory + !! + !! Deallocate memory that stores the gwfpackages array + !< subroutine gwtfmi_deallocate_gwfpackages(this) ! -- modules ! -- dummy - class(GwtFmiType) :: this + class(TspFmiType) :: this ! -- local integer(I4B) :: n ! @@ -915,8 +966,8 @@ subroutine gwtfmi_deallocate_gwfpackages(this) call this%gwfpackages(n)%da() end do ! - ! -- return + ! -- Return return end subroutine gwtfmi_deallocate_gwfpackages -end module GwtFmiModule +end module TspFmiModule diff --git a/src/meson.build b/src/meson.build index e86adab8e17..e6f28bdb05e 100644 --- a/src/meson.build +++ b/src/meson.build @@ -98,7 +98,6 @@ modflow_sources = files( 'Model' / 'GroundWaterTransport' / 'gwt1disv1idm.f90', 'Model' / 'GroundWaterTransport' / 'gwt1dsp1.f90', 'Model' / 'GroundWaterTransport' / 'gwt1dsp1idm.f90', - 'Model' / 'GroundWaterTransport' / 'gwt1fmi1.f90', 'Model' / 'GroundWaterTransport' / 'gwt1ic1.f90', 'Model' / 'GroundWaterTransport' / 'gwt1idm.f90', 'Model' / 'GroundWaterTransport' / 'gwt1ist1.f90', @@ -134,6 +133,7 @@ modflow_sources = files( 'Model' / 'ModelUtilities' / 'Xt3dAlgorithm.f90', 'Model' / 'ModelUtilities' / 'Xt3dInterface.f90', 'Model' / 'TransportModel' / 'tsp1.f90', + 'Model' / 'TransportModel' / 'tsp1fmi1.f90', 'Model' / 'BaseModel.f90', 'Model' / 'ExplicitModel.f90', 'Model' / 'NumericalModel.f90',