refactor(swf): minor swf cleanup (#2084)

* refactor(swf): minor swf cleanup

* more ruffing

autotest/test_chf_dfw_loop.py

@@ -258,6 +258,16 @@ def build_models(idx, test):
         ],
     )
 
+    # note: for specifying zero-based reach number, put reach number in tuple
+    fname = f"{name}.flw.obs.csv"
+    flw_obs = {
+        fname: [
+            ("INFLOW1", "FLW", (0,)),
+            ("INFLOW5", "FLW", (4,)),
+        ],
+        # "digits": 10,
+    }
+
     flwlist = [
         [(0,), "reach1"],
         [(4,), "reach5"],
@@ -267,6 +277,7 @@ def build_models(idx, test):
         maxbound=len(flwlist),
         print_input=True,
         print_flows=True,
+        observations=flw_obs,
         stress_period_data=flwlist,
     )
 
@@ -317,65 +328,112 @@ def make_plot(test):
 
     name = test.name
     ws = test.workspace
-    mfsim = flopy.mf6.MFSimulation.load(sim_ws=ws)
+    mfsim = test.sims[0]
+
     chf = mfsim.get_model(name)
 
-    fpth = test.workspace / f"{name}.obs.csv"
-    obsvals = np.genfromtxt(fpth, names=True, delimiter=",")
+    # observations for flow package
+    output = chf.obs[0].output
+    obs_names = output.obs_names
+    print("Obs output names: ", obs_names)
+    print("Observations", output.methods())
+    flw_obsvals = output.obs().get_data()
+    print(flw_obsvals)
 
-    fig = plt.figure(figsize=(10, 10))
+    fig = plt.figure(figsize=(4, 4))
     ax = fig.add_subplot(1, 1, 1)
-    for irch in [1, 4, 15, 18]:
+    colors = ["blue", "red"]
+    sec2hour = 1.0 / 60.0 / 60.0
+    for i, irch in enumerate([1, 5]):
         ax.plot(
-            obsvals["time"],
+            flw_obsvals["totim"] * sec2hour,
+            flw_obsvals[f"INFLOW{irch}"],
+            color=colors[i],
+            marker="o",
+            mfc="none",
+            mec="none",
+            lw=0.5,
+            label=f"Inflow reach {irch}",
+        )
+    plt.xlabel("time, in hours")
+    plt.ylabel("inflow, in cubic meters per second")
+    plt.legend()
+    ax.set_xlim(0.0, 25.0)
+    fname = ws / "loop_network_inflow.png"
+    plt.savefig(fname)
+
+    # observations for model
+    output = chf.obs[1].output
+    obs_names = output.obs_names
+    print("Obs output names: ", obs_names)
+    print("observations", output.methods())
+    obsvals = output.obs().get_data()
+
+    fig = plt.figure(figsize=(6, 4))
+    ax = fig.add_subplot(1, 1, 1)
+    colors = ["b", "g", "r", "c"]
+    for i, irch in enumerate([1, 4, 15, 18]):
+        ax.plot(
+            obsvals["totim"] * sec2hour,
             obsvals[f"REACH{irch}"],
             marker="o",
             mfc="none",
-            mec="k",
+            mec=colors[i],
             lw=0.0,
             label=f"MF6 reach {irch}",
         )
         ax.plot(
-            obsvals["time"],
+            obsvals["totim"] * sec2hour,
             answer[f"STAGE00000000{irch:02d}"],
-            "k-",
+            f"{colors[i]}-",
             label=f"SWR Reach {irch}",
         )
     ax.set_xscale("log")
-    plt.xlabel("time, in seconds")
+    plt.xlabel("time, in hours")
     plt.ylabel("stage, in meters")
     plt.legend()
-    fname = ws / f"{name}.obs.1.png"
+    fname = ws / "loop_network_stage.png"
     plt.savefig(fname)
 
-    fig = plt.figure(figsize=(10, 10))
+    fig = plt.figure(figsize=(6, 4))
     ax = fig.add_subplot(1, 1, 1)
+    sign = -1.0
     ax.plot(
-        obsvals["time"],
-        obsvals["FLOW45"],
+        obsvals["totim"] * sec2hour,
+        obsvals["FLOW45"] * sign,
         marker="o",
         mfc="none",
         mec="b",
         lw=0.0,
-        label="MF6 Gauge 4",
+        label="MF6 Location 4",
     )
     ax.plot(
-        obsvals["time"],
-        obsvals["FLOW56"],
+        obsvals["totim"] * sec2hour,
+        obsvals["FLOW56"] * sign,
         marker="o",
         mfc="none",
         mec="g",
         lw=0.0,
-        label="MF6 Gauge 5",
+        label="MF6 Location 5",
+    )
+    ax.plot(
+        answer_flow["TOTIME"] * sec2hour,
+        answer_flow["FLOW45"] * sign,
+        "b-",
+        label="SWR Location 4",
+    )
+    ax.plot(
+        answer_flow["TOTIME"] * sec2hour,
+        answer_flow["FLOW56"] * sign,
+        "g-",
+        label="SWR Location 5",
     )
-    ax.plot(answer_flow["TOTIME"], answer_flow["FLOW45"], "b-", label="SWR Gauge 4")
-    ax.plot(answer_flow["TOTIME"], answer_flow["FLOW56"], "g-", label="SWR Gauge 5")
     # ax.plot(obsvals["time"], answer["STAGE0000000014"], marker="o", mfc="none", mec="k", lw=0., label="swr")  # noqa
     ax.set_xscale("log")
-    plt.xlabel("time, in seconds")
-    plt.ylabel("flow, in cubic meters per second")
+    plt.xlabel("time, in hours")
+    plt.ylabel("reach outflow, in cubic meters per second")
     plt.legend()
-    fname = ws / f"{name}.obs.2.png"
+    fname = ws / "loop_network_flow.png"
     plt.savefig(fname)
 
     fig = plt.figure(figsize=(6, 6))
@@ -388,6 +446,26 @@ def make_plot(test):
     ax.plot(vertices["xv"], vertices["yv"], "bo")
     for iv, x, y in vertices:
         ax.text(x, y, f"{iv + 1}")
+    cell1d = chf.disv1d.cell1d.get_data()
+    print(cell1d.dtype)
+    for row in cell1d:
+        ic = row["icell1d"]
+        ncvert = row["ncvert"]
+        if ncvert == 2:
+            n0 = row["icvert_0"]
+            n1 = row["icvert_1"]
+            x0 = vertices["xv"][n0]
+            x1 = vertices["xv"][n1]
+            y0 = vertices["yv"][n0]
+            y1 = vertices["yv"][n1]
+            xm = 0.5 * (x0 + x1)
+            ym = 0.5 * (y0 + y1)
+        elif ncvert == 3:
+            n0 = row["icvert_1"]
+            xm = vertices["xv"][n0] + 150
+            ym = vertices["yv"][n0] + 150
+        ax.text(xm, ym, f"{ic + 1}", color="red")
+    # raise Exception()
     fname = ws / "grid.png"
     plt.savefig(fname)
 

doc/mf6io/swf/disv1d.tex

@@ -6,7 +6,7 @@ \subsubsection{Structure of Blocks}
 \lstinputlisting[style=blockdefinition]{./mf6ivar/tex/swf-disv1d-dimensions.dat}
 \lstinputlisting[style=blockdefinition]{./mf6ivar/tex/swf-disv1d-griddata.dat}
 \lstinputlisting[style=blockdefinition]{./mf6ivar/tex/swf-disv1d-vertices.dat}
-\lstinputlisting[style=blockdefinition]{./mf6ivar/tex/swf-disv1d-cell2d.dat}
+\lstinputlisting[style=blockdefinition]{./mf6ivar/tex/swf-disv1d-cell1d.dat}
 
 \vspace{5mm}
 \subsubsection{Explanation of Variables}

doc/mf6io/swf/swf.tex

@@ -66,8 +66,4 @@ \subsection{Critical Depth Boundary (CDB) Package}
 \subsection{Observation (OBS) Utility for a SWF Model}
 \input{swf/swf-obs}
 
-\newpage
-\subsection{Inflow (FLW) Package}
-\input{swf/flw}
-
 

src/Model/GroundWaterFlow/gwf-buy.f90

@@ -22,31 +22,31 @@ module GwfBuyModule
   public :: buy_cr
 
   type :: ConcentrationPointer
-    real(DP), dimension(:), pointer :: conc => null() ! pointer to concentration array
-    integer(I4B), dimension(:), pointer :: icbund => null() ! store pointer to gwt ibound array
+    real(DP), dimension(:), pointer :: conc => null() !< pointer to concentration array
+    integer(I4B), dimension(:), pointer :: icbund => null() !< store pointer to gwt ibound array
   end type ConcentrationPointer
 
   type, extends(NumericalPackageType) :: GwfBuyType
-    type(GwfNpfType), pointer :: npf => null() ! npf object
-    integer(I4B), pointer :: ioutdense => null() ! unit number for saving density
-    integer(I4B), pointer :: iform => null() ! formulation: 0 freshwater head, 1 hh rhs, 2 hydraulic head
-    integer(I4B), pointer :: ireadelev => null() ! if 1 then elev has been allocated and filled
-    integer(I4B), pointer :: ireadconcbuy => null() ! if 1 then dense has been read from this buy input file
-    integer(I4B), pointer :: iconcset => null() ! if 1 then conc is pointed to a gwt model%x
-    real(DP), pointer :: denseref => null() ! reference fluid density
-    real(DP), dimension(:), pointer, contiguous :: dense => null() ! density
-    real(DP), dimension(:), pointer, contiguous :: concbuy => null() ! concentration array if specified in buy package
-    real(DP), dimension(:), pointer, contiguous :: elev => null() ! cell center elevation (optional; if not specified, then use (top+bot)/2)
-    integer(I4B), dimension(:), pointer :: ibound => null() ! store pointer to ibound
+    type(GwfNpfType), pointer :: npf => null() !< npf object
+    integer(I4B), pointer :: ioutdense => null() !< unit number for saving density
+    integer(I4B), pointer :: iform => null() !< formulation: 0 freshwater head, 1 hh rhs, 2 hydraulic head
+    integer(I4B), pointer :: ireadelev => null() !< if 1 then elev has been allocated and filled
+    integer(I4B), pointer :: ireadconcbuy => null() !< if 1 then dense has been read from this buy input file
+    integer(I4B), pointer :: iconcset => null() !< if 1 then conc is pointed to a gwt model%x
+    real(DP), pointer :: denseref => null() !< reference fluid density
+    real(DP), dimension(:), pointer, contiguous :: dense => null() !< density
+    real(DP), dimension(:), pointer, contiguous :: concbuy => null() !< concentration array if specified in buy package
+    real(DP), dimension(:), pointer, contiguous :: elev => null() !< cell center elevation (optional; if not specified, then use (top+bot)/2)
+    integer(I4B), dimension(:), pointer :: ibound => null() !< store pointer to ibound
 
-    integer(I4B), pointer :: nrhospecies => null() ! number of species used in equation of state to calculate density
-    real(DP), dimension(:), pointer, contiguous :: drhodc => null() ! change in density with change in concentration
-    real(DP), dimension(:), pointer, contiguous :: crhoref => null() ! reference concentration used in equation of state
-    real(DP), dimension(:), pointer, contiguous :: ctemp => null() ! temporary array of size (nrhospec) to pass to calcdens
-    character(len=LENMODELNAME), dimension(:), allocatable :: cmodelname ! names of gwt models used in equation of state
-    character(len=LENAUXNAME), dimension(:), allocatable :: cauxspeciesname ! names of gwt models used in equation of state
+    integer(I4B), pointer :: nrhospecies => null() !< number of species used in equation of state to calculate density
+    real(DP), dimension(:), pointer, contiguous :: drhodc => null() !< change in density with change in concentration
+    real(DP), dimension(:), pointer, contiguous :: crhoref => null() !< reference concentration used in equation of state
+    real(DP), dimension(:), pointer, contiguous :: ctemp => null() !< temporary array of size (nrhospec) to pass to calcdens
+    character(len=LENMODELNAME), dimension(:), allocatable :: cmodelname !< names of gwt models used in equation of state
+    character(len=LENAUXNAME), dimension(:), allocatable :: cauxspeciesname !< names of gwt models used in equation of state
 
-    type(ConcentrationPointer), allocatable, dimension(:) :: modelconc ! concentration pointer for each transport model
+    type(ConcentrationPointer), allocatable, dimension(:) :: modelconc !< concentration pointer for each transport model
 
   contains
     procedure :: buy_df

src/Model/SurfaceWaterFlow/swf-cdb.f90

@@ -15,12 +15,9 @@ module SwfCdbModule
   use BndModule, only: BndType
   use BndExtModule, only: BndExtType
   use ObsModule, only: DefaultObsIdProcessor
-  use SmoothingModule, only: sQSaturation, sQSaturationDerivative
   use ObserveModule, only: ObserveType
   use TimeSeriesLinkModule, only: TimeSeriesLinkType, &
                                   GetTimeSeriesLinkFromList
-  use BlockParserModule, only: BlockParserType
-  use InputOutputModule, only: GetUnit, openfile
   use MatrixBaseModule
   use BaseDisModule, only: DisBaseType
   use SwfCxsModule, only: SwfCxsType

src/Model/SurfaceWaterFlow/swf-cxs.f90

@@ -8,8 +8,8 @@ module SwfCxsModule
   use ConstantsModule, only: LENMEMPATH, DZERO, DTWOTHIRDS
   use MemoryHelperModule, only: create_mem_path
   use MemoryManagerModule, only: mem_allocate
-  use SimVariablesModule, only: errmsg, warnmsg
-  use SimModule, only: count_errors, store_error, store_error_unit
+  use SimVariablesModule, only: errmsg
+  use SimModule, only: store_error
   use NumericalPackageModule, only: NumericalPackageType
   use BaseDisModule, only: DisBaseType
 

src/Model/SurfaceWaterFlow/swf-dfw.f90

@@ -5,21 +5,20 @@
 module SwfDfwModule
 
   use KindModule, only: DP, I4B, LGP
-  use ConstantsModule, only: LENMEMPATH, LENVARNAME, LINELENGTH, &
-                             DZERO, DHALF, DONE, DTWO, DTHREE, &
-                             DNODATA, DEM5, DTWOTHIRDS, DP9, DONETHIRD, &
+  use ConstantsModule, only: LENMEMPATH, LINELENGTH, &
+                             DZERO, DHALF, DONE, DTWO, &
+                             DTWOTHIRDS, DP9, DONETHIRD, &
                              DPREC, DEM10
   use MemoryHelperModule, only: create_mem_path
-  use MemoryManagerModule, only: mem_reallocate
-  use MemoryManagerModule, only: mem_allocate, mem_setptr, get_isize
-  use SimVariablesModule, only: errmsg, warnmsg
+  use MemoryManagerModule, only: mem_allocate, mem_setptr, get_isize, &
+                                 mem_reallocate
+  use SimVariablesModule, only: errmsg
   use SimModule, only: count_errors, store_error, store_error_unit, &
                        store_error_filename
   use NumericalPackageModule, only: NumericalPackageType
   use BaseDisModule, only: DisBaseType
   use SwfCxsModule, only: SwfCxsType
   use ObsModule, only: ObsType, obs_cr
-  use ObsModule, only: DefaultObsIdProcessor
   use ObserveModule, only: ObserveType
   use MatrixBaseModule
 

src/Model/SurfaceWaterFlow/swf-flw.f90

@@ -7,20 +7,15 @@
 module SwfFlwModule
   ! -- modules
   use KindModule, only: DP, I4B
-  use ConstantsModule, only: DZERO, DEM1, DONE, LENFTYPE, DNODATA, &
-                             LINELENGTH
+  use ConstantsModule, only: DZERO, LENFTYPE, DNODATA
   use SimVariablesModule, only: errmsg
   use SimModule, only: store_error, store_error_filename
-  use MemoryHelperModule, only: create_mem_path
   use BndModule, only: BndType
   use BndExtModule, only: BndExtType
   use ObsModule, only: DefaultObsIdProcessor
-  use SmoothingModule, only: sQSaturation, sQSaturationDerivative
   use ObserveModule, only: ObserveType
   use TimeSeriesLinkModule, only: TimeSeriesLinkType, &
                                   GetTimeSeriesLinkFromList
-  use BlockParserModule, only: BlockParserType
-  use InputOutputModule, only: GetUnit, openfile
   use MatrixBaseModule
   !
   implicit none

src/Model/SurfaceWaterFlow/swf-ic.f90

@@ -3,7 +3,6 @@ module SwfIcModule
   use KindModule, only: DP, I4B, LGP
   use ConstantsModule, only: LINELENGTH
   use NumericalPackageModule, only: NumericalPackageType
-  use BlockParserModule, only: BlockParserType
   use BaseDisModule, only: DisBaseType
 
   implicit none
@@ -71,7 +70,6 @@ end subroutine ic_cr
   !<
   subroutine ic_load(this)
     ! -- modules
-    use BaseDisModule, only: DisBaseType
     ! -- dummy
     class(SwfIcType) :: this
     !

src/Model/SurfaceWaterFlow/swf-obs.f90

@@ -1,7 +1,7 @@
 module SwfObsModule
 
   use KindModule, only: DP, I4B
-  use ConstantsModule, only: LINELENGTH, MAXOBSTYPES
+  use ConstantsModule, only: LINELENGTH
   use BaseDisModule, only: DisBaseType
   use SwfIcModule, only: SwfIcType
   use ObserveModule, only: ObserveType

src/Model/SurfaceWaterFlow/swf-oc.f90

@@ -2,7 +2,6 @@ module SwfOcModule
 
   use BaseDisModule, only: DisBaseType
   use KindModule, only: DP, I4B
-  use ConstantsModule, only: LENMODELNAME
   use OutputControlModule, only: OutputControlType
   use OutputControlDataModule, only: OutputControlDataType, ocd_cr
 

src/Model/SurfaceWaterFlow/swf-sto.f90

@@ -7,14 +7,12 @@
 module SwfStoModule
 
   use KindModule, only: DP, I4B, LGP
-  use ConstantsModule, only: DZERO, DEM6, DEM4, DHALF, DONE, DTWO, &
-                             LENBUDTXT, LINELENGTH, LENMEMPATH
+  use ConstantsModule, only: DZERO, LENBUDTXT, LINELENGTH, LENMEMPATH
   use MemoryHelperModule, only: create_mem_path
   use SimVariablesModule, only: errmsg
-  use SimModule, only: store_error, store_error_filename, count_errors
+  use SimModule, only: store_error, store_error_filename
   use BaseDisModule, only: DisBaseType
   use NumericalPackageModule, only: NumericalPackageType
-  use InputOutputModule, only: GetUnit, openfile
   use MatrixBaseModule
   use Disv1dModule, only: Disv1dType
   use SwfCxsModule, only: SwfCxsType
@@ -576,7 +574,6 @@ end subroutine allocate_arrays
   !<
   subroutine source_options(this)
     ! -- modules
-    use ConstantsModule, only: LENMEMPATH
     use MemoryManagerExtModule, only: mem_set_value
     use SourceCommonModule, only: filein_fname
     use SwfStoInputModule, only: SwfStoParamFoundType