diff --git a/autotest/test_prt_exg.py b/autotest/test_prt_exg.py index 04caca2b396..6ebb0a54937 100644 --- a/autotest/test_prt_exg.py +++ b/autotest/test_prt_exg.py @@ -8,8 +8,15 @@ Results are compared against a MODPATH 7 model. -This test includes two cases, one which gives -boundnames to particles and one which does not. +This test includes four cases: one which gives +boundnames to particles, one which does not, a +third in which the flow model has a uniformly +active idomain array while the tracking model +does not, and a final case in which flow and +tracking model have different IDOMAIN arrays, +both non-uniform, where the active region is +the same size but consists of different cells. +Both latter cases should be caught as errors. """ from pathlib import Path @@ -30,7 +37,7 @@ ) simname = "prtexg01" -cases = [simname, f"{simname}bnms"] +cases = [simname, f"{simname}bnms", f"{simname}idmu", f"{simname}idmn"] def get_model_name(idx, mdl): @@ -41,18 +48,31 @@ def build_mf6_sim(idx, test): # create simulation name = cases[idx] sim = FlopyReadmeCase.get_gwf_sim(name, test.workspace, test.targets["mf6"]) + gwf = sim.get_model() # create prt model prt_name = get_model_name(idx, "prt") prt = flopy.mf6.ModflowPrt(sim, modelname=prt_name) # create prt discretization + idomain = np.ones( + (FlopyReadmeCase.nlay, FlopyReadmeCase.nrow, FlopyReadmeCase.ncol) + ) + if "idm" in name: + idomain[-1, -1, -1] = 0 + if "idmn" in name: + # TODO alter flow model's idomain + gwf_idomain = idomain.copy() + gwf_idomain[-1, -1, -1] = 1 + gwf_idomain[0, 0, 0] = 0 + gwf.dis.idomain = gwf_idomain flopy.mf6.modflow.mfgwfdis.ModflowGwfdis( prt, pname="dis", nlay=FlopyReadmeCase.nlay, nrow=FlopyReadmeCase.nrow, ncol=FlopyReadmeCase.ncol, + idomain=idomain, ) # create mip package @@ -155,7 +175,7 @@ def build_models(idx, test): gwf_name = get_model_name(idx, "gwf") gwf = mf6sim.get_model(gwf_name) mp7sim = build_mp7_sim(idx, test.workspace / "mp7", test.targets["mp7"], gwf) - return mf6sim, mp7sim + return mf6sim, None if "idm" in test.name else mp7sim def check_output(idx, test): @@ -178,6 +198,9 @@ def check_output(idx, test): # extract model grid mg = gwf.modelgrid + if "idm" in name: + return + # check mf6 output files exist gwf_budget_file = f"{gwf_name}.bud" gwf_head_file = f"{gwf_name}.hds" @@ -362,5 +385,6 @@ def test_mf6model(idx, name, function_tmpdir, targets, plot): plot=lambda t: plot_output(idx, t) if plot else None, targets=targets, compare=None, + xfail="idm" in name, ) test.run() diff --git a/src/Exchange/exg-gwfgwe.f90 b/src/Exchange/exg-gwfgwe.f90 index 156abe01a74..536a57b232f 100644 --- a/src/Exchange/exg-gwfgwe.f90 +++ b/src/Exchange/exg-gwfgwe.f90 @@ -181,6 +181,9 @@ end subroutine exg_df subroutine exg_ar(this) ! -- modules use MemoryManagerModule, only: mem_checkin + use DisModule, only: DisType + use DisvModule, only: DisvType + use DisuModule, only: DisuType ! -- dummy class(GwfGweExchangeType) :: this ! -- local @@ -194,6 +197,11 @@ subroutine exg_ar(this) & GWF Model has ', i0, ' user nodes and ', i0, ' reduced nodes.& & GWE Model has ', i0, ' user nodes and ', i0, ' reduced nodes.& & Ensure discretization packages, including IDOMAIN, are identical.')" + character(len=*), parameter :: fmtidomerr = & + "('GWF and GWE Models do not have the same discretization for exchange& + & ',a,'.& + & GWF Model and GWE Model have different IDOMAIN arrays.& + & Ensure discretization packages, including IDOMAIN, are identical.')" ! ! -- set gwfmodel mb => GetBaseModelFromList(basemodellist, this%m1_idx) @@ -220,6 +228,34 @@ subroutine exg_ar(this) call store_error(errmsg, terminate=.TRUE.) end if ! + ! -- Make sure idomains are identical + select type (gwfdis => gwfmodel%dis) + type is (DisType) + select type (gwedis => gwemodel%dis) + type is (DisType) + if (.not. all(gwfdis%idomain == gwedis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + type is (DisvType) + select type (gwedis => gwemodel%dis) + type is (DisvType) + if (.not. all(gwfdis%idomain == gwedis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + type is (DisuType) + select type (gwedis => gwemodel%dis) + type is (DisuType) + if (.not. all(gwfdis%idomain == gwedis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + end select + ! ! -- setup pointers to gwf variables allocated in gwf_ar gwemodel%fmi%gwfhead => gwfmodel%x call mem_checkin(gwemodel%fmi%gwfhead, & diff --git a/src/Exchange/exg-gwfgwt.f90 b/src/Exchange/exg-gwfgwt.f90 index d96c1db6ec3..aa198e7a100 100644 --- a/src/Exchange/exg-gwfgwt.f90 +++ b/src/Exchange/exg-gwfgwt.f90 @@ -184,6 +184,9 @@ end subroutine exg_df subroutine exg_ar(this) ! -- modules use MemoryManagerModule, only: mem_checkin + use DisModule, only: DisType + use DisvModule, only: DisvType + use DisuModule, only: DisuType ! -- dummy class(GwfGwtExchangeType) :: this ! -- local @@ -197,6 +200,11 @@ subroutine exg_ar(this) & GWF Model has ', i0, ' user nodes and ', i0, ' reduced nodes.& & GWT Model has ', i0, ' user nodes and ', i0, ' reduced nodes.& & Ensure discretization packages, including IDOMAIN, are identical.')" + character(len=*), parameter :: fmtidomerr = & + "('GWF and GWT Models do not have the same discretization for exchange& + & ',a,'.& + & GWF Model and GWT Model have different IDOMAIN arrays.& + & Ensure discretization packages, including IDOMAIN, are identical.')" ! ! -- set gwfmodel mb => GetBaseModelFromList(basemodellist, this%m1_idx) @@ -223,6 +231,34 @@ subroutine exg_ar(this) call store_error(errmsg, terminate=.TRUE.) end if ! + ! -- Make sure idomains are identical + select type (gwfdis => gwfmodel%dis) + type is (DisType) + select type (gwtdis => gwtmodel%dis) + type is (DisType) + if (.not. all(gwfdis%idomain == gwtdis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + type is (DisvType) + select type (gwtdis => gwtmodel%dis) + type is (DisvType) + if (.not. all(gwfdis%idomain == gwtdis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + type is (DisuType) + select type (gwtdis => gwtmodel%dis) + type is (DisuType) + if (.not. all(gwfdis%idomain == gwtdis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + end select + ! ! -- setup pointers to gwf variables allocated in gwf_ar gwtmodel%fmi%gwfhead => gwfmodel%x call mem_checkin(gwtmodel%fmi%gwfhead, & diff --git a/src/Exchange/exg-gwfprt.f90 b/src/Exchange/exg-gwfprt.f90 index 87393123634..61929f13693 100644 --- a/src/Exchange/exg-gwfprt.f90 +++ b/src/Exchange/exg-gwfprt.f90 @@ -177,6 +177,9 @@ end subroutine exg_df subroutine exg_ar(this) ! -- modules use MemoryManagerModule, only: mem_checkin + use DisModule, only: DisType + use DisvModule, only: DisvType + use DisuModule, only: DisuType ! -- dummy class(GwfPrtExchangeType) :: this ! -- local @@ -190,6 +193,11 @@ subroutine exg_ar(this) & GWF Model has ', i0, ' user nodes and ', i0, ' reduced nodes.& & PRT Model has ', i0, ' user nodes and ', i0, ' reduced nodes.& & Ensure discretization packages, including IDOMAIN, are identical.')" + character(len=*), parameter :: fmtidomerr = & + "('GWF and PRT Models do not have the same discretization for exchange& + & ',a,'.& + & GWF Model and PRT Model have different IDOMAIN arrays.& + & Ensure discretization packages, including IDOMAIN, are identical.')" ! ! -- set gwfmodel mb => GetBaseModelFromList(basemodellist, this%m1id) @@ -216,6 +224,34 @@ subroutine exg_ar(this) call store_error(errmsg, terminate=.TRUE.) end if ! + ! -- Make sure idomains are identical + select type (gwfdis => gwfmodel%dis) + type is (DisType) + select type (prtdis => prtmodel%dis) + type is (DisType) + if (.not. all(gwfdis%idomain == prtdis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + type is (DisvType) + select type (prtdis => prtmodel%dis) + type is (DisvType) + if (.not. all(gwfdis%idomain == prtdis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + type is (DisuType) + select type (prtdis => prtmodel%dis) + type is (DisuType) + if (.not. all(gwfdis%idomain == prtdis%idomain)) then + write (errmsg, fmtidomerr) trim(this%name) + call store_error(errmsg, terminate=.TRUE.) + end if + end select + end select + ! ! -- setup pointers to gwf variables allocated in gwf_ar prtmodel%fmi%gwfhead => gwfmodel%x call mem_checkin(prtmodel%fmi%gwfhead, & diff --git a/src/Model/ModelUtilities/FlowModelInterface.f90 b/src/Model/ModelUtilities/FlowModelInterface.f90 index 60a4729e257..7437722f307 100644 --- a/src/Model/ModelUtilities/FlowModelInterface.f90 +++ b/src/Model/ModelUtilities/FlowModelInterface.f90 @@ -23,7 +23,7 @@ module FlowModelInterfaceModule logical, pointer :: flows_from_file => null() !< if .false., then flows come from GWF through GWF-Model exg type(ListType), pointer :: gwfbndlist => null() !< list of gwf stress packages integer(I4B), pointer :: iflowsupdated => null() !< flows were updated for this time step - integer(I4B), dimension(:), pointer, contiguous :: ibound => null() !< pointer to Model ibound + integer(I4B), dimension(:), pointer, contiguous :: ibound => null() !< pointer to model ibound real(DP), dimension(:), pointer, contiguous :: gwfflowja => null() !< pointer to the GWF flowja array real(DP), dimension(:, :), pointer, contiguous :: gwfspdis => null() !< pointer to npf specific discharge array real(DP), dimension(:), pointer, contiguous :: gwfhead => null() !< pointer to the GWF head array