Analytic Element Domain Boundary Conditions for Site-Scale Groundwater Flow Modeling Los Angeles Basin

Physics-based groundwater flow modeling is a useful tool for the design and optimization of pump-and-treat systems for groundwater site cleanup. Numerical methods like finite differences and finite elements, and hybrid analytic elements, require boundary conditions (BC) to be assigned to the outer domain of the grid, mesh, or line elements. These outer BC do not always correspond with hydrogeologic features. Common practice in model setup is to either: (1) extend the model domain boundary outward such that introduced artificial outer BCs (e.g., first type head specified, second type flux specified) do not have undue influence on near-field scale simulations; or (2) assign outer BCs to capture the effective far-field influence (e.g., third type head-dependent flux). Groundwater flow modeling options for assigning BCs were demonstrated for the extensively documented Dual Site Superfund cleanup in Torrance, California. The existing MODFLOW models for the Dual Site scale and the Los Angeles basin scale document the current hydrogeologic conceptual site model. Simplified analytic element AnAqSim models at the LA Basin scale, West Coast Subbasin scale, and Dual Site scale, were used for mapping near-field domain velocity vector fields and pathline envelopes. The pump-treat-inject system demonstrated hydraulic containment and showed pathline envelopes relatively insensitive to BC choices. However, the nearfield domain boundary groundwater flow fields were sensitive to BC choices. The Los Angeles basin case study demonstrated the use of analytic element groundwater modeling for testing stress dependent boundaries during site pump-treat-inject design.

Groundwater flow model BCs are classified as three main types: (Type 1) Dirichlet specified head including constant head; (Type 2) Neumann specified flux including no-flow; and (Type 3) Robin head dependent flux.Jazayeri and Werner (2019) frame the boundary problem for groundwater flow systems within a closed region Ω with boundary Γ and x representing temporal or spatial dimensions, and ϕ (representing hydraulic head) and n (representing the normal direction) and a (representing a non-zero coefficient).BCs are referenced by three main types as summarized in Table S1.
Table S1.Groundwater model boundary condition types (BC) after Jazayeri, Werner (2019).The Cauchy (Type 4) BC has Dirichlet and Neumann BC independently specified along the same boundary but is rarely used in real world applications.Jazayeri and Werner (2019) document inconsistencies in the literature and found most references to Cauchy (Type 4) BC are misidentified.
Analytical solutions tend to have unbounded (infinite) domains.
A seepage face is an example of a switching BC between Dirichlet and Neumann.
Stress-dependency is of primary concern wherever the model boundaries differ from the natural system boundaries, or the natural boundaries that may extend beyond the boundaries of the model (ASTM 2016).If sensitivity tests indicate stress dependency of a specified head or flux boundary, then the head-dependent flux boundary is suggested.
The practical use of the Robin head dependent flux (Type 3) BC is implemented in MODFLOW as a general head boundary (GHB).The head dependent flux at the boundary cell is defined by a far-field constant head (effectively infinite source of water, hypothetical or real) a set distance away and the material of transfer has an assigned effective conductance (Fig. S1a).The head dependent normal flux (Type 3) BC is represented in AnAqSim as an outer domain line element with an estimated effective outside domain conductance and head (Fig. S1b).

Final 26Apr2023
Disclaimer: Supporting Information is not peer reviewed.

Final 26Apr2023
Disclaimer: Supporting Information is not peer reviewed.
The pathlines associated with the Dual Site pump-treat-inject system are visualized using the MODFLOW solution for groundwater flow and the MODPATH model for reverse particle tracking from the extraction wells (Fig. S3).Attention focused on the Gage aquifer layer that is the target for the injection wells and the model suggests that the injection wells and extraction wells form a closed pathline envelope in the Gage aquifer.The model files come from ICF/Sundance (2019).

Final 26Apr2023
Disclaimer: Supporting Information is not peer reviewed.

Final 26Apr2023
Disclaimer: Supporting Information is not peer reviewed.The Los Angeles Basin including the Central, West Coast, and Orange County Basins, as represented in the single level steady GFLOW groundwater model, has outer boundary conditions having mostly geohydrological correspondence: no-flow boundaries associated with rock hills; specified head boundaries associated with the ocean shorelines; specified flux focused through the LA Narrows and Whittier Narrows.The layout of elements is shown in Fig. S4.The calibration of GFLOW for "predevelopment" conditions are shown in Fig. S5 and S6.In Fig. S5, the basemap is the map of head contours reported by Mendenhall (1905).The shaded areas were artesian at that time prior to over pumping.GFLOW test points were associated with the 10, 20, 30 foot contours in the West Coast Basin.GFLOW shows triangles centered on test points --red triangles with tip pointing down show model prediction less than observed, green triangles with tip pointing up show model prediction greater than observed.Ideally the distribution of triangles in space and in size should be random to minimize bias.
Figure S5.The GFLOW head contours for "predevelopment" conditions are shown as dashed lines with 10 ft intervals.The GFLOW basemap is the map of observations and head contours reported by Mendenhall (1905).GFLOW test points were associated with the 10, 20, and 30 foot contours in the West Coast Basin.GFLOW shows triangles centered on test points -red triangles with tip pointing down show model prediction less than observed, green triangles with tip pointing up show model prediction greater than observed.

Final 26Apr2023
Disclaimer: Supporting Information is not peer reviewed. In

Figure S1 .
Figure S1.Conceptualization of the Robin (Type 3) head dependent flux BC: (a) presented in the finite difference model MODFLOW for boundary cells, afterReilly (2001) and WaterlooHydrogeologic.com; and (b)  presented in the analytic element model AnAqSim as line boundary element, after AnAqSim User Guide.

Figure S2 .
Figure S2.Summary of calibration of the MODFLOW Dual Site model, Layer 9, Gage Aquifer: (a) Simulated water level contours and monitoring wells; (b) calibration statistics for all layers; (c) scatter plot of hydraulic heads all layers.AfterCH2M Hill 2008.

Figure S3 .
Figure S3.MODFLOW/MODPATH Dual Site capture zones based on tracelines for the extraction wells.The extraction wells are mostly screened in Layer 9 Gage aquifer.Extraction wells are also screen in Layer 2 Upper Bellflower aquifer (UBA) and Layers 3,4 Bellflower aquifer (BF).The injection wells are screened in the Gage aquifer.

Figure S4 .
Figure S4.GFLOW layout of the "predevelopment" Los Angeles Coastal Plain including Los Angeles County and OrangeCounty, California USA.The basemap is the USGS 1 arc second DEM.The analytic elements include constant head line-sinks (h=0.0 ft) for the Santa Monica Bay and the San Pedro Bay; specified head line-sinks (linearly variable) for the LA River, San Gabriel River, Santa Ana River, Ballona Creek, LA Narrows, Whittier Narrows; specified discharge line-sinks for mountain front recharge associated with Yorba Linda subbasin, Tustin Hills, Irvine subbasin, La Habra subbasin, Santa Monica Hills, Hollywood Hills, Elysian Hills, Repeto Hills, Merced Hills ; horizontal barrier line elements (no-flow) surrounding the basin; horizontal barriers line elements (leakage) for the Newport-Inglewood fault zones; area elements for recharge in the West Coast basin.The Dual Site Montrose-Del Amo Superfund site near Torrance is shown.For "predevelopment" conditions no pumping wells are explicitly introduced.
Figure S6, the GFLOW calibration statistics are shown for the test points.Parameterization is based on MODFLOW (Reichard et al. 2003) including effective horizontal hydraulic conductivity for West Coast Basin, and mountain front recharge discharge specified linesinks.The simple conceptual model, averaged steady flow, single layer with horizontal base elevation, single hydraulic conductivity and West Coast basin calibrated areal recharge gives minimized mean residual 0.0 ft and normalized RMSE 28.9%.The dashed lines represent the eastern artificial boundary of the USGS Paulinski et al. (2021) MODFLOW-USG.Notice that the head contours are perpendicular to this line under pre-development conditions suggesting justification of an artificial no-flow BC.Superimposing a flux condition to the OC line would reflect the influence of well extractions and barrier well injections.

Figure S6 .
Figure S6.GFLOW calibration to "predevelopment" West Coast basin test points based onMendenhall (1905) head contours  (10 ft, 20 ft, 30 ft).Parameterization based on MODFLOW(Reichard et al. 2003) including effective horizontal hydraulic conductivity for West Coast Basin, and mountain front recharge discharge specified linesinks.The statistics include the minimized mean residual 0.0 ft and normalized RMSE 28.9%.The dashed line represents the location of the Paulinski et al. (2021) Orange County MODFLOW model boundary.

Figure S7 .Figure S8 .Figure S9 .
Figure S7.Calibration of the AnAqSim model for the Dual Site scale.The bias in the Dual Site scale AnAqSim modeled heads in the Gage aquifer (Level 5) (green test points values being model too low; red being model too high) is likely due to the simplification of domain levels with horizontal bottoms.

Figure S10 .
Figure S10.Contours of averaged hydraulic heads in productive aquifers for Spring and Fall of 2006 in the West Coast Basin and Central Basin, Los Angeles County (WRD, 2007).The outer boundary of the MODFLOW grid for the Dual Site model is shown in red.The highlighted yellow area shows the influence of the fall extractions and the Newport Inglewood uplift as a hydraulic barrier.

Table S3 .
AnAqSim Dual Site model well parameterization.

Table S6 .
AnAqSim LA Basin spatially variable area sources/sinks by domain

Table S7 .
AnAqSim LA Basin spatially variable area sources/sinks by polygon.