From 5ccb1164789310a3b6cdf3162512a77809fec6c3 Mon Sep 17 00:00:00 2001
From: mjreno <mjreno@IGSAAA071L01144.gs.doi.net>
Date: Mon, 27 Oct 2025 16:32:49 -0400
Subject: [PATCH] add chdg baseline

---
 flopy4/mf6/gwf/__init__.py |  4 +++-
 flopy4/mf6/gwf/chdg.py     | 48 ++++++++++++++++++++++++++++++++++++++
 test/test_component.py     | 44 +++++++++++++++++++++++++++++++++-
 3 files changed, 94 insertions(+), 2 deletions(-)
 create mode 100644 flopy4/mf6/gwf/chdg.py

diff --git a/flopy4/mf6/gwf/__init__.py b/flopy4/mf6/gwf/__init__.py
index d2051a7f..0f4495ca 100644
--- a/flopy4/mf6/gwf/__init__.py
+++ b/flopy4/mf6/gwf/__init__.py
@@ -8,6 +8,7 @@
 from xattree import xattree
 
 from flopy4.mf6.gwf.chd import Chd
+from flopy4.mf6.gwf.chdg import Chdg
 from flopy4.mf6.gwf.dis import Dis
 from flopy4.mf6.gwf.drn import Drn
 from flopy4.mf6.gwf.ic import Ic
@@ -19,7 +20,7 @@
 from flopy4.mf6.utils import open_cbc, open_hds
 from flopy4.utils import to_path
 
-__all__ = ["Gwf", "Chd", "Dis", "Drn", "Ic", "Npf", "Oc", "Wel"]
+__all__ = ["Gwf", "Chd", "Chdg", "Dis", "Drn", "Ic", "Npf", "Oc", "Wel"]
 
 
 def convert_grid(value):
@@ -76,6 +77,7 @@ def budget(self):
     oc: Oc = field(block="packages")
     npf: Npf = field(block="packages")
     chd: list[Chd] = field(block="packages")
+    chdg: list[Chdg] = field(block="packages")
     wel: list[Wel] = field(block="packages")
     drn: list[Drn] = field(block="packages")
     output: Output = attrs.field(
diff --git a/flopy4/mf6/gwf/chdg.py b/flopy4/mf6/gwf/chdg.py
new file mode 100644
index 00000000..e3278b13
--- /dev/null
+++ b/flopy4/mf6/gwf/chdg.py
@@ -0,0 +1,48 @@
+from pathlib import Path
+from typing import ClassVar, Optional
+
+import numpy as np
+from attrs import Converter
+from numpy.typing import NDArray
+from xattree import xattree
+
+from flopy4.mf6.converter import structure_array
+from flopy4.mf6.package import Package
+from flopy4.mf6.spec import array, field, path
+from flopy4.utils import to_path
+
+
+@xattree
+class Chdg(Package):
+    multi_package: ClassVar[bool] = True
+    auxiliary: Optional[list[str]] = array(block="options", default=None)
+    auxmultname: Optional[str] = field(block="options", default=None)
+    print_input: bool = field(block="options", default=False)
+    print_flows: bool = field(block="options", default=False)
+    save_flows: bool = field(block="options", default=False)
+    ts_filerecord: Optional[Path] = path(
+        block="options", default=None, converter=to_path, inout="filein"
+    )
+    obs_filerecord: Optional[Path] = path(
+        block="options", default=None, converter=to_path, inout="fileout"
+    )
+    dev_no_newton: bool = field(default=False, block="options")
+    maxbound: Optional[int] = field(block="dimensions", default=None, init=False)
+    head: Optional[NDArray[np.float64]] = array(
+        block="period",
+        dims=(
+            "nper",
+            "nodes",
+        ),
+        default=None,
+        converter=Converter(structure_array, takes_self=True, takes_field=True),
+    )
+    aux: Optional[NDArray[np.float64]] = array(
+        block="period",
+        dims=(
+            "nper",
+            "nodes",
+        ),
+        default=None,
+        converter=Converter(structure_array, takes_self=True, takes_field=True),
+    )
diff --git a/test/test_component.py b/test/test_component.py
index 0be0d31a..2e971d32 100644
--- a/test/test_component.py
+++ b/test/test_component.py
@@ -8,11 +8,13 @@
 
 from flopy4.mf6.component import COMPONENTS
 from flopy4.mf6.constants import FILL_DNODATA, LENBOUNDNAME
-from flopy4.mf6.gwf import Chd, Dis, Gwf, Ic, Npf, Oc
+from flopy4.mf6.gwf import Chd, Chdg, Dis, Gwf, Ic, Npf, Oc
 from flopy4.mf6.ims import Ims
 from flopy4.mf6.simulation import Simulation
 from flopy4.mf6.tdis import Tdis
 
+DNODATA = 3.0e30
+
 
 def test_registry():
     assert COMPONENTS["simulation"] is Simulation
@@ -394,6 +396,46 @@ def test_quickstart(function_tmpdir):
     sim.run()
 
 
+def test_quickstart_array(function_tmpdir):
+    sim_name = "quickstart"
+    gwf_name = "mymodel"
+    nlay = 1
+    nrow = 10
+    ncol = 10
+    time = ModelTime(perlen=[1.0], nstp=[1], tsmult=[1.0])
+    ims = Ims(filename="mymodel.ims", models=[gwf_name])
+    dis = Dis(
+        nlay=nlay,
+        nrow=nrow,
+        ncol=ncol,
+        top=1.0,
+        botm=0.0,
+    )
+    sim = Simulation(
+        tdis=time,
+        workspace=function_tmpdir,
+        name=sim_name,
+        solutions={"ims": ims},
+    )
+    gwf = Gwf(parent=sim, dis=dis, name=gwf_name)
+    ic = Ic(parent=gwf)
+    oc = Oc(
+        parent=gwf,
+        budget_file=f"{gwf_name}.bud",
+        head_file=f"{gwf_name}.hds",
+        save_head=["all"],
+        save_budget=["all"],
+    )
+    npf = Npf(parent=gwf, icelltype=0, k=1.0)
+    head = {0: np.full((nlay, nrow, ncol), DNODATA, dtype=float)}
+    head[0][0, 0, 0] = 1.0
+    head[0][0, 9, 9] = 0.0
+    chd = Chdg(parent=gwf, head=head)
+
+    sim.write()
+    sim.run()
+
+
 def test_write_ascii(function_tmpdir):
     sim_name = "sim"
     gwf_name = "gwf"