refactor(gwe-est): remove packagedata block and put values in options block (#1882)

* refactor(gwe-est): remove packagedata block and put values in options block

* spelling mistake

* more spelling fixes

* addresses #1882 (comment) and #1882 (comment)

* right, don't forget the fortran

* fix bad unit specification in dfn file

* Update doc/ReleaseNotes/develop.tex

Co-authored-by: langevin-usgs <[email protected]>

* update keywords in mf6io est example tex file.  Addresses #1882 (comment)

* cleanup latent heat of vaporization; address #1882 (comment)

* forgot to run mf6ivar.py

* another unit blunder

* apparently I biffed develop.tex during a recent rebase

---------

Co-authored-by: langevin-usgs <[email protected]>

Author: emorway-usgs

autotest/test_gwe_cnd.py

@@ -307,10 +307,11 @@ def build_models(idx, test):
     flopy.mf6.ModflowGweest(
         gwe,
         save_flows=True,
+        heat_capacity_water=cpw,
+        density_water=rhow,
         porosity=prsity,
         cps=760.0,
         rhos=1500.0,
-        packagedata=[cpw, rhow, lhv],
         filename="{}.est".format(gwename),
     )
 

autotest/test_gwe_drycell_cnd0.py

@@ -289,9 +289,11 @@ def build_models(idx, test):
         gwe1,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=760.0,
         rhos=1500.0,
-        packagedata=[cpw, rhow, lhv],
         pname="EST-2",
         filename="{}.est".format(gwename1),
     )
@@ -403,7 +405,9 @@ def check_output(idx, test):
         "Pass through cell should not be as warm as its neighbor to "
         "the left"
     )
-    assert np.all(np.round(conc1[:, 0, 0, 3] - conc1[:, 0, 0, 2], 8) <= 0), msg5
+    assert np.all(
+        np.round(conc1[:, 0, 0, 3] - conc1[:, 0, 0, 2], 8) <= 0
+    ), msg5
 
 
 # - No need to change any code below

autotest/test_gwe_drycell_cnd1.py

@@ -329,9 +329,11 @@ def build_models(idx, test):
         gwe1,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=760.0,
         rhos=1500.0,
-        packagedata=[cpw, rhow, lhv],
         pname="EST-2",
         filename="{}.est".format(gwename1),
     )

autotest/test_gwe_drycell_cnd2.py

@@ -487,9 +487,11 @@ def build_models(idx, test):
         gwe,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=cps,
         rhos=rhos,
-        packagedata=[cpw, rhow, lhv],
         pname="EST-3",
         filename="{}.est".format(gwename),
     )

autotest/test_gwe_drycell_cnd4.py

@@ -348,9 +348,11 @@ def add_gwe_model(sim, gwename):
         gwe,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=cps,
         rhos=rhos,
-        packagedata=[cpw, rhow, lhv],
         pname="EST",
         filename="{}.est".format(gwename),
     )

autotest/test_gwe_esl01.py

@@ -269,9 +269,11 @@ def build_models(idx, test):
         gwe,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=cps,
         rhos=rhos,
-        packagedata=[cpw, rhow, lhv],
         pname="EST-2",
         filename="{}.est".format(gwename),
     )

autotest/test_gwe_esl02.py

@@ -268,9 +268,11 @@ def build_models(idx, test):
         gwe,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=cps,
         rhos=rhos,
-        packagedata=[cpw, rhow, lhv],
         pname="EST-1",
         filename="{}.est".format(gwename),
     )

autotest/test_gwe_esl_analyt_sln.py

@@ -144,6 +144,7 @@
 # Eqn 7-3: Bulk thermal diffusivity
 D = K_t_bulk / rho_C_bulk
 
+
 # Energy input to boundary (q_x term in the documentation)
 def calc_ener_input(primer_val):
     ener_add_rate = delr * delc * delz * rho_C_bulk * primer_val
@@ -334,9 +335,11 @@ def build_models(idx, test, ener_input):
     flopy.mf6.ModflowGweest(
         gwe,
         porosity=theta,
+        heat_capacity_water=Cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=Cps,
         rhos=rhos,
-        packagedata=[Cpw, rhow, lhv],
         pname="EST-1",
         filename="{}.est".format(gwename),
     )
@@ -412,18 +415,18 @@ def eq7_24(x, t, l, D, T_0, ener_add_rate):
     x_hat = x / l  # Dimensionless distance
 
     # Compute corresponding t_hat term
-    t_hat = D * t / l ** 2  # Dimensionless time
+    t_hat = D * t / l**2  # Dimensionless time
 
     # Solve equation 7-24
-    term1 = (1 / 2) * (x_hat ** 2 - 1 / 3)
+    term1 = (1 / 2) * (x_hat**2 - 1 / 3)
     summation_terms = [
         ((-1) ** n)
-        / n ** 2
-        * math.exp(-1 * n ** 2 * math.pi ** 2 * t_hat)
+        / n**2
+        * math.exp(-1 * n**2 * math.pi**2 * t_hat)
         * math.cos(n * math.pi * x_hat)
         for n in np.arange(1, 1000)
     ]
-    term2 = 2 / (math.pi ** 2) * np.sum(summation_terms)
+    term2 = 2 / (math.pi**2) * np.sum(summation_terms)
     T = T_0 + ener_add_rate * l / K_t_bulk * (t_hat + term1 - term2)
 
     return T
@@ -434,17 +437,17 @@ def eq7_25(x, t, l, D, T_0, ener_add_rate):
     x_hat = x / l  # Dimensionless distance
 
     # Compute corresponding t_hat term
-    t_hat = D * t / l ** 2  # Dimensionless time
+    t_hat = D * t / l**2  # Dimensionless time
 
     # Solve equation 7-25
     summation_terms = [
         ((-1) ** n)
         / (2 * n + 1) ** 2
-        * math.exp(-1 * (2 * n + 1) ** 2 * math.pi ** 2 * t_hat / 4)
+        * math.exp(-1 * (2 * n + 1) ** 2 * math.pi**2 * t_hat / 4)
         * math.sin((2 * n + 1) * math.pi * x_hat / 2)
         for n in np.arange(0, 1000)
     ]
-    term1 = (8 / math.pi ** 2) * np.sum(summation_terms)
+    term1 = (8 / math.pi**2) * np.sum(summation_terms)
 
     T = T_0 + ener_add_rate * l / K_t_bulk * (x_hat - term1)
 
@@ -456,19 +459,19 @@ def eq7_26(x, t, el, D, T_0, ener_add_rate):
     x_hat = x / el  # Dimensionless distance
 
     # Compute corresponding t_hat term
-    t_hat = D * t / el ** 2  # Dimensionless time
+    t_hat = D * t / el**2  # Dimensionless time
 
     # Solve equation 7-26
     term1 = x_hat * (1 - x_hat)
     summation_terms = [
         1
         / (2 * n + 1) ** 3
-        * math.exp(-1 * (2 * n + 1) ** 2 * math.pi ** 2 * t_hat)
+        * math.exp(-1 * (2 * n + 1) ** 2 * math.pi**2 * t_hat)
         * math.sin((2 * n + 1) * math.pi * x_hat)
         for n in np.arange(0, 1000)
     ]
-    term2 = (8 / math.pi ** 3) * np.sum(summation_terms)
-    T = T_0 + 0.5 * ener_add_rate * el ** 2 / K_t_bulk * (term1 - term2)
+    term2 = (8 / math.pi**3) * np.sum(summation_terms)
+    T = T_0 + 0.5 * ener_add_rate * el**2 / K_t_bulk * (term1 - term2)
 
     return T
 

autotest/test_gwe_lke_conduction.py

@@ -575,9 +575,11 @@ def build_models(idx, test):
     flopy.mf6.ModflowGweest(
         gwe,
         porosity=porosity,
+        heat_capacity_water=Cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=Cps,
         rhos=rhos,
-        packagedata=[Cpw, rhow, lhv],
         pname="MST-1",
         filename=f"{gwename}.mst",
     )

autotest/test_gwe_mve.py

@@ -567,9 +567,11 @@ def build_mf6_model(idx, ws):
         gwe,
         save_flows=True,
         porosity=prsity,
+        heat_capacity_water=cpw,
+        density_water=rhow,
+        latent_heat_vaporization=lhv,
         cps=cps,
         rhos=rhos,
-        packagedata=[cpw, rhow, lhv],
         pname="EST",
         filename=f"{gwename}.est",
     )