My simulated water table is above the ground surface — what went wrong?

Almost certainly the 10×-rainfall pitfall — a recharge unit mismatch. The most common case: you entered recharge in inches per year into a field that expects meters per day (or vice versa). The platform accepts the value and simulates with it, so the model runs — but the resulting recharge flux is 10× or 100× or 365× too large, and the water table is pushed above the surface. Check your recharge value and its units in Ch. 14 §14.7. Rule of thumb: in the mid-latitudes, realistic mean annual recharge is usually 5–15 inches per year (~125–380 mm/year, ~0.34–1.04 mm/day, ~0.00034–0.00104 m/day). If your entered value is orders of magnitude outside this range, the units are probably wrong.

I changed Manning's n and nothing happened — why?

Manning's n is always stored on polylines (streams) and zones (lakes), but it is only active when the Watershed Solver is on. If you're running the groundwater solver alone, Manning's n is inert — the value sits in the model file but doesn't influence anything. To activate overland-flow physics that use Manning's n, enable the Watershed Solver (Ch. 16). This is a silent-failure case because the UI doesn't warn you that your change had no effect; the parser stores the value regardless of whether the solver will ever read it.

My Tóth-flow cross-section looks beautiful but the fluxes are unrealistic — is the model wrong?

The model is numerically correct but physically over-constrained. The prescribed-water-table approach (Ch. 10 §10.2) finds a head field consistent with Darcy's law and your drawn water table, but there's no constraint tying the implied flux to rainfall recharge. A steep water table paired with high K produces enormous simulated fluxes that no climate could supply. The geometry of flow is real; the magnitude is only realistic when the drawn water table and the K field are mutually compatible. Sanity check: compute (domain area × realistic recharge rate) as an upper-bound flux, and compare against your simulated discharge. If simulated flux exceeds this by more than ~2×, reduce K, flatten the water table, or switch to the 3D-slice approach (Ch. 10 §10.2 Approach 2) where the water table emerges as a solution consistent with rainfall.

How do I know if my submodel boundary is too close to the area of interest?

Apply the 10-cell buffer rule: submodel boundaries should be at least 10 cells away from any feature of interest inside the submodel (Ch. 13 §13.3). Closer than that, boundary effects contaminate the results you care about — heads get pinned to parent-model values, drawdowns from wells inside the submodel hit the boundary and rebound, plumes exit unrealistically. Symptom: results near the submodel edge look 'too clean' or over-constrained compared to results in the interior. Fix: either expand the submodel so the features of interest have a 10-cell buffer on all sides, or refine the parent model instead of using a submodel.

My stochastic K-field simulation looks almost identical to the uniform-K case — why isn't heterogeneity showing up?

Your grid is not resolving the correlation scale of the random field. If the correlation length Lambda is smaller than a few cell widths, the random field is effectively averaged by the discretization — the cell-to-cell K variations that theory says should produce macrodispersion are being filtered out by the grid. Rule of thumb: at least 5–10 cells must fit within one correlation length for heterogeneity to be visible in the solution. If you have a 10 km domain with 100 cells (100 m cells) and you want to see heterogeneity at a 50 m correlation scale, your cells are too big. Fix: increase NX/NY, or add a refined submodel over the area where heterogeneity matters (Ch. 13). See Ch. 19 §19.3 for full random-field setup.

What does 'structure-first, calibration-second' mean, and why does it matter?

Calibration can compensate for parameter errors but not for structural errors in the conceptual model. If your model has the wrong layer geometry, the wrong boundary-condition types, or missing features (a stream that should be there but isn't), no amount of parameter tuning will give you a physically meaningful model — you'll just be finding the parameter values that best match head observations given a broken conceptual structure. Always verify the conceptual model (geometry, boundaries, features, layering) first, and only then adjust parameters. The most common failure mode is head-only calibration that produces a 'well-calibrated' model with systematic flux errors hidden in the residuals. See Ch. 18 §18.1 and Ch. 25 (errors) for full treatment.

My head residuals are tiny but the mass balance report shows large errors — what's going on?

One of the most insidious failure modes in groundwater modeling. Head-match is the conventional calibration metric but it is blind to mass-balance errors. Several mechanisms can produce excellent heads alongside bad mass balance: an overly relaxed solver tolerance that lets iterates almost-converge cell-by-cell while accumulating budget errors globally; poorly-posed boundary conditions with two BCs in tension at the same location, producing spurious flux into or out of the domain; free-surface artifacts in transient runs where cells change saturation state; or time-stepping schemes in coupled lake-aquifer runs that drift from conservation over long periods. The danger is the head-match deceives the modeler (and client, and reviewer) into trusting the model — then mass-balance-dependent quantities like sustainable yield or plume mass are exactly what the model is most wrong about. Fix: always read the global mass balance report. If steady-state discrepancy exceeds ~1% or cumulative transient error exceeds ~2–3%, investigate before trusting derived quantities. See Ch. 22 §22.4.2 for the full pitfall treatment and Ch. 8 for reading mass balance as diagnostic conversation.

My model is stable at one K value but goes divergent or produces dry cells at a slightly different K — why is it so sensitive?

Geologic-geometry-driven ill-posedness. When the aquifer bottom (bedrock, confining unit) is locally steep, or sediment thickness varies sharply (thin sediments pinching out against rising bedrock, bedrock valley walls, glacial terrain with rapidly-varying till thickness), the saturated thickness of adjacent cells can differ by orders of magnitude. Transmissivity T = K × saturated thickness, so neighboring cells have T values differing by factors of 10 to 1000, producing a stiff, poorly-conditioned linear system. On steep bedrock slopes, small errors in computed head can push head below cell bottom and the cell dries up — whether that drying is physical or numerical depends on details that are hard to diagnose by inspection. The upshot is a narrow window of numerically stable K values around the physically correct K. Calibrating by convergence picks K for the wrong reason. Fixes: smooth spurious bedrock geometry (raster artifacts, single-borehole bull's-eyes); refine the grid in the transition zone; add sublayers near the pinchout; verify required transmissivity is geologically plausible; check whether dried cells match water-well, wetland, or spring evidence; calibrate against flux observations, not stability. See Ch. 22 §22.3.7 for the full pitfall treatment and Ch. 5, Ch. 7 §7.4, and Ch. 10 for the underlying mechanics.

Should I add more processes and more zones to make my groundwater model more accurate?

No — this is one of the most common and costly misconceptions in applied groundwater modeling. More processes and more parameters do not generally improve accuracy. Every added process brings parameters; every parameter brings uncertainty; the uncertainty compounds. A model with four zones can be well-constrained by twenty observations; the same observations cannot constrain forty zones, and calibration will find parameters that fit data without physical meaning. The right complexity depends on the modeling objective. If the goal is to understand processes, adding complexity is reasonable but should be done incrementally — start simple, add one process at a time using synthetic mode, observe the effect. If the goal is to predict the future (sustainable yield, wellhead protection, plume migration), use the simplest conceptual model that can achieve the goal. Keep it simple. If a constant-K, constant-recharge model explains the observations, don't introduce zones. When you add zones, they should be supported by independent evidence (geologic maps, pumping tests, borehole logs) — not by the desire to match local data the simpler model missed. Zones that exist only to fit observations are overfitting; they have tight fit but poor predictive skill. Each feature needs a one-line justification: if you can't write it, the feature probably shouldn't be in the model. See Ch. 22 §22.4.7 for the full treatment of the complexity fallacy, and Ch. 1 §1.3 for the refine-progressively philosophy.

What boundary condition does IGW-NET apply by default when I draw a box or polygon, and how should I think about it?

By default, IGW-NET applies a no-flow boundary condition to the domain edge — no water crosses the boundary. This is consequential: the domain water balance is determined entirely by internal sources and sinks (recharge, wells, rivers, lakes); the aquifer outside the box is treated as if it does not exist. For a domain that coincides with a true physical no-flow boundary (bedrock ridge, groundwater divide, impermeable contact), this default is correct. For a domain whose edges cut arbitrarily through flowing aquifer — almost every small-domain or site-scale model drawn on an arbitrary box — the default is physically wrong, and boundary flux can be as large as the internal sources and sinks. The recommended IGW-NET workflow is to make the boundary hydraulically remote: draw a large parent domain that puts multiple strong internal features (rivers, lakes, groundwater divides) between the boundary and the area of interest, then use a submodel (Chapter 13) to zoom into the area of interest at site-scale resolution. When the boundary is hydraulically remote, its choice becomes irrelevant to the area of interest because boundary flux is absorbed by the surface-water features and divides in between before it can propagate to the interior. The Global Base Model makes drawing large parent domains essentially free (data is pre-populated), and hierarchical submodels keep computation cheap (fine resolution only where you need it). Classical groundwater modeling agonized over hand-picking boundaries on flow lines or hydraulic features because large domains were expensive to set up and slow to solve; neither constraint applies in modern IGW-NET workflows. Verify remoteness by expanding the parent domain and confirming the area-of-interest results don't change — the sensitivity test is trivial with global data. See Ch. 22 §22.2.1 for the full treatment and Chapter 13 for hierarchical modeling.

A narrow stream in my model shows huge flux — becoming a massive source of water when it should be losing. Why?

Almost certainly a DEM-as-stage resolution failure. The DEM-as-stage mechanism assigns river stage from the DEM elevation at the model cell containing the river polyline. This works for wide rivers in gentle topography — the DEM cell center sits near the river surface. But for narrow rivers in variable topography (steep valleys, cliffs, terraces, alluvial fans), a model cell may be hundreds of meters wide while the stream is tens of meters wide, and the DEM cell center may sit on the valley wall or nearby upland rather than on the river. The model then imposes a river stage that is the upland elevation — sometimes tens of meters higher than the real river — and drives spurious recharge from this impossibly-high 'river' into the aquifer. Diagnostic: inspect per-stream mass balance (not just global), compare assigned stage values to USGS gage data or surveyed cross-sections, check the DEM value at the river cell against the actual river elevation. Fixes: refine the grid where narrow streams meet variable topography (submodel, Ch. 13); override DEM-as-stage with explicit stage values on affected reaches; use a higher-resolution DEM (LiDAR) for stream-stage lookup; or reduce narrow tributaries from Level 2 (stream polyline with DEM-as-stage) back to Level 1 (DEM-as-drain), which doesn't use stage values. See Ch. 22 §22.2.4 and Ch. 14 for the full treatment.

How should I lay out vertical layers in my groundwater model — along geologic units, or in strictly horizontal slices?

For almost all applied groundwater modeling, layers conforming to geologic units is the right pragmatic choice — this is the approach MODFLOW and most modern groundwater codes (including IGW-NET) use. Strict rectilinear discretization with uniform horizontal layers is more mathematically rigorous in the classical finite-difference sense, but it would require many horizontal layers to resolve the K structure of a real aquifer because the K contrasts concentrate at geologic boundaries that are usually oblique to horizontal. Geology-based layering resolves the same K structure with a handful of model layers, and what you give up in formal mathematical rigor you recover many times over in applied accuracy — because the error from having too few layers to represent a geologic contrast is usually far larger than the error from non-strictly-rectilinear cells. Given the horizontal-to-vertical scale disparity in most aquifers (flow within a layer is approximately horizontal), the pragmatic choice is clearly correct. However, geology-based layering assumes that geologic layers are fundamentally approximately horizontal at real scale. This assumption breaks down when geology is locally complex — when a unit pinches out to near-zero thickness, when two units converge and the layer between becomes locally singular, when a regionally 20 m-thick layer has a spot where it drops to 20 cm. Near-zero-thickness cells produce near-infinite vertical conductance and near-zero horizontal transmissivity, creating local singularities that compromise solution quality in that neighborhood. Mitigations: use IGW-NET's minimum-layer-thickness protection (any layer thinner than a configurable percentage of maximum thickness is automatically set to that minimum), add more sublayers within geologic units where complexity demands extra resolution, simplify the geologic representation locally if the complexity doesn't affect your answer, and validate results in locally-complex areas specifically. See Ch. 22 §22.3.1 for full treatment and Ch. 10 for the layering framework.

I imported SwaNET or SWAT recharge into my IGW-NET model and calibration is still off — what's wrong?

Check projection alignment first. SwaNET (and the related SWAT watershed model) require UTM internally because watershed physics (Manning's n, overland flow, flow accumulation) is only meaningful in a projected coordinate system with consistent linear distance units. IGW-NET, by contrast, is flexible about projections — it can work in geographic coordinates (lat/lon), in various UTM zones, or in other projections. This flexibility becomes a silent failure mode when SwaNET output is imported into an IGW-NET model that isn't in the matching UTM zone. The import appears to succeed — raster values go into raster cells — but the cells are geographically offset from where SwaNET computed them. Recharge highs land on wrong parts of the map; calibration can't match because the recharge field is spatially wrong regardless of its magnitudes. Fixes: always build IGW-NET in the same UTM zone SwaNET uses for the region (enable the UTM-Only checkbox in IGW-NET before drawing the domain); verify alignment visually by overlaying the two models before transferring any data; spot-check imported recharge values at specific known locations against SwaNET-reported values for the same geographic points; don't try to convert an existing non-UTM IGW-NET model after the fact — build a new UTM-projected model. When SwaNET and IGW-NET are talking to each other, SwaNET sets the projection and IGW-NET conforms. See Ch. 22 §22.2.7 and Ch. 16 for the full treatment.

Can I use IGW-NET to model karst aquifers, or aquifers with caves and major fractures?

It depends on the scale of your question. IGW-NET is built on Darcy's law, which assumes flow through a continuous porous medium where velocity is linearly proportional to hydraulic gradient. This works for sediments, sedimentary rocks, weathered bedrock, and most common aquifer materials. Darcy's law breaks down in karst (conduits and caves where flow is turbulent), major fractures (where Reynolds number exceeds the Darcy threshold), and fractured crystalline rock at small scale (where flow happens through sparse discrete fractures). However — and this is the subtlety — Darcy applicability is scale-dependent. A karst aquifer that plainly violates Darcy at the conduit scale may behave like an equivalent porous medium (EPM) when observed over kilometers to tens of kilometers, as the collective behavior of many conduits, matrix blocks, and fractures averages toward something diffusive. Whether the EPM representation in IGW-NET is useful depends entirely on the question: for regional water-balance in karst, it may be perfectly adequate; for predicting the arrival time of a contaminant in a specific well in karst, it is useless because the contaminant travels through a few conduits in hours to days rather than diffusive Darcy flow over months to years. Identifying karst indicators (large-discharge springs, rapid water-level response to storms, turbid water after rain events, sinkholes, cave records, dye-tracer velocities greater than ~10 m/day) is the first step. If the question is at a scale where individual conduits or fractures matter, use a discrete-fracture-network code or a conduit-flow model (MODFLOW-CFP, CFPv2) rather than IGW-NET. If the question is at a scale where EPM is defensible, IGW-NET is appropriate — with explicit acknowledgment in reporting that it is a coarse-resolution descriptor of averaged behavior. See Ch. 22 §22.4.8 for full treatment.

Ch. 22 — Common Pitfalls & Gotchas

How to use this chapter

These are modeling errors, not software bugs — that's what makes them hard. Every pitfall in this chapter is about how the model is built or how the software is being used, not about anything broken in the platform. The simulation runs cleanly. The solver converges. The numbers that come out are what the mathematics says they should be, given the inputs and choices you supplied. The problem is that the inputs or choices don't represent the physics you actually intended. A head field that matches observations beautifully while hiding a mass-balance error of 30%. A plume that spreads sensibly but is being blurred by numerical dispersion rather than physical dispersion. A cross-section that reproduces Tóth's classical circulation patterns but implies fluxes that no rainfall could supply. These are conceptual errors in the modeling — wrong boundary choice for the question, wrong physics applicability for the scale, wrong discretization philosophy for the geology, wrong calibration target for what needs to be predicted. There is nothing to “debug” in the code; what needs correction is the thinking that went into the model. That is why these errors are harder than outright failures — the software gives no warning because, from the software's point of view, nothing is wrong. The detection has to come from the modeler cross-checking results against independent physics (water balance, expected flux magnitudes, analytical-solution comparisons, mass-in vs. mass-out at plume scale) and against modeling judgment about what the numbers should look like for the problem being solved. This chapter is therefore not a debugging guide; it is a guide to modeling well — to recognizing where modeling decisions have gone silently wrong, and to correcting them at the level where they live.
If your simulation fails to run or produces an error message, go to Chapter 25 (Error Messages) first. This chapter is for simulations that do run but produce wrong or suspicious results.
If you know roughly what kind of problem you have, jump to the relevant section (unit issues → §22.1, boundary issues → §22.2, grid issues → §22.3, calibration issues → §22.4, workflow or integration issues → §22.5).
If you don't know where the problem lies, scan the "Symptom" line at the start of each pitfall. Pitfalls are organized so the symptom is always first — if a symptom matches what you're seeing, read the cause and fix.
Every entry cross-references the chapter that covers the mechanism in depth. This chapter is an index and diagnostic aid; the underlying material lives in Parts I–III.
The eleven most dangerous pitfalls — four foundational traps and seven “looks right, is wrong” traps:
1. Using IGW-NET where Darcy's law doesn't apply (§22.4.8) — applying the platform to karst, major-fracture, or cave systems without recognizing that Darcy-based continuum models are appropriate only at scales where the equivalent-porous-medium approximation is defensible. If the question is at the scale where individual conduits or fractures matter, IGW-NET is the wrong tool regardless of what K values you choose. The most fundamental pitfall — if applicability is wrong, nothing else matters.
2. Geology-based vertical discretization breaking down in complex geology (§22.3.1) — most groundwater platforms including IGW-NET use vertical layers that conform to geologic units, not strictly horizontal rectilinear slices. This is the pragmatic choice that made MODFLOW popular — dramatically fewer layers needed to resolve aquifer structure, efficiency gain over strict rectilinear discretization far outweighing the loss of formal mathematical rigor. But the approach relies on an assumption that geologic layers are fundamentally approximately horizontal at real scale. When geology becomes complex — locally singular, pinching out, or very thin — the assumption breaks down and local numerical trouble results. Use IGW-NET's minimum-layer-thickness protection and know where this assumption limits you.
3. The silent no-flow boundary default (§22.2.1) — IGW-NET applies a no-flow condition to the domain edge by default. For small domains cutting arbitrarily through flowing aquifer, this is physically wrong and boundary flux can dominate the interior. The IGW-NET resolution — go hydraulically remote and use submodels, rather than agonizing over flow-line picking — fundamentally changes how boundaries should be chosen in the modern workflow. Traps modelers coming from classical practice who are used to tight hand-picked boundaries and don't realize the workflow has changed with global data and hierarchical modeling.
4. The 10×-rainfall recharge unit mismatch (§22.1.1) — a model with water tables above ground, wrong by orders of magnitude
5. Prescribed-water-table flux unrealism (§22.2.2) — a 2D cross-section that reproduces Tóth's classical nested-cell flow patterns but implies fluxes that no rainfall could supply, because the numerical prescribed-water-table BC constrains head but not flux
6. Head-only calibration hiding flux error (§22.4.1) — a model that fits observed heads beautifully while missing baseflow, capture, or pumping response
7. Mass balance error despite excellent head match (§22.4.2) — related to #6 but its own distinct signature: heads look great, the global water budget doesn't close
8. Numerical dispersion (dilution) swamping physical dispersion in transport (§22.5.3) — a plume shape that looks clean and sensible but is dominated by grid-induced smearing; conservative for predicting plume extent (overstates reach) but severely non-conservative for peak concentration at a receptor (understates the real value) — the same model run gives oppositely-biased answers to “how far?” versus “how much?”
9. Geologic-geometry-driven ill-posedness (§22.3.7) — steep bedrock slopes or thin sediments pinching out, where the model is physically reasonable but only numerically stable in a narrow K window, and small K errors produce unphysical dry cells or divergence
10. DEM-as-stage resolution failure for narrow streams (§22.2.4) — a narrow stream in variable topography produces huge spurious mass-balance flux because DEM cells don't resolve the river geometry; a losing stream becomes a major water source when the DEM cell center sits on nearby uplands rather than the river surface. One of the most common causes of stream-aquifer mass-balance disasters.
11. The complexity fallacy (§22.4.7) — the assumption that simulating more processes makes a model more accurate. For applied prediction work, the opposite is usually true: parsimony beats complexity, unsupported zones become overfitting, and the simplest conceptual model that can achieve the goal is the right one. This one traps experienced modelers more than beginners.
If you're new to the platform, read at least these eleven. If you only have time for two, read §22.4.8 (wrong model for the physics) and §22.2.1 (wrong boundary for the question) — no calibration procedure can recover from either.
A map of the pitfall hierarchy — modeling errors form layers, and an error at a lower-numbered layer cannot be repaired by work at a higher-numbered layer:
1. Applicability (§22.4.8) — is IGW-NET the right tool for the physics and scale of the question? Darcy-based continuum modeling is appropriate where an equivalent-porous-medium representation is defensible; not appropriate at the scale where karst conduits, caves, or major fractures dominate.
2. Modeling philosophy (§22.4.7) — is the model's complexity matched to the question? For prediction, parsimony; for understanding, incremental process addition. Complexity added without supporting evidence becomes overfitting.
3. Discretization philosophy (§22.3.1) — does the vertical discretization represent the geology well? Geology-based layer-conforming discretization is the pragmatic choice almost universally, but breaks down in locally-complex geology where layers become singular, pinching, or very thin.
4. Boundary placement (§22.2.1) — where is the domain boundary drawn, and what does the platform do with it? A hydraulically-remote boundary with submodel refinement is the modern IGW-NET workflow; the default is silent no-flow, which is physically wrong for most small domains.
5. Input data (§22.1, §22.2.5, §22.2.7) — are the input values, units, and spatial data correct and correctly interpreted? DEM resolution matched to grid, coupled-model projections matched to IGW-NET's projection, recharge units verified.
6. Calibration and solver (§22.4) — does calibration target the right observations (heads AND fluxes), and is the solver converging for the right reasons (physics, not numerical convenience)?
7. Numerics and execution (§22.3, §22.5) — is the grid adequate for the physics, is the transport scheme appropriate, are workflow interactions understood?
An error at Layer 1 cannot be fixed by perfect execution at Layers 3–7. Diagnose from the top down: applicability first, then philosophy, then structure, then data, then numerics. Most frustrating diagnostic sessions are attempts to fix a Layer-1 or Layer-4 modeling problem by adjusting Layer-6 parameters; the calibration converges to plausible-looking nonsense because the structural error is being absorbed into K values that have no physical meaning. These are modeling errors, not bugs — they live in the thinking, and they have to be corrected there.
How to verify each pitfall in your own model — every entry includes a “How to check” recipe. This chapter is about using IGW-NET, not about abstract modeling wisdom — every pitfall includes a concrete How to check line pointing to the specific dialog, attribute panel, info display, or menu path in the IGW-NET interface where you can verify whether the pitfall is live in your own model. Is Manning's n stored on your streams but inactive because the Watershed Solver is off (§22.1.2)? Does the projection info in your DomainAttr panel match the UTM zone that SwaNET used for its recharge output (§22.2.7)? Is your submodel boundary close enough to your area of interest to violate the 10-cell buffer (§22.3.3)? Is your transport scheme the default upwinding when you need sharper peaks (§22.5.3)? Each of these can be verified in a few clicks. The combination of interface-based diagnosis and modeling understanding for interpretation is where the chapter delivers its value — the interface shows you what the model is configured to do; the modeling discussion tells you whether that is what you actually want, and why it matters for your question.

22.1 Unit and Value Traps

The commonest category. A value is numerically entered correctly but its units mismatch what the field expects, or the value is interpreted by a field whose semantics differ from intuition. The simulation runs; the results are wrong. Symptoms are usually orders-of-magnitude errors in output quantities — a water table above the surface, a recharge flux ten times too large, a K value that produces either a flat pond or a dry aquifer.

22.1.1 The 10×-rainfall pitfall (the most common pitfall in IGW-NET)

Symptom: simulated water table is above the ground surface across broad areas of the domain, or mass balance shows recharge flux 10× / 100× / 365× what the climate can plausibly supply. This is the single most frequent user error on the platform.

Cause: unit mismatch on the recharge field. The recharge input can accept different units (inches/year vs mm/year vs m/day vs m/year) depending on where and how it's entered — Data Center default, user-entered constant, zone override, imported raster. If you type 10 into a field expecting m/day but meaning inches/year, you've just entered a recharge rate roughly 1000× larger than intended. The platform accepts the value and simulates.

Fix: always confirm units on the recharge field explicitly. Rule of thumb: in mid-latitudes, realistic mean annual recharge is 5–15 inches/year (≈125–380 mm/year, ≈0.34–1.04 mm/day, ≈0.00034–0.00104 m/day). If your entered value is orders of magnitude outside this range, the units are wrong. When uncertain, use the Data Center defaults first — they carry correct units automatically. For the strongest protection against this pitfall entirely, couple IGW-NET with an external recharge or watershed model (USGS INFIL, USDA SWAT, or the SwaNET platform within MAGNET4WATER) — these produce physically-derived, properly-unit'd spatial recharge fields from climate, soil, land-use, and DEM inputs, bypassing the manual-entry step where unit mistakes happen. See Ch. 14 §14.7 for the canonical treatment, Ch. 16 for the coupled-watershed workflow, and Ch. 23 for how Data Center units are handled.

22.1.2 Manning's n stored but silently inactive

Symptom: you adjust Manning's n on a polyline (stream) or zone (lake) and the simulation result is unchanged. The UI lets you set the value, the model file stores it, but nothing happens.

Cause: Manning's n is universal — stored on every polyline and zone regardless of whether the Watershed Solver is on. But Manning's n only does anything when the Watershed Solver is active. The parser stores it whether it'll be used or not. If you're running groundwater alone (Watershed Solver off), Manning's n is inert.

Fix: enable the Watershed Solver if you want Manning's-n-driven overland flow. See Ch. 16. If you don't need overland flow, the Manning's n value simply doesn't matter — leave it at default.

How to check: open the Stream Polyline attribute dialog and look at the Manning's n value for each stream feature. It will usually show a value (the default for that stream order, or a user-entered value) — this is the stored value. Then open the Watershed Model Solver dialog to see whether the Watershed Solver is enabled. If Manning's n is stored on streams but the Watershed Solver is off, Manning's is silently inactive — the stored value is documentation only. To activate overland-flow physics that use Manning's n, enable the Watershed Solver (Ch. 16).

22.1.3 DrnL (DEM-as-drain leakance) auto-upgrade from 1 to 4 values

Symptom: an older model opened in the current platform shows different surface-drainage behavior than it used to, or the DrnL field has four values where older documentation describes one.

Cause: the DrnL field format was expanded from a single leakance value to four values (leakance, specific yield, specific storage, porosity) to support transient and unsaturated-zone mass-balance protection. The parser auto-upgrades old 1-value files to the 4-value format, using sensible defaults for the new fields. The behavior is equivalent for steady-state flow but can differ for transient runs.

Fix: no action needed for most cases — the auto-upgrade produces sensible defaults. For transient work where the DrnL storage terms matter, verify the upgraded values in the parser output. See Ch. 5 §5.8 for DrnL semantics and Ch. 14 §14.2.1 for the surface-drainage context.

How to check: open the DrnL (surface-drainage leakance) dialog. Newer IGW-NET versions show four values per entry: leakance, specific yield, specific storage, and porosity. If you opened an older model file that had only a single leakance value, the platform auto-upgrades to the 4-value format using sensible defaults for the three new fields. For steady-state runs this is invisible. For transient runs, open the DrnL dialog and replace those defaults with values that match your actual aquifer material — accepting the auto-upgrade defaults silently can bias transient results.

22.1.4 K source interpretation ambiguity

Symptom: you set a K value in the attributes dialog but the simulation uses a different value; or a random-field K looks the same as uniform K even though you enabled randomization.

Cause: the K field supports multiple sources (constant, Data Center raster, scattered points from wells, T-PROGS geostatistical realization, random field) and the interpretation depends on which source is selected. Entering a new constant value doesn't override a K-source that's currently set to "Data Center" — the constant is stored but the raster takes precedence. Enabling random-field generation on a zone that's inheriting global K doesn't automatically activate randomization on the zone; you have to flip the K source for the zone to "random field."

Fix: check which K source is active for the zone or domain you're configuring. The source selector is in the attributes dialog, and the parser field .hydraulic_conductivity encodes both the source type and the values. See Ch. 5 for the full customization ladder, Ch. 17 for T-PROGS sources, Ch. 19 for random-field sources.

How to check: open the K Attribute dialog for the zone or domain you're configuring. The source selector at the top (“Constant,” “Data Center,” “Random Field,” “T-PROGS realization,” “Import”) determines which K value is actually used — the other fields may be populated but inert. Entering a new constant value does not override a K source currently set to “Data Center”; the constant is stored, but the raster takes precedence. To enable randomization on a zone that inherits global K, you must switch the zone's K source to “Random Field.” The Info panel for the current model gives an at-a-glance summary of which source is active for K, BotE, and recharge across the whole domain.

22.2 Boundary-Condition Traps

Boundary conditions are where the conceptual model meets the numerics. Most BC-related pitfalls come from implicit assumptions about what a given BC type actually constrains — head, flux, or both — and what the platform does when those assumptions don't match the real physics of the site. This is the largest section in the chapter because BCs are where the most consequential modeling decisions live. Three entries carry the weight and reward careful reading: §22.2.1 — the silent no-flow boundary default, which is the foundational BC choice that precedes all the others and which the modern IGW-NET workflow (hydraulically-remote boundary + submodel) fundamentally changes; §22.2.5 — DEM resolution choice and artifacts, which bundles together the “higher is always better” misconception, the LiDAR-through-DataNET architecture that resolves the data-volume tradeoff, the modeler-vigilance discipline for DEM artifacts, and the flagship application to headwater GW-dependent ecosystems; and §22.2.7 — the SwaNET/SWAT recharge handshake, which catches the silent projection-mismatch failure that can make a spatially-variable recharge field land in the wrong place on the map. The remaining entries cover specific traps that occur once the major BC decisions have been made correctly.

22.2.1 The silent no-flow boundary default — and the hydraulically-remote strategy that resolves it

Symptom: a small-domain model (a site, a parcel, a short stream reach) produces heads and flow directions near the boundary that don't match observations or common sense. Pumping tests compress against an invisible wall. Regional head gradients that should pass through the model area are truncated at the domain edge. A stream that should carry flow out of the domain shows no baseflow leaving. In larger domains the symptom can be milder — just a broad area of subtly wrong results near the boundary — but the underlying cause is the same: the model boundary is doing something physical the modeler didn't intend.

Cause: when a user draws a box or polygon to define the model domain, IGW-NET assigns a no-flow boundary condition to the domain edge by default. No flow crosses the boundary; no water enters or leaves through it. This is a consequential assumption: no-flow means the domain's water balance is determined entirely by internal sources and sinks (recharge, wells, rivers, lakes) — the aquifer beyond the box is treated as if it doesn't exist. For a domain that happens to coincide with a physical no-flow boundary (a true groundwater divide, a bedrock ridge, an impermeable contact), this default is correct. For a domain whose edges cut arbitrarily through flowing aquifer — which is the case for almost every small-domain or site-scale model drawn on an arbitrary box — the default is physically wrong, and boundary flux can be as large as or larger than the internal sources and sinks the modeler is actually trying to simulate. The silent nature of this is the trap: the user is rarely told that a default has been applied; they simply see a simulation that runs, and the boundary condition that dominates their results is a choice they never explicitly made.

The historical agony of boundary selection. Traditional groundwater modeling has dealt with this by hand-picking the boundary to coincide with hydrologically meaningful features: a flow line (where no-flow is correct by symmetry), a surface-water body (where a head boundary can be justified), a groundwater divide (where no-flow is correct). This is hard. Flow lines shift with pumping and seasonal recharge. Surface-water levels change. Divides are hypotheses, not observations. Inferring a defensible boundary from sparse head data is an exercise in educated guessing, and the chosen boundary often becomes the dominant source of uncertainty in the model. Entire chapters of classical groundwater texts are devoted to this problem, and it remains one of the most time-consuming steps in setting up a professional site-scale model the traditional way.

The IGW-NET resolution — go hydraulically remote, use submodels to focus. IGW-NET's hierarchical modeling framework, combined with the Global Base Model, turns this problem on its head. Rather than agonizing over boundary placement, the recommended workflow is:

Draw a large box around your area of interest — large enough that the boundaries are hydraulically remote, meaning far enough away that whatever happens at the boundary cannot propagate to the area of interest over the simulation timescale.
Accept the default no-flow boundary — because when the boundary is hydraulically remote and the large domain contains multiple strong internal sources and sinks (rivers, lakes, groundwater divides, pumping centers), the no-flow choice becomes irrelevant to the area of interest. Any spurious flux at the boundary gets absorbed by nearby surface-water features or divides long before it could affect the interior.
Use a submodel to zoom into the area of interest (Chapter 13) with the resolution and detail that actual analysis requires. The parent model supplies physically-realistic BCs to the submodel automatically, derived from the solved head field — no hand-picking, no flow-line hunting.

This works because of two features that weren't available in classical modeling. First, global base data — the pre-assembled DEM, stream network, aquifer properties, and recharge fields that cover Earth — make drawing a regional-scale domain essentially free. You aren't choosing between a tight domain that's fast to build and a large domain that takes weeks of data collection; both are equally easy. Second, hierarchical submodels mean you don't pay a computational cost for the large parent domain, because high-resolution work happens in the submodel. The parent solves at regional scale and coarse resolution; the submodel inherits BCs and resolves the area of interest at site scale. Both run fast.

The rule of thumb for “hydraulically remote”: the boundary should be far enough from the area of interest that several strong features (rivers, lakes, divides) lie between them. A good working diagnostic is that the boundary should enclose the full drainage context of the area of interest — typically a whole subwatershed, not a parcel-scale box. Modern global data makes this trivial to set up.

Fix:

Do not draw small boxes on arbitrary boundaries. A small-domain model with no-flow edges cutting through flowing aquifer is almost always wrong in a hidden way.
Go hydraulically remote by default. Draw the parent domain large enough that multiple strong internal features (rivers, lakes, major divides) lie between the boundary and the area of interest. Use the Global Base Model to populate this large domain without effort.
Use submodels for focus (Chapter 13). The 10-cell buffer rule in §22.3.3 still applies at the submodel scale, but the parent-model boundary doesn't need the same careful treatment because by construction it is remote from the area of interest.
Test remoteness by sensitivity analysis. Because the Global Base Model makes domain expansion trivial, it is easy to rerun with a larger parent domain and confirm that results in the area of interest don't change. If they do change, the original boundary wasn't remote enough — expand further. If they don't, the boundary is adequately remote and you can trust the submodel results. This sensitivity check is one of the clearest practical validations of model setup available in modern groundwater modeling.
Only hand-pick boundaries when there's strong physical reason to. If a known-impermeable bedrock contact or a reliable surface-water body passes close to the area of interest and you're confident it represents a true BC, using it can reduce the domain size — but in practice, the remote + submodel approach is usually simpler, more defensible, and more transparent.

The underlying principle: hydraulically-remote boundaries are not a compromise — they are the correct BC strategy when global data and hierarchical modeling are available. Classical modeling demanded tight, hand-picked boundaries because collecting data for a large domain was prohibitively expensive and solving a large model was slow. Neither constraint applies anymore. Modelers who still agonize over boundary-line hunting when a remote-box-with-submodel approach is available are solving a problem that the tool has already solved for them. See Chapter 13 for the hierarchical framework, Chapter 5 for the Global Base Model details, and Chapter 4 for domain-drawing mechanics.

How to check: look at the Domain outline on the map — is it a small box cutting across flowing aquifer, or a large parent domain enclosing the full drainage context? If it is a small box and the model is not marked as a submodel (no parent in the model hierarchy), the default no-flow BC is silently in effect along all edges. Two fixes: enlarge the parent domain until multiple strong features (rivers, lakes, divides) sit between the boundary and the area of interest, or convert to a submodel (Ch. 13) that inherits BCs from a larger parent model. The model Info panel tells you whether you're looking at a parent-level model or a submodel.

22.2.2 Prescribed-water-table flux unrealism — a trap in numerical reimplementations of Tóth-style boundary conditions

Symptom: a 2D vertical-profile synthetic model produces a beautiful Tóth regional-flow pattern, but the total flux moving through the system is one or two orders of magnitude larger than rainfall could supply. Head distribution looks physically sensible; the flux doesn't.

Cause: József Tóth's classical 1962–63 analysis of regional groundwater flow derived its elegant nested-cell circulation pattern from exactly this boundary condition — a prescribed water table serving as the upper boundary, with head specified along a smooth curve representing topography. Tóth's original solution is analytical, presented as a long, complex infinite-series expression, and it correctly predicts the head field that satisfies Darcy's law under that BC. The solution is a landmark in hydrogeology. The trap is not in Tóth's formulation — it is in the numerical reimplementation of a Tóth-style prescribed-water-table BC, where modelers (especially those learning the method by reproducing Tóth's classical patterns in 2D vertical profiles) may not realize that the prescribed-water-table approach (Ch. 10 §10.2, Approach 1) constrains head along the water-table line (head = y) but does not constrain flux. The numerical solver finds whatever flow field satisfies Darcy's law given the drawn water table and the assigned K. If K is too large for the water-table slope, the implied flux is enormous — and there is no constraint tying it to rainfall or any other physical source. Tóth's original analytical treatment carried the same mathematical feature, but in his usage the water-table shape was chosen to represent a plausible real landscape and the K values were treated as hypothetical — the infinite-series solution was interpreted for what it demonstrated about flow-system structure, not for what it implied about fluxes. A modern numerical model can accidentally read fluxes from the solution as if they were physical, and get a very wrong answer.

Fix: perform a flux-realism sanity check. Compute (domain horizontal area × plausible recharge rate) as an upper-bound total flux the climate could supply. If simulated flux exceeds this by more than ~2×, either reduce K, flatten the drawn water table, or switch to the 3D-slice approach (Ch. 10 §10.2, Approach 2), where recharge is the input and the water table emerges as the solution — guaranteeing water-balance consistency by construction.

How to check: in the Layering dialog, confirm which water-table formulation is active — the prescribed-water-table approach (Approach 1 in Ch. 10 §10.2) constrains head along the water-table line but not flux. If you are using this approach, do a flux-realism sanity check from the Flow Summary after simulation: compare total simulated through-flow against (domain area × plausible recharge rate). If simulated flux exceeds this upper bound by more than ~2×, reduce K in the K Attribute dialog, flatten the drawn water table, or switch to the 3D-slice approach (Approach 2) where recharge is the input and the water table emerges as a solution consistent with rainfall.

22.2.3 Stream-order default mismatch with desired SW level

Symptom: streams you drew appear in the simulation with different behavior than you expected — some are acting like drains, others like prescribed-head, and the boundary response isn't what your conceptual model called for.

Cause: the platform auto-populates stream representation based on Strahler stream order — small streams default to Level 1 (DEM-as-drain), medium streams to Level 2 (explicit stream polylines with one-way or two-way coupling), large streams and lakes to Level 3 (coupled lake-aquifer with ODE water balance). The defaults are usually correct, but if your conceptual model has strong regional flow driven by a mid-order stream that happened to auto-populate as Level 1, you'll see flux behavior that doesn't match reality.

Fix: review the auto-populated stream representation on the UI and override by feature where your conceptual model requires a different level. See Ch. 14 for the three-level framework and override mechanics. For coupled lake-aquifer work specifically, Ch. 15. See also §22.2.4 for the DEM-as-stage narrow-stream resolution trap, and §22.2.5 for the broader DEM resolution-and-artifacts pitfall.

How to check: the Stream Polyline dialog shows each stream feature's order, width, Manning's n, flow direction, and stage. Auto-population assigns stream order from hydrography data — review the auto-populated orders against your conceptual model, and override per feature where the desired SW representation level differs. Also verify from the Watershed Model Solver dialog or the Info panel that streams are active globally; if they're disabled at the global level, none of the per-polyline settings matter.

22.2.4 DEM-as-stage failure for narrow streams in variable topography

Symptom: a stream that should be gaining or losing modestly appears in the simulation as an enormous source or sink of water — mass balance shows huge stream-aquifer exchange flux, local heads are pinned to values tens of meters different from the broader water table, and what should have been a subtle hydrologic feature dominates the model's water budget. Often the affected reach is a narrow creek, a stream in a steep valley, or a river crossing variable topography such as cliffs, terraces, or alluvial fans. A losing stream appears to be recharging the aquifer at rainfall-exceeding rates, or a gaining stream appears to be drying the aquifer far faster than it could sustain.

Cause: the DEM-as-stage mechanism assigns river stage from the DEM elevation at the model cell containing the river polyline. This works beautifully for wide rivers in gentle topography — the DEM cell center is close to the river surface and the assigned stage is physically sensible. But for narrow rivers in variable topography, the assumption breaks down: a narrow river may be tens of meters wide, but a model cell may be hundreds of meters wide, and the DEM cell center (which provides the stage elevation) may sit on the valley wall, the cliff top, a terrace, or adjacent high ground rather than on the actual river surface. The model then imposes a river stage that is the nearby upland elevation, sometimes tens of meters higher than the real river. The prescribed-head mechanism (or the head-dependent leakage mechanism) then drives flux from this impossibly-high "river" into the aquifer, producing the giant spurious source term.

This is one of the most common causes of mass-balance disasters in site-scale or regional models that include narrow tributaries in variable terrain. It is particularly insidious because the river polyline is drawn correctly, the DEM is correct, the stream order is correct, and the feature is being represented at the right SW level — the error lives entirely in the resolution mismatch between model cells and the actual river geometry. Inspecting individual river-cell stage values against independent knowledge (USGS gage data, surveyed cross-sections, photos of the channel) will reveal the problem immediately, but few modelers think to check.

Fix: several mitigations, usable individually or in combination:

Refine the grid where narrow streams meet variable topography. A submodel (Ch. 13) with cell size comparable to the stream width resolves the river-cell geometry so the DEM value at that cell center is actually the river surface, not the valley wall.
Override DEM-as-stage with explicit stage values on affected reaches. For reaches where the DEM-cell-center doesn't correspond to the river surface, provide the stage directly — from USGS gage data, from stream-gage cross-sections, or from a manually-specified elevation that represents the channel. Ch. 14 covers the override mechanics in the stream-polyline attributes.
Use a higher-resolution DEM where available. A 1 m LiDAR DEM resolves narrow channels that a 10 m or 30 m regional DEM cannot. If the model is regional and a high-resolution DEM isn't feasible for the whole area, use the high-resolution DEM just for the stream-stage lookup (a workflow Ch. 23 supports).
Reduce the stream representation level for narrow tributaries in steep terrain. If the small tributary isn't driving the regional flow pattern, dropping it from Level 2 (explicit stream polyline with DEM-as-stage) back to Level 1 (DEM-as-drain) eliminates the problem entirely — Level 1 doesn't use stage values, so resolution mismatch can't produce spurious flux. Ch. 14 §14.2 covers when this trade-off is appropriate.
Always inspect per-feature mass balance. The global mass balance may look fine if a losing stream's spurious source is counterbalanced by legitimate discharge elsewhere. Per-stream mass balance is where this pitfall surfaces cleanly — a narrow stream contributing megaliters/day when it should contribute liters/day is immediately diagnostic.

The underlying principle: DEM-as-stage is an excellent choice when the DEM resolution, model-cell size, and river width are all commensurate. For wide rivers in gentle topography, they usually are. For narrow rivers in steep or variable topography, they often are not — and the user has to verify rather than assume. See Ch. 14 for the DEM-as-stage canonical treatment and Ch. 5 for the DEM customization ladder.

How to check: compare three numbers visible in the IGW-NET interface. (1) The Stream Polyline dialog shows width for each stream feature. (2) The Data Center DEM configuration shows the DEM resolution in use (90 m / 30 m / 10 m / 1 m LiDAR). (3) The model Info panel shows model cell size. If stream width is smaller than either DEM resolution or cell size, this pitfall is live — the DEM value assigned at the river cell is almost certainly not the river surface. Mitigate by refining via submodel (Ch. 13), by entering explicit stage values on affected reaches in the Stream Polyline dialog, or by reducing the affected polylines from Level 2 (stream polyline with DEM-as-stage) to Level 1 (DEM-as-drain, which doesn't use stage values).

22.2.5 DEM resolution choice and artifacts — “higher resolution is always better” is not quite true

Symptom: one of several patterns, all stemming from DEM-choice misconceptions or uninspected DEM artifacts. (a) A beginner assumes that because DEM is such a critical dataset, the finest-resolution available must always be the right choice, and tries to specify LiDAR 1 m for every domain regardless of scale — not realizing that for most regional models a coarser DEM captures exactly the same information because the grid cell is much larger than any LiDAR pixel. (b) A simulated water table converges to a local low that doesn't correspond to any real topographic feature — investigation reveals a hole or pit in the DEM, typically a raster processing artifact or an unmapped gravel mine. (c) The model produces physically unreasonable surface-drainage patterns along a specific reach — the DEM has a systematic bias there (cloud shadow in the source imagery, snowpack during acquisition, vegetation-canopy top rather than ground surface). (d) A site-scale model of small surface-water features (headwater streams, wetlands, small channels) fails to resolve them adequately because a 30 m DEM was used where LiDAR would have been far more informative — the opposite mistake from (a).

Cause: three distinct issues to untangle.

First, the resolution-choice principle. DEM is arguably the single most important dataset that transformed hydrologic and groundwater modeling — it determines the top of the aquifer, the drainage network, river stages (when used as-stage), and the surface-drainage relief. Multiple DEM resolutions are available across North America: 90 m (SRTM, global coverage), 30 m (NED, continental US, ASTER globally), 10 m (NED in many US regions), and 1 m or finer (LiDAR, now flown across most of North America and served from USGS). A natural beginner assumption is that higher resolution is always better. But DEM resolution should match model grid resolution, not just availability. A regional model with 500 m cells gains nothing by using a 1 m LiDAR DEM; the DEM value assigned to each 500 m cell is an average or sampled value regardless of input resolution, and the extra detail is discarded in the assignment. A 30 m or even 90 m DEM is entirely appropriate for a regional model.

Second, IGW-NET's architecture resolves the data-volume concern that would otherwise limit LiDAR use. This is a subtle point worth understanding clearly, because a practitioner familiar with other modeling platforms may carry over the wrong mental model. LiDAR DEM data for a large region is genuinely massive — tens of gigabytes to terabytes for a state-scale area, and effectively impossible to move as a local file. IGW-NET's interface shows 1 m LiDAR as a selectable DEM option, but the actual 1 m DEM is not stored in the platform or moved to the user's session. Instead, LiDAR is delivered dynamically through the DataNET services hub, which accesses USGS WMS/WCS services providing LiDAR coverage for all of North America. When a modeler selects the 1 m LiDAR option, DataNET automatically subsamples the LiDAR data to match the model's domain size: for a large area, the subsampled resolution is coarser (appropriate to the model's grid); for a small area, the original 1 m resolution is preserved. The user's model file stays light — effectively immobile data, served to the model on demand. This architecture means the tradeoff between LiDAR's detail and LiDAR's data volume has been resolved at the platform level, and the practical barrier to using LiDAR for local-scale modeling — the sheer size of the dataset — simply does not apply in IGW-NET. This is one of the quiet but consequential capabilities of the platform: hierarchical modeling coupled with LiDAR through DataNET enables modeling across scales with detail where it matters, massive underlying datasets, and a system controlled by physics rather than by data-logistics limitations.

Third, the DEM-artifacts issue. Even the best DEM is not perfect. Any DEM will contain locally-spurious values — holes or pits where drainage processing introduced artifacts, systematic bias from acquisition conditions (cloud cover, snow, vegetation canopy mistaken for ground surface), outdated data where landscape has changed (mining, construction, land subsidence, river channel migration), or straightforwardly incorrect values in rare cells. These are usually small in spatial extent but can have disproportionate effect where they occur: a hole in the DEM can act as a spurious drain in the model, causing the simulated water table to converge to it and producing a local low that doesn't correspond to any real feature. A systematic elevation bias along a river reach can make a losing stream appear gaining or vice versa. A gravel-mining pit or quarry that exists in reality but not yet in the DEM (or vice versa) creates genuine confusion — is the model's behavior wrong because of a DEM error, or is it right and the DEM simply doesn't yet reflect recent human activity? Either case requires the modeler to recognize the anomaly, not to accept it silently. Groundwater modeling has not reached the point where it can be used as a calculator with blind trust in inputs — a vigilant user is still the first-order quality control.

The scope-of-impact insight makes artifact diagnosis tractable: the DEM directly affects the simulation only where the water table reaches the surface. In most of the domain, the water table sits below ground by several meters to hundreds of meters, and the DEM value at that cell enters the solution only through (a) the top-of-aquifer geometry for layer thickness calculations and (b) DEM-as-drain leakance if enabled. A small local DEM error in such an area typically has negligible effect on the solution. The places where DEM errors do matter are specifically (a) where rivers, lakes, and wetlands exist — because those are the places where the water table is at or near the surface, and the DEM value feeds directly into surface-water boundary conditions (§22.2.4), (b) where DEM-as-drain is active and surface drainage is a meaningful fraction of the water balance, and (c) where the water table emerges naturally in the solution — discharge zones, seeps, springs. A DEM error in any of these locations can drive the simulation wrong; a DEM error in a dry upland far from these features usually doesn't matter.

Fix:

Match DEM resolution to model grid resolution, and trust the platform's defaults. For regional models with coarse cells (hundreds of meters to kilometers), the default DEM selection from the Data Center (typically 30 m or 90 m) is appropriate. For site-scale models with tens-of-meters cells, 10 m DEM is usually the right default. The platform generally selects reasonable resolution by grid cell size — accept the default unless you have specific reason to override.
Use LiDAR through DataNET for local-scale and small-feature work. LiDAR is genuinely valuable — often transformative — for site-scale modeling that must resolve small surface-water features: headwater streams, small wetlands, narrow channels in variable topography, urban drainage features, and the stream-stage detail that resolves the §22.2.4 narrow-stream pitfall. To use LiDAR, select the 1 m LiDAR option in the DEM configuration; DataNET handles the rest — streaming data from USGS WMS/WCS services, automatically subsampling for the domain size, and serving the result to the model without putting a massive file in your session. The few extra clicks to opt into LiDAR are well worth it for local-scale work where small features matter. See Ch. 23 for DataNET DEM options.
Combine LiDAR with hierarchical modeling. The most powerful pattern is a regional parent model with standard-resolution DEM and a fine-resolution submodel using LiDAR over the area where small features matter. Because IGW-NET resolves the data-volume concern automatically — subsampling for the regional domain, keeping original LiDAR resolution for the submodel — this combination is practical in a way it is not in platforms that require manual DEM preparation. Hierarchical modeling + LiDAR + DataNET is a distinctive IGW-NET capability worth exploiting when the problem calls for it. See Ch. 13 for the nested modeling framework.
Inspect the model output for DEM-driven artifacts. After the first simulation, examine the water-table map and the surface-drainage pattern. Does the water table converge to any unexpected local minimum? Is there a sharp drainage pattern that doesn't correspond to a real feature? Either can signal a DEM issue. If you see a local low that has no corresponding surface-water feature, cross-check the DEM elevation at that location against aerial imagery or a topographic map.
Understand before you patch. A local low in the model can be: (a) a real DEM artifact (processing hole, cloud shadow, canopy bias), (b) a real landscape feature not yet in the DEM (gravel mine, new quarry, subsidence), (c) a real hydrologic feature (closed depression, sinkhole that actually drains the water table) that the model is correctly capturing, or (d) a computational artifact in the simulation itself. The response depends on which — patch a DEM hole, update the landscape, accept the real feature, or fix the simulation. Do not silently patch every local low; that can hide real information.
Patch DEM artifacts locally when appropriate — the impact is usually contained. Because DEM only directly affects the solution where the water table is at the surface, a local DEM patch in an upland area of the domain (where the water table sits meters or tens of meters below the surface) has minimal effect on the overall solution. Where DEM patching is needed near a river, lake, or wetland, do it carefully: the patch propagates into the SW boundary condition and can be noticed in head or flux results. See Ch. 5 for DEM customization options.
Accept that current groundwater modeling requires modeler vigilance. The field has not reached the point where a model can be built, run, and trusted as a calculator. A model with a subtly-wrong DEM can produce a subtly-wrong answer that calibration will compensate for (pushing K or recharge values to unphysical ranges to fit heads around the DEM error) without ever flagging the underlying problem. First-principles inspection of the results — does the water-table map make physical sense? do the drainage patterns match what you'd expect? — is the quality-control step that catches these issues. See §22.4.3 on structure-first and §22.4.1 on head-only-calibration hazards.

The underlying principle: DEM is the most important single dataset in groundwater modeling, precisely because it connects the model to real landscape. Using it well means choosing a resolution that matches the question (the platform's defaults usually do this correctly), exploiting the LiDAR-through-DataNET capability for local-scale work where small features matter, inspecting the results for artifacts, understanding the cause before patching, and retaining vigilance — the work is not done when the simulation finishes. The combination of hierarchical modeling and streamed LiDAR through DataNET is distinctive to IGW-NET and expands what a single user can accomplish in groundwater modeling: massive underlying datasets, modeling across scales, detail where it matters, system control by physics rather than by data logistics — with a model file that stays light because the data remains virtually immobile on the server. See Ch. 5 for DEM customization and auto-by-grid selection, Ch. 23 for DataNET DEM options including LiDAR delivery, Ch. 13 for hierarchical modeling, and Ch. 14 for how DEM feeds surface-water boundary conditions.

Where LiDAR matters most — headwater groundwater-dependent ecosystems

The most consequential application of the LiDAR-through-DataNET capability is modeling of groundwater flow systems in headwater areas — the upper reaches of watersheds where small streams begin, where groundwater discharges as seeps and springs, and where groundwater-dependent ecosystems (fens, hanging bogs, coldwater headwater stream reaches, wetland complexes) concentrate. These are fragile landscapes, often where a small piece of land supports disproportionate biodiversity. They are also landscapes under active competing pressure: bottled-water extraction industries are economically drawn to the high-quality groundwater available in headwater areas, while ecological-conservation communities — land trusts, government agencies, academic groups — spend substantial resources documenting GW-dependent features and the ecosystems that depend on them in order to protect what the extraction pressure threatens.

Groundwater modeling in such areas has historically been limited by DEM resolution. Even a 10 m DEM is often insufficient to resolve the topographic features that control groundwater discharge in headwater terrain: small valley floors a few cells wide, subtle breaks in slope that define seep locations, narrow wetland margins where vegetation and hydrology transition sharply, headwater stream channels only a few meters wide. A model that cannot see these features cannot accurately represent where groundwater emerges, which wetlands are GW-dependent versus surface-fed, or how pumping at one scale affects discharge at another. Modelers working on these questions with classical workflows faced an unwelcome choice: model coarsely and miss the physics that matters, or prepare LiDAR manually through a days-to-weeks data-management exercise for each site.

1 m LiDAR combined with IGW-NET's hierarchical modeling and DataNET streaming changes the math. A regional parent at standard DEM captures the watershed-scale flow system; a local submodel at LiDAR over the headwater feature resolves the topography that controls discharge. The LiDAR data never moves as a file — it streams from USGS WMS/WCS through DataNET, subsampled appropriately for the regional parent and preserved at original resolution for the local submodel. The model file stays light. What used to require weeks of data-preparation infrastructure now takes a few extra clicks. For the modelers, land trusts, agencies, and researchers working to understand and protect headwater groundwater-dependent ecosystems, this is a genuine shift in what is practically feasible — and for communities balancing extraction pressure against ecosystem protection, it is the kind of modeling capability that can support defensible, quantitative decisions.

How to check: open the Data Center DEM dialog and note the selected DEM resolution. Compare it against the model cell size shown in the Info panel. If DEM resolution is much finer than cell size (for example 1 m LiDAR in a model with 500 m cells), you are not benefiting from the extra detail — either let the platform auto-select DEM by grid, or drop to a coarser DEM to save query volume. If the model is small-scale and needs to resolve narrow streams or small wetlands, switch to the 1 m LiDAR option; DataNET's auto-subsampling handles the data volume transparently. If a custom DEM has been uploaded to override the Data Center default, the DEM dialog will show the imported file — review that source against expected coverage and quality.

22.2.6 Pumping near a surface-water feature — two cases, two treatments

Symptom: a pumping well placed near a stream or lake produces either (a) an implausibly flat drawdown that barely extends beyond the well, or (b) implausibly rapid lake depletion.

Cause: a pumping well drawing from an aquifer near surface water has two distinct physical cases depending on whether the SW body is an effectively infinite source or a small lake with a limited water budget. If the SW is represented as prescribed-head (Level 2 with high-K streambed), the well can draw unlimited water from the SW — drawdown stays shallow. If the SW is a lake with two-way ODE coupling (Level 3), the well can deplete the lake faster than the watershed can refill it, producing unrealistic rapid decline.

Fix: match the SW representation level to the physical reality. For a well pumping from a large river or lake with strong recharge from its watershed, prescribed-head behavior is realistic. For a well pumping from a small isolated lake, the lake's water budget must be represented explicitly. See Ch. 14 and Ch. 15.

How to check: open the Well properties for the well near the surface-water feature and note its pumping rate. Then open the properties of the nearby SW feature: if it is a coupled lake (Level 3), the Lake dialog shows the lake's water-balance components (precipitation, runoff inflow, GW exchange, direct outflow, pumping). Confirm the Ins and Outs close to a plausible balance at the well's pumping rate — if pumping rate exceeds plausible Ins, lake depletion will be unrealistically fast. If the feature is a prescribed-head stream or lake (Level 2), pumping treats it as an infinite source; check whether that is physically reasonable for the water body you're modeling.

22.2.7 Projection mismatch — the SwaNET/SWAT recharge handshake

Symptom: one of several failure modes, depending on which side of the handshake breaks down. (a) IGW-NET and SwaNET models of the same region don't align spatially when overlaid. (b) SwaNET-computed recharge is successfully transferred to IGW-NET as an imported raster, but the recharge appears in the wrong location — offset by hundreds of meters to several kilometers, or geographically sensible but with wrong units because of the projection misinterpretation. (c) The simulation runs and the global water budget looks roughly plausible, but per-cell recharge values don't match their SwaNET source — areas that SwaNET modeled as high recharge (e.g., forested uplands) show up in IGW-NET as low recharge, and vice versa. (d) In the worst case, the apparent symptom is indistinguishable from ordinary calibration difficulty: you can't get baseflow or heads to match observations, and you spend weeks adjusting K and recharge multipliers when the real problem is that your recharge field is in the wrong place on the map.

Cause: a fundamental asymmetry between the two platforms. IGW-NET is flexible about projections — it can work in geographic coordinates (lat/lon), in various UTM zones, or in other regional projected coordinate systems, because groundwater modeling at the site-to-regional scale can be meaningful in any of these. SwaNET, by contrast, uses UTM internally because watershed modeling depends on derivatives of the DEM (slopes, flow accumulation, flow direction) that are only physically meaningful in a projected coordinate system with consistent linear units. This means when the two models are working on the same landscape, SwaNET's native projection is UTM and IGW-NET must adapt to match. If IGW-NET is in a different projection (or in geographic coordinates) and SwaNET output is imported, the import may succeed numerically — raster values go into raster cells — but the cells themselves are misaligned with the IGW-NET model grid. The platform often does not warn about this because each individual step is technically valid: SwaNET produces a valid raster, IGW-NET accepts a valid raster, and the projection mismatch only manifests when the two are combined into a coherent physical picture. This is why the failure is so commonly silent.

The problem is particularly consequential for SwaNET-computed recharge because recharge is one of the most sensitive inputs in a groundwater model — a spatially-variable recharge field from SwaNET (or from a companion SWAT run) can be the single most valuable piece of input for calibration, but only if it lands in the right geographic location on the IGW-NET grid. A 1 km offset between recharge high and low zones can mean the difference between a model that matches baseflow observations and one that systematically misses them regardless of how you tune K.

Fix: always build IGW-NET in the same UTM zone that SwaNET is using for the same region. Specific steps:

Enable UTM-only at model setup. In IGW-NET, enable the UTM-Only checkbox before drawing the domain. This forces IGW-NET to use UTM coordinates from the start, preventing later projection conflicts. See Ch. 16 §16.6.3 for the UTM-Only checkbox location and behavior.
Match UTM zones explicitly. North American projects split across UTM zones require a specific choice: either use a single zone and accept distortion near the edge, or use the zone appropriate for your study area's center. If SwaNET has already been run for the watershed, use the same UTM zone IGW-NET that SwaNET used. If you're setting up both from scratch, pick the zone before drawing either.
Verify alignment by overlay before transferring data. After setting up both models, display them on the same map and confirm the SwaNET watershed boundary, stream network, and land-use polygons align with the IGW-NET domain, features, and observation points. If they don't align visually, transferring recharge will fail silently in the data domain — fix the projection before importing anything.
Spot-check recharge values after import. Pick a few known locations (a specific well, a stream gage, a land-use transition) and confirm the imported SwaNET-recharge value at that IGW-NET cell matches the SwaNET-reported value for that same geographic location. If the values disagree, alignment is off even if the visual overlay looked right — often a datum-level (NAD83 vs WGS84) or units-level issue on top of the projection mismatch.
When retrofitting an existing IGW-NET model to accept SwaNET recharge, don't try to convert the existing non-UTM model to UTM after the fact. Build a new UTM-projected IGW-NET model aligned to SwaNET and transfer features explicitly. Post-hoc projection conversion in a complex model is error-prone and usually worse than starting fresh in the right projection.

The underlying principle: SwaNET's UTM requirement is a physical one, not a software choice — watershed physics (Manning's n, slope-driven overland flow, flow accumulation) requires a consistent linear-distance coordinate system, and UTM is the standard choice. IGW-NET's flexibility about projections is a feature in standalone work but becomes a trap in cross-platform workflows. When SwaNET and IGW-NET are talking to each other, SwaNET sets the projection and IGW-NET conforms. See Ch. 16 for the full SwaNET-handshake treatment and the two-platform principle.

How to check: open DomainAttr → Miscellaneous and confirm the UTM Only checkbox is enabled, and that the displayed UTM zone and datum match whatever SwaNET is using for the same region (for example, Zone 16N NAD83 for southern Michigan, Zone 12N for Arizona, and so on). If UTM Only is not enabled — or if the zone doesn't match — SwaNET-exported recharge imported via the From SW Model Output recharge-source option will silently land in the wrong geographic location. The safest workflow is to enable UTM Only before drawing the domain; converting a non-UTM model after the fact is error-prone. Verify alignment visually by overlaying the imported recharge field against the IGW-NET domain boundary before trusting results.

22.3 Grid and Discretization Traps

The grid is how continuous physics meets finite computation. When the grid's resolution doesn't match the scale of the processes you're trying to simulate, the solver faithfully computes something — but it's not the answer to your question. A related set of traps comes from the interaction between geologic geometry and the grid: steep bedrock slopes, sharp pinchouts of sediment, and rapid thickness variations can create numerically ill-posed configurations where a physically reasonable setup produces dry cells, divergence, or hypersensitive K calibration. These traps are particularly easy to miss because the simulation runs cleanly — or, worse, runs cleanly at one K value and goes unstable at another equally plausible value. The first entry in this section (§22.3.1) addresses a foundational question that precedes all the others: how should vertical discretization be laid out in the first place — along geologic layers (the pragmatic MODFLOW approach) or in strictly horizontal rectilinear slices (the classical finite-difference approach) — and what goes wrong when the geology-based choice meets complex geology.

22.3.1 Geology-based vertical discretization — when the pragmatic choice breaks down

Symptom: a model with vertical layers aligned to geologic units runs cleanly across most of the domain but develops one or more of the following in locations where geology is complex: (a) cells that go dry or produce spurious heads where a geologic layer pinches out to near-zero thickness, (b) numerical instability — non-convergence, wildly oscillating iterations, or sensitivity to small K changes — concentrated in zones where a layer becomes very thin, (c) transmissivity contrasts across adjacent cells so large that the solver behaves as if approaching a singularity, or (d) head results that are physically reasonable at the regional scale but nonsensical in the specific areas where two geologic units converge or where a unit thins against bedrock.

Cause: this pitfall concerns the philosophy of vertical discretization, which is a choice that precedes almost every other discretization decision. Two philosophies have historically competed in groundwater modeling:

Strict rectilinear finite difference. Classical finite-difference theory assumes a rectilinear grid — horizontal layers of uniform or slowly-varying thickness, with cell faces aligned to coordinate axes. This is mathematically most rigorous: the finite-difference approximation of Darcy's equation converges reliably, error analysis is clean, and the numerical properties are well understood. Applied to an aquifer with realistic stratigraphy, however, strict rectilinear discretization requires many horizontal layers to adequately resolve the geologic structure — the K contrasts that matter for groundwater flow are concentrated at geologic boundaries, and if those boundaries are oblique to horizontal, many thin rectilinear layers are needed before any one of them cleanly represents a given geologic unit. The number of layers required grows quickly with structural complexity, and the computational cost grows with it.

Pragmatic geology-based layering. MODFLOW and most modern groundwater codes instead allow layer-conforming vertical discretization — model layer boundaries follow the tops and bottoms of geologic units rather than strict horizontal planes. This is not rigorously rectilinear in the classical sense, but it is enormously more efficient: the K structure that actually controls flow is resolved with a handful of model layers instead of dozens, and the solver addresses the physics at the natural scale of the geology. What this approach gives up in formal mathematical rigor, it recovers many times over in applied accuracy — because the error from having too few layers to represent a geologic contrast is usually far larger than the error from having non-strictly-rectilinear cells. For practical groundwater modeling, especially given the horizontal-to-vertical scale disparity in most aquifers (where flow within a layer is approximately horizontal), the pragmatic choice is clearly correct and is why MODFLOW is so widely used. IGW-NET follows this same pragmatic approach (Ch. 10).

The pragmatic choice comes with its own failure modes, and this is where the pitfall lives. Geology-based layers work cleanly when the layers are fundamentally approximately horizontal at real scale — that is, when layer thickness and layer dip vary smoothly enough that each model cell has a well-defined, non-tiny thickness and the horizontal-flow-within-layer assumption remains valid. When geology becomes complex — when the natural layering loses clean structure, when a unit pinches out to near-zero thickness over a short lateral distance, when two units converge and the layer between them becomes locally "singular" (effectively vanishing), or when a layer that averages 20 m thickness regionally has a spot where its thickness drops to 20 cm — the geology-based approach runs into numerical difficulties. Near-zero-thickness cells produce near-infinite vertical conductance and near-zero horizontal transmissivity, which the finite-difference matrix cannot reasonably represent. The system becomes singular or near-singular locally; the solver struggles; results in that neighborhood can be nonsensical even if the rest of the model is fine.

This is the pitfall: a modeler who has grown accustomed to the efficiency of geology-based layering applies it to a geologically complex area where the layers are not fundamentally approximately horizontal at real scale, and the resulting model has local singularities that compromise solution quality. This is distinct from the geometry-driven ill-posedness pitfall in §22.3.7 (bedrock slopes, pinchouts against bedrock) but closely related — §22.3.7 is about what happens at the bottom of the modeled aquifer, while §22.3.1 is about what happens between geologic layers within the aquifer stack.

Fix:

Ensure layers are fundamentally approximately horizontal at real scale before adopting geology-based discretization. Look at cross-sections of the geologic model before finalizing the layer structure. Does each layer maintain a reasonable minimum thickness everywhere in the domain? Are the layers locally coherent, or do they become fragmented, pinching out, or crossing each other in ways that geology-based layering cannot cleanly represent? If the geology does not support this assumption, the geology-based approach is not appropriate — consider using either more sublayers within a single geologic unit (to give the solver more flexibility), or reverting locally to a rectilinear approach in the problematic area.
Use IGW-NET's minimum-layer-thickness protection. Because near-zero-thickness layers are such a common failure mode, IGW-NET includes an option to enforce a minimum layer thickness — by default, any layer that would be thinner than a configurable percentage of the maximum layer thickness is automatically set equal to that minimum. This prevents the local singularities that would otherwise destabilize the solver while preserving the overall geology-based layer structure. The default percentage is usually appropriate; increase it if the solver is still struggling, or decrease it if the default is over-regularizing away real thinness that needs to be modeled. See Ch. 10 for the minimum-thickness configuration.
Add more sublayers within geologic units where the geology is locally complex. If a geologic layer locally becomes very thin or has complex internal structure, adding sublayers within that geologic unit gives the solver more vertical resolution where it is needed, without abandoning the geology-based approach elsewhere. This is a middle-ground between pure geology-based and pure rectilinear — each geologic layer is divided into a user-specified number of model sublayers, so the model retains the geology-conforming advantage while providing extra resolution where local complexity demands it. See Ch. 10 §10.2.3 for the “add sublayers when” list.
Simplify the geologic model in problematic areas if the complexity does not change the answer. If a layer pinches out locally in a way that is causing numerical trouble but the pinchout is not central to the question being asked (for example, it happens tens of kilometers from the area of interest), consider simplifying the geologic representation in that area rather than fighting the discretization. Geologic realism in the model should be held to the standard of what affects the answer to your question (see §22.4.7 on the complexity fallacy), not to the standard of matching every known detail of the geology.
Validate results in locally-complex areas specifically. The rest of the model may be numerically sound while a specific zone has hidden numerical trouble. Make a habit of examining heads and flux outputs in the neighborhood of known complex geologic features — pinchouts, unit convergences, thin zones — and look for signs of local nonsense: oscillating heads between adjacent cells, vertical flow patterns inconsistent with the regional flow system, mass-balance errors concentrated in those cells. If the local geology is complex and results look physically odd there, the discretization is the likely cause.

The underlying principle: geology-based vertical discretization is the right pragmatic choice for almost all applied groundwater modeling because the efficiency gain over strict rectilinear discretization far outweighs the loss of formal mathematical rigor — but it rests on an assumption that geologic layers are fundamentally approximately horizontal at real scale. When that assumption breaks down because geology is locally complex, the modeler must either adjust the discretization (add sublayers, enforce minimum thickness, simplify the geologic representation locally) or accept that results in the complex area need extra scrutiny. See Ch. 10 for the full layering framework, Ch. 7 §7.4 for the water-table-as-top formulation, and §22.3.7 for the closely-related geometry-driven ill-posedness that occurs at the aquifer bottom rather than between layers.

How to check: open the Layering dialog and review the configured layers and their sublayer counts. Open a cross-section view that cuts through areas of known geologic complexity (pinchouts, unit convergences, bedrock valleys). Look for layers that become very thin or vanish locally — these are the locations where geology-based discretization can become locally singular. In the same dialog, verify the Minimum Layer Thickness option is enabled (this auto-sets any layer thinner than a configurable percentage of the maximum thickness to that minimum, preventing near-zero-thickness cells from destabilizing the solver). Also check the BotE source in the Data Center BotE dialog — if bottom elevation came from a noisy raster, smoothing it before import often resolves local numerical pathologies.

22.3.2 Random-field correlation scale vs grid resolution mismatch

Symptom: a stochastic K-field simulation looks almost identical to the uniform-K case. Macrodispersion doesn't emerge. The plume spreads like advection–dispersion in a homogeneous medium, with no fingering, no channeling.

Cause: the grid is not resolving the correlation scale of the random field. If the correlation length Lambda is smaller than a few cell widths, the random field is effectively averaged by discretization — cell-to-cell K variations get filtered out, and the heterogeneity that theory says should produce macrodispersion isn't visible in the solution.

Fix: ensure at least 5–10 cells fit within one correlation length. For a 10 km domain with Lambda = 50 m, you need cells ≤ 10 m, which means NX × NY around 1000 × 1000. If that's too expensive globally, add a refined submodel over the area where heterogeneity matters (Ch. 13). See Ch. 19 §19.3 for full random-field setup and Ch. 19 §19.7.4 for the macrodispersion treatment.

How to check: open the K Attribute dialog and confirm the K source is set to “Random Field.” Then open the Random Field configuration (Ch. 19) and note the correlation scale. Compare it against the model cell size shown in the Info panel. The correlation scale should span several cell widths for the random K field to express meaningful heterogeneity — if correlation scale is smaller than cell size, the field is effectively averaged to uniform-looking K within each cell, defeating the purpose of the stochastic representation.

22.3.3 The 10-cell buffer violation in submodels

Symptom: results near the submodel edge look "too clean" or over-constrained. Wells' drawdown inside the submodel hits the boundary and rebounds unphysically. Plumes exit the submodel boundary and disappear.

Cause: the submodel boundary is too close to a feature of interest. Parent-model heads are imposed on the boundary; if the boundary is within the zone of influence of a well, a pumping-induced cone of depression, or a plume, the parent's relatively coarse heads contaminate the submodel's fine-scale solution.

Fix: the 10-cell buffer rule — submodel boundaries must be at least 10 cells away from any feature of interest. Expand the submodel until the buffer is satisfied, or refine the parent instead of using a submodel if the feature of interest is too central to separate. See Ch. 13 §13.3.

How to check: in the submodel view, measure the distance from the submodel boundary to your area of interest in model cells (count cells on the map, or divide real-distance by cell size from the Info panel). If the distance is fewer than 10 cells, boundary inheritance artifacts may propagate into the area of interest. Fix: enlarge the submodel outward until 10+ cells of buffer sit between the boundary and the area of interest, or relocate the area of interest toward the submodel center.

22.3.4 Sublayer count too low for vertical processes

Symptom: a partially-penetrating well produces the same drawdown pattern it would in a fully-penetrating configuration. A DNAPL source doesn't show the expected downward sinking. Leakage from a leaky river is uniformly distributed vertically instead of concentrated near the bed.

Cause: the conceptual layer has only one computational sublayer — the model is effectively 2D vertically, so vertical head gradients and vertical transport are averaged out. The physics that requires resolved vertical gradients is invisible in the solution.

Fix: add sublayers where vertical resolution matters. For partially-penetrating wells, 3–5 sublayers resolve the screen interval. For DNAPL sinking or leaky-river bed dynamics, more may be needed. Every added sublayer raises cost. See Ch. 10 §10.2.3 for the "add sublayers when" list and Ch. 7 §7.4 for cost implications.

How to check: open the Layering dialog and note the sublayer count per geologic layer. If your question involves vertical flow physics — plume sinking through dense zones, partial-penetration wells, well-screen depth effects, vertical head or concentration gradients, 3D particle tracking — one sublayer per geologic layer is usually not enough; 3–5 sublayers are a typical minimum. The sublayer count applies to the computational (ijk) grid produced at SIMULATE time; open a vertical cross-section after simulation to confirm you have the vertical resolution the question requires.

22.3.5 Horizontal refinement without matching vertical refinement

Symptom: a refined submodel resolves horizontal flow patterns beautifully but the multi-layer behavior at the refined scale still looks like the parent's coarse vertical representation.

Cause: refining NX/NY in the submodel without refining the sublayer count leaves the vertical resolution unchanged. The solver sees a fine horizontal grid and a coarse vertical grid — high aspect-ratio cells that misrepresent vertical physics.

Fix: horizontal refinement usually requires vertical refinement. If you refine cells horizontally by a factor of 4, consider increasing sublayers to keep cell aspect ratios reasonable. See Ch. 13 and Ch. 7.

How to check: compare the Grid settings (horizontal cell size and number of cells) against the Layering sublayer counts. If you increased horizontal resolution (e.g., reduced cell size by 5×) without proportionally increasing sublayers, the refinement is isotropic-horizontal only — vertical gradients that mattered at coarse resolution are still invisible to the refined model. Refine both dimensions together for problems where 3D physics matters.

22.3.6 Submodel-as-just-zoomed-in misconception

Symptom: you build a submodel expecting it to show "the same thing the parent shows, but with more detail," and instead the submodel produces quantitatively different results — different heads, different flow directions near some features.

Cause: a submodel is not a zoom — it's a separate simulation with BCs inherited from the parent. The submodel resolves physics the parent couldn't see (heterogeneity, small features, local pumping). Where those effects matter, the submodel's answer is different from the parent's, and the submodel is the correct one.

Fix: this is not a bug; it's the whole point of nested modeling. The refined submodel is revealing scale-dependent features the parent averaged over. See Ch. 13 for the full hierarchy-as-framework treatment.

How to check: open the submodel's Info panel and compare cell size, sublayer count, K source, recharge source, and stream representation against the parent model. A submodel that has only finer cell size but otherwise identical parameters to its parent is under-using the submodel mechanism — the submodel should also have refined K (often from local borehole data via T-PROGS or via imported raster), refined recharge where meaningful, and possibly finer sublayers. The per-property source selector in each attribute dialog is how you override inherited values locally.

22.3.7 Geologic-geometry-driven ill-posedness (thin, steep sediments over rising bedrock)

Symptom: a physically reasonable model setup produces one of several numerical pathologies: (a) dry cells appearing along valley edges or bedrock highs where the physical aquifer is actually wet but thin; (b) the solver failing to converge or needing many more iterations than a comparable flatter-bedrock setup; (c) divergent solutions where cells oscillate between wet and dry across iterations; (d) a solution that is stable at one K value but goes divergent with a modest (factor of 2 or 3) change in K, even though both K values are physically plausible. The model isn't wrong as-specified — it's close to the edge of what the numerical formulation can handle.

Cause: this is one of the most important interactions between geology and numerics in groundwater modeling. When the aquifer-bottom (bedrock, confining unit, or geologic contact) is locally steep or the sediment thickness varies sharply — a common situation near bedrock valley walls, buried channel edges, glacial terrain with rapidly-varying till thickness, or thin sediments pinching out against rising bedrock — the saturated thickness of cells along that geometry varies by orders of magnitude over short distances. In the finite-difference formulation, flow between adjacent cells is proportional to transmissivity T = K × saturated thickness, so neighboring cells can have transmissivities differing by factors of 10–1000, producing a stiff system with very poor conditioning. On top of that, the water-table-as-top formulation (Ch. 7 §7.4) requires the computed head to stay above the bedrock bottom for the cell to remain wet; on a steep bedrock slope, even small errors in head (a few cm) can push computed head below cell bottom, and the cell dries up. Whether that drying is physical (the cell really would be unsaturated in reality) or numerical (the real cell is wet but the discretization forces it dry) depends on geometry, K, recharge, and lateral flow in a way that is not easy to diagnose by inspection.

This creates a dangerous sensitivity. A setup with 2 m of sediment over a 10° bedrock slope under 8 in/yr recharge might be stable at K = 1 m/day (enough conductivity to move recharge laterally along the slope without over-building the water table) but go divergent at K = 0.3 m/day (water can't drain fast enough, water table tries to rise above surface, numerics break down) — or also at K = 10 m/day (water drains too fast, saturated thickness in the thin region becomes near-zero, cells dry). The physically correct answer exists somewhere in this range, but the numerically stable answer may occupy only a narrow window. Calibrating K "by convergence" is the wrong criterion — you can easily calibrate to a numerically convenient but physically meaningless K.

Fix: several mitigations, in order of effort and typical effectiveness:

Smooth the bedrock geometry locally where extreme gradients don't correspond to real geologic knowledge. If your bedrock raster has artifacts — DEM noise, interpolation bull's-eyes around a single borehole — smooth them. The real geology at the valley edge is rarely as sharp as a raster cell's worth of elevation change.
Refine the grid in the transition zone. A refined submodel (Ch. 13) over the bedrock-edge region resolves the geometry better, reduces the transmissivity-contrast stiffness, and gives the solver more cells over which to resolve the water-table gradient.
Use sublayers near the transition. Multiple thin sublayers of the same material give the solver more degrees of freedom to resolve vertical head gradients in the thin-sediment region, which can prevent spurious drying (Ch. 10).
Verify K against the geometry. A useful sanity check: compute the lateral transmissivity required to pass the recharge from the recharge area to the discharge boundary through the thin-sediment region. If the required T is physically unreasonable (K orders of magnitude outside geologic plausibility), the geometry is wrong or the recharge is wrong — no K value will make the model both physically correct and numerically stable.
Check for genuine vs numerical drying. If cells dry up in a region where independent evidence (water-well records, wetland presence, springs) indicates the aquifer is actually wet, the drying is numerical. If the region also lacks physical evidence of saturation, it may be genuinely dry and the model is correctly representing it.
Avoid calibrating K to stability. Auto-cal can easily pick K values that happen to land in a numerically stable window but are physically wrong. Calibrate against independent flux observations (Ch. 18) to pin K to physical reality, and fix the geometry or grid separately if stability is still an issue.

See Ch. 5 for the bedrock-elevation customization ladder (including smoothing options), Ch. 7 §7.4 for the water-table-as-top linearize-within/iterate-between formulation, Ch. 10 for layering and sublayers, and Ch. 24 for solver behavior near ill-posed geometries.

How to check: open a cross-section cutting through the geologic feature of concern — bedrock valley wall, pinchout against rising bedrock, thin sediments over steep basement. Look at the BotE (bottom-elevation) surface: are there rapid elevation changes over just a few cells, single-cell spikes, or sharp discontinuities that suggest a noisy raster? If so, smoothing the BotE raster before import often resolves the numerical pathology. Then in the K Attribute dialog, vary K by a factor of 2–3 around your best estimate and re-simulate — if the solution becomes divergent, that sensitivity is the signature of the narrow stability window described above. Calibrate against flux observations (baseflow from stream gages), not against solver convergence.

22.4 Calibration and Solver Traps

Calibration can compensate for parameter errors but cannot fix structural errors in the conceptual model. Most calibration pitfalls arise when parameter adjustment is used to paper over a problem that lives at the structural level. Solver-tuning pitfalls are a related family — when a solver is struggling, the instinct to tighten tolerances or bump iterations often masks an underlying well-posedness issue. Two fundamental-level pitfalls run through this section and deserve particular attention: a philosophy trap — the assumption that more processes simulated, more parameters, more zones produce a more accurate model, when for applied prediction work the opposite is usually true (§22.4.7) — and an applicability trap — the use of IGW-NET for aquifer systems where Darcy's law does not apply at the scale of the question being asked (§22.4.8). Both are easy to miss because they occur before any parameter choice has been made, and no calibration procedure can recover from them.

22.4.1 Head-only calibration hiding flux error

Symptom: calibration converges cleanly with small head residuals everywhere — the "fit" looks excellent — but the simulated flow system doesn't match what's known about the site: pumping wells produce less drawdown than observed, a stream shows less baseflow than gaged, a spring's discharge is wrong by an order of magnitude.

Cause: heads are insensitive to K-recharge tradeoffs over broad ranges. For a given hydraulic gradient, you can fit the heads with high K + high recharge, or with low K + low recharge; both produce the same head field but wildly different fluxes. Head-only calibration picks some point in this continuum with no mechanism to pick the physically correct one.

Fix: always calibrate against flux observations in addition to heads — stream baseflow, well-pumping drawdown measured independently, spring discharge. Even one flux observation constrains the K-recharge non-uniqueness enormously. See Ch. 18, particularly §18.2 on the non-uniqueness problem and §18.5 on multi-observation calibration.

How to check: open the Calibration dialog and review the list of observation types being fit. If only head observations are registered, flux information (stream baseflow, lake levels, pumping drawdown) is not constraining the calibration — parameters can adjust freely to fit heads while producing wildly wrong fluxes. Add flux observations (baseflow from USGS stream gages, measured pumping response, lake-level time series) before trusting the calibrated model. See Ch. 18 for observation-type recommendations.

22.4.2 Massive mass-balance errors despite excellent head match

Symptom: head residuals are tiny across all observation points — the model appears extraordinarily well-calibrated — but the global mass balance report shows a discrepancy of several percent or tens of percent. Inflow and outflow totals don't match the recharge you specified. A plot of cumulative error over time drifts steadily upward. The model looks right cell-by-cell but fails the most basic conservation-of-mass test.

Cause: this is one of the most insidious failure modes in groundwater modeling because the conventional calibration metric (head-residual statistics) is blind to it. Several mechanisms can produce it, often in combination: (1) an overly relaxed solver tolerance that lets the iterative solution "almost-converge" cell-by-cell while accumulating residual budget errors globally; (2) a poorly-posed boundary-condition configuration where two BCs are in tension at the same location and the solver finds a head solution that minimizes local head residuals while producing spurious flux into or out of the domain; (3) a free-surface problem where cells becoming wet or dry during transient simulation introduce mass-balance artifacts the solver cannot fully correct; (4) in coupled lake-aquifer or coupled watershed runs, time-step sizes or sub-stepping schemes that drift from conservation over long simulation periods. The danger is that the head-match deceives the modeler (and the client, and the reviewer) into trusting the model; the mass-balance error is then used as the basis for decisions about water budget, sustainable yield, or plume mass — all of which are the quantities the model is most wrong about.

Fix: always read the mass-balance report; head-match alone is never sufficient validation. IGW-NET reports global and per-feature mass balance for every simulation — make inspecting it part of every workflow, not an optional check. If the global discrepancy is above ~1% for a steady-state run or above ~2–3% cumulative for a transient run, investigate before trusting any derived quantity. Common fixes depend on the cause: tighten solver tolerance in Ch. 24 (but note §22.4.4 — if tighter tolerance doesn't close the balance, the problem is structural); resolve conflicting BCs at the same location; for transient runs with coupled SW, reduce time step and enable sub-stepping in Ch. 15; for unconfined problems with dry-cell risk, use the water-table-as-top formulation rather than fully-3D unconfined. See Ch. 8 for reading the mass balance as diagnostic conversation and Ch. 25 for diagnostic signatures of structural budget errors.

How to check: open the Flow Summary or Water Balance panel after simulation. The panel lists Ins (recharge, lateral inflow, prescribed-head inflow, lake injection) and Outs (baseflow to streams, lake discharge, pumping, evapotranspiration, prescribed-head outflow). The difference between Ins and Outs is the mass-balance error — should be <1% of the total. If larger, inspect each BC type for asymmetric configuration: constant-head inflow with no balancing outflow, pumping with no replenishment, lake coupling with inconsistent Ins/Outs. Review pumping wells in the Wells dialog, recharge source in the Data Center, and lake water-balance components in the Lake dialog for the most common asymmetries.

22.4.3 Structure-first principle violation

Symptom: calibration is producing K values that are physically implausible (orders of magnitude off from any geologic estimate), or parameter ranges are unstable across similar runs.

Cause: the conceptual model has a structural error — wrong layering, missing feature, wrong boundary type — and the auto-cal is trying to compensate by pushing parameters to implausible values. Calibration cannot fix structural errors; it can only hide them temporarily.

Fix: stop and re-examine the conceptual model. Does the layering match the site hydrostratigraphy? Are all significant features represented? Are boundaries the correct type? Only after structure is sound should parameter calibration proceed. See Ch. 18 §18.1 for the structure-first principle and Ch. 25 for diagnostic signatures of structural error.

How to check: open the model Info panel and read through what the model says about its own structure — K source (Constant / Data Center / Random / T-PROGS / Import), BotE source, recharge source, stream representation level, lake coupling level, solver choice. For each entry, ask: is this the structure I chose deliberately, or a default I haven't inspected and justified? If the model inherits most of its structure from defaults, structure hasn't actually been chosen yet — parameter calibration on top will tune values against an unexamined structural framework.

22.4.4 Overfitting in auto-cal

Symptom: calibration produces tight head-residual statistics at observation points but the simulated flow field between observations looks physically strange — oscillating contours, implausible gradients, sharp transitions that don't correspond to any geologic feature.

Cause: too many parameters relative to observation data. Auto-cal is exploiting degrees of freedom to fit noise in the observations, producing a model that fits the data points but is physically wrong between them.

Fix: start simple — calibrate K and recharge only (the big two), and only add parameter complexity if initial calibration can't match observations. Regularize if needed (see Ch. 18). Signs of overfitting: residuals near zero at obs points but contours that "dip" unphysically between them.

How to check: count the adjustable parameters in your calibration (one per zone for zoned K, one per recharge zone, etc.) and compare against the number of independent observations. Zones are defined in the Zone dialog — each zone with its own calibrated K is a degree of freedom. If the model has 20 K zones and the calibration fits 15 head observations, the system is under-determined and automatic calibration will fit noise. Reduce zone count, simplify zoning to match data density, or add observations. See Ch. 18 for the parameter-vs-observation rule.

22.4.5 Solver tuning masking a well-posedness problem

Symptom: the solver fails to converge, you tighten tolerance or increase iterations, it "converges" but the results look nonsensical; or you change solver type (PCG → SIP → SOR → Hybrid SOR) hoping one will work.

Cause: the problem is usually not the solver; it's the problem specification. Isolated subdomains with no hydraulic connection, overlapping constant-head features with conflicting values, domains that are insufficiently constrained by BCs (resulting in a null space in the solution) — these all produce ill-posed systems that no solver can fix.

Fix: when solver tuning doesn't work quickly, stop and look at the problem structure. Does every cell have a path to a BC? Are there conflicting BCs at the same location? Is the system over-constrained or under-constrained? Is the bedrock geometry creating a stiff transmissivity contrast (see §22.3.7)? Is a geologic layer becoming locally singular or too thin (see §22.3.1)? See Ch. 24 for the "when numerical adjustments don't work, the problem is structural" principle.

How to check: open the Solver Options dialog. If iteration counts are unusually high, convergence tolerances have been loosened, or damping has been increased well beyond defaults, the solver may be masking a well-posedness problem rather than solving the right one. Before tightening settings further, step back and inspect problem structure: are BCs consistent, does every cell have a path to a BC, is bedrock geometry smooth (§22.3.7), are geologic layers approximately horizontal at real scale (§22.3.1)? Solver settings are the last place to fix a structural problem.

22.4.6 Transient initial-condition mismatch

Symptom: a transient simulation shows a large artificial transient in the first few time steps — heads swing dramatically, flux shows a spike — before settling into the expected behavior.

Cause: the initial condition (starting heads) is not consistent with the initial stresses. If the transient simulation starts with uniform heads but the boundary conditions and stresses would produce a very non-uniform steady state, the solver spends the first time steps "relaxing" the initial condition toward the appropriate state.

Fix: run a steady-state simulation first to produce a physically consistent starting head field, then use those heads as the initial condition for the transient. This is the default recommendation in Ch. 9. The exception: when the initial state you care about is genuinely non-equilibrium (e.g., just after a large rainfall event), you may want an explicit non-steady-state initial condition.

How to check: in the Simulation dialog, review the initial-condition configuration for transient runs. Options include uniform head, output from a prior steady-state run, or an imported head raster. If the IC is uniform head but the aquifer has significant head gradients, the transient simulation will spend real simulated time adjusting from the unphysical IC — results during that adjustment period are not responses to the real transient signal. The best practice is to run a steady-state spin-up first, then use its output head field as the transient IC; both options are selectable in the Simulation dialog.

22.4.7 “More processes simulated equals more accurate” — the complexity fallacy

Symptom: a modeler who has access to a capable platform reaches for every available feature — transient dynamics, T-PROGS geostatistics, reactive transport, stochastic Monte Carlo, density-dependent flow, coupled surface water, dozens of calibration zones — for a project whose objective could be answered by a much simpler model. The resulting model is expensive to run, hard to analyze, impossible to present cleanly, and produces predictions with wider (not narrower) uncertainty than a simpler treatment would. Often the modeler can't clearly explain why any given piece of complexity is there, only that “more realistic” must be better.

Cause: this is the complexity fallacy, and it is one of the most widespread and costly pitfalls in professional groundwater modeling. The assumption that more processes simulated makes the model more accurate is wrong in both theory and practice. Every process added brings parameters; every parameter brings uncertainty; the uncertainty compounds. A model with four zones can be well-constrained by twenty observations; the same twenty observations cannot constrain a model with forty zones, and auto-cal will find parameter values that fit the data without physical meaning (see §22.4.4 on overfitting). More fundamentally, modeling is a goal-directed activity, not a physics-reproduction exercise. The right level of complexity depends entirely on what question the model is answering:

If the goal is to understand processes — how heterogeneity affects transport, what vertical circulation looks like under a given water-table shape, how coupled lake-aquifer dynamics work in principle — then adding complexity can be appropriate. But even here, do not add it all at once. Start simple, add one process at a time, observe the effect, understand what changed. The real-time modeling capability of IGW-NET (§1.3.2, synthetic mode) is designed for this incremental workflow.
If the goal is to predict the future — sustainable yield, wellhead protection area, plume arrival time, drawdown from proposed pumping — use the simplest conceptual model that can achieve the goal. This is often called KISS (Keep It Simple, Stupid) or the principle of parsimony, and it is the most important modeling virtue for applied work. If you can explain the available observations with a constant-K, constant-recharge model, do not introduce zones. If you add zones, they should be supported by clear independent evidence (geologic maps, pumping test results, borehole logs) — not by the desire to match local observations that the homogeneous model missed. Local patches that exist only to force the model through a data point are the calibration equivalent of overfitting a polynomial through noisy data: tight fit, poor predictive skill.

The parsimony principle has a sharper formulation that is particularly useful: zones are hypotheses about subsurface heterogeneity, and each hypothesis needs independent evidence. If your conceptual model requires a low-K zone to explain a high head at one observation, ask what geologic evidence — a mapped clay unit, a documented till layer, a consistent pattern across multiple wells — supports it. If the only “evidence” is that the head matches better when you add the zone, you have not explained the observation; you have fitted it. The distinction matters because hypotheses supported by independent evidence generalize to predictions (the model is more trustworthy out-of-sample); hypotheses introduced solely to fit observations do not.

Fix: treat simplicity as a first-class modeling virtue and make complexity a justified addition. Practical steps:

State the modeling objective explicitly before building. What question is this model answering? What decision will be made from the answer? What level of precision does that decision actually need? Models built for regulatory precedent, litigation, or high-stakes infrastructure decisions have different complexity budgets than models built for classroom demonstrations or first-cut scoping.
Start with the simplest conceptual model plausibly consistent with the objective. For most site-scale prediction questions, a homogeneous or 2–3 zone model with steady-state flow and advective particle tracking is a legitimate starting point, even if the site is known to be heterogeneous in detail. Add complexity only when the simple model demonstrably cannot answer the question.
Use synthetic mode incrementally when exploring complex processes. Add one complexity at a time (a low-K zone, then heterogeneity, then transport, then coupling) and observe its effect on the solution. This builds intuition and helps you recognize which complexities actually matter for your question (Ch. 1 §1.3.2).
Require independent evidence for each added feature. Before adding a zone, a new BC, a new process, ask: what evidence (beyond the fit to existing observations) supports this? If the answer is “none,” the feature is overfitting, not modeling.
Evaluate the simpler model honestly before replacing it. If a homogeneous-K model produces a 5-foot head residual at one well while the aquifer is characterized by 50 feet of regional head variation, that 5-foot residual may be acceptable given the goal. Not every residual needs to be explained — the ones that matter are the ones that would change the answer to your question.
Document the complexity budget. For every feature in the model, a one-line justification: “zone A exists because borehole logs show till from 20–45 ft depth across three wells.” “Transient dynamics are included because the prediction is for seasonal pumping.” If you can't write the justification, the feature probably shouldn't be in the model.

See Ch. 1 §1.3 for the refine-progressively philosophy, Ch. 5 for the customization ladder (defaults first, refine where evidence supports it), Ch. 18 for structure-first and overfitting, and Ch. 27 for how different professional roles (consulting, research, teaching) carry different complexity budgets.

How to check: open the model Info panel and scan the list of active processes and properties. For each enabled feature — retardation, decay, half-life, sorption, random fields, T-PROGS, multi-level surface-water, coupled lakes, variable density, reactive species — ask: is there independent data supporting this, or is it stored-as-default because I never disabled it? Each active process is a place where complexity was added, and parsimony demands that each have an explicit justification tied to the question being asked. The Biochemistry dialog, the Concentration Source dialog, and the surface-water-level settings are the most common places where unjustified complexity accumulates.

22.4.8 Using IGW-NET where Darcy's law does not apply — the karst/fracture/cave trap

Symptom: a modeler applies IGW-NET to an aquifer system dominated by karst features, large fractures, conduits, or caves, and finds that (a) no parameter set can reproduce the observed flow behavior — pumping tests show responses that IGW-NET cannot simulate regardless of K values, (b) tracer-test velocities are tens or hundreds of times faster than anything a continuum Darcy model can produce, (c) water-level response to rainfall events shows timescales (minutes to hours) inconsistent with diffusive Darcy flow, or (d) the model calibrates reasonably to steady-state heads but fails catastrophically on any transient or transport prediction. In some cases the model runs and produces numerical results that look like a valid simulation, but the physics underneath is wrong — the model is computing Darcy flow for a system that is actually dominated by conduit or turbulent flow.

Cause: IGW-NET is built on Darcy's law, which assumes flow through a continuous porous medium where velocity is linearly proportional to hydraulic gradient and the medium can be characterized by a representative elementary volume (REV) much larger than any individual pore but much smaller than the scale of heterogeneity being simulated. Darcy's law is the foundation of almost all quantitative groundwater hydrology because it works superbly for sediments, sedimentary rocks, weathered bedrock, and many other geologic materials. But it breaks down in specific contexts:

Karst aquifers — carbonate rocks (limestone, dolomite) where dissolution has produced conduits, caves, and enlarged fractures through which flow is often turbulent rather than laminar. Single conduits can carry more water than kilometers of surrounding matrix. No single K value represents the system.
Cave systems and large fractures — any flow through a discrete geometric feature large enough that the Reynolds number exceeds the Darcy threshold (Re ≳ 10) produces inertial or turbulent flow. Q is no longer proportional to gradient.
Fractured crystalline rock at site scale — granite, basalt, or metamorphic rock where transmissivity is concentrated in a sparse network of discrete fractures. At very small scales a discrete-fracture-network model is the appropriate tool; continuum models are not meaningful when only a handful of fractures exist.
Very-close-to-well turbulent flow — even in Darcy media, the region within a few well-radii of a pumping well can experience non-Darcy behavior if pumping rates are high enough; this is usually handled with well-loss coefficients rather than flagged as a modeling limitation.

Critically, however, Darcy's law applicability is scale-dependent, and this is where the trap becomes more subtle. A system that is plainly non-Darcy at one scale can be approximately Darcy at another. A karst aquifer with individual conduits that plainly violate Darcy's law at the scale of the conduit itself may behave like an equivalent porous medium when observed over a large enough region — kilometers to tens of kilometers — where the collective behavior of many conduits, matrix blocks, and fractures averages toward something that looks diffusive. Whether this equivalent-porous-medium (EPM) representation is useful depends entirely on what question you're asking. For a regional water-balance study, an EPM model of a karst aquifer may be perfectly adequate. For predicting the arrival time of a contaminant injected upgradient of a municipal well in karst, the EPM representation is useless because the contaminant follows a few conduits in hours to days, not the diffusive Darcy flow the model will compute over months to years. The issue is not whether Darcy's law “really” applies — it doesn't, at the conduit scale — but whether the averaged behavior at the scale of your question can be captured by an equivalent diffusive description.

Fix: before building an IGW-NET model, verify the aquifer system is appropriate for a Darcy-based continuum representation at the scale of your question:

Identify the aquifer material and structure from geologic information. Is the system dominated by porous sediments or sedimentary rocks? Fractured but small-aperture bedrock? Karst with mapped caves or sinkholes? The USGS Water Atlas, state geologic surveys, and regional hydrogeologic reports will tell you whether you're dealing with materials that typically support continuum modeling.
Look for karst indicators. Springs with large discharge (hundreds to thousands of gpm), rapid water-level response to storms, turbid water after rain events, sinkholes on the landscape, cave exploration records, or documented losing streams that disappear underground are diagnostic. Dye tracer studies showing transport velocities greater than about 10 m/day are strongly suggestive of non-Darcy flow.
Evaluate the question-to-scale match. If the system has known non-Darcy features but your question is at a scale much larger than the individual features, an EPM model may still be appropriate — with the explicit caveat that the model is a coarse-resolution descriptor of an averaged behavior, not a predictor of local flow paths. State this limitation explicitly in any reporting. If the question is at a scale comparable to or smaller than the non-Darcy features (wellhead protection in karst, cave-hydrology studies, contaminant transport along a known fracture zone), IGW-NET is not the right tool regardless of how the parameters are chosen.
Use appropriate tools for non-Darcy problems. Discrete-fracture-network (DFN) codes for fractured-rock systems at site scale; conduit-flow models (such as CFPv2, MODFLOW-CFP, or SWMM-based approaches) for karst where the conduit network is mapped; hybrid codes that couple Darcy matrix flow with conduit flow for intermediate cases. IGW-NET's role in these contexts is to answer the regional-scale questions where the EPM approximation is defensible, while explicitly acknowledging what it cannot address.
When EPM modeling of karst is defensible, calibrate transport separately from heads. The head field may be reasonably represented by an EPM; the transport behavior almost never is. Calibrate K and recharge to head and baseflow observations at the regional scale, but use particle tracking and tracer-test observations to understand the travel-time distribution — and present both together so the user of the model knows what the EPM representation covers and what it doesn't.

The honest summary: IGW-NET applies where Darcy's equation applies, but whether Darcy's equation applies is a question about scale, not about geology alone. For a porous sedimentary aquifer at any practical scale, Darcy's law holds cleanly and IGW-NET is the right tool. For a karst aquifer being studied at the regional water-balance scale, an EPM representation in IGW-NET may be defensible — but only as an averaged descriptor, with explicit acknowledgment of the limitation. For a karst, cave, or major-fracture system being studied at the scale where individual features matter, IGW-NET is not appropriate and a different modeling approach is required. Recognizing this boundary is part of the modeling judgment that separates competent applied work from misuse of the tool. See Ch. 1 §1.2 for the platform's explicit capabilities and limits, and Ch. 5 for the characterizations appropriate to different aquifer types.

How to check (framing, not detection): the IGW-NET interface cannot tell you whether Darcy's law applies — that is a question about the physical aquifer at the scale of your question, not about the model configuration. But the interface shows you what the model is claiming about the aquifer: open the K Attribute dialog to see K values and their source. If K is a single lumped value over a region known to have karst, major fractures, or conduit flow, and your question is at the scale where those features individually matter, the model is implicitly treating the aquifer as an equivalent porous medium — whether that assumption is defensible is a modeling judgment, not an interface setting. This pitfall lives in the conceptual model, not in any dialog; the interface shows the model's claims, and you compare those against what you know about the real aquifer.

22.5 Workflow and Platform-Integration Traps

Pitfalls that come from how the platform's features interact rather than from any single feature. Feature-precedence rules, engine selection, workflow interactions between chapters — these are easy to trip over when you're combining capabilities that each work fine in isolation. One entry in this section (§22.5.3 — numerical dispersion / numerical dilution) deserves particular attention because it lives at the boundary between workflow and transport physics: it is intrinsic to cell-averaged transport methods, it cannot be eliminated through grid refinement alone, and it has an asymmetric effect on predictions — conservative for plume extent but severely non-conservative for peak concentration at a receptor. For regulatory or health-risk work where peak concentration matters, this trap can quietly understate real risk while appearing to produce a physically-reasonable plume picture. The dual-simulation diagnostic described in that entry is worth adopting as routine practice for transport modeling.

22.5.1 DEM-as-drain vs DEM-as-stage confusion

Symptom: when you inspect a model's surface-water behavior, you can't tell whether streams are being represented at Level 1 (DEM-as-drain, surface drainage as a pressure-relief valve) or Level 2 (DEM-as-stage, where the DEM provides default stream stages for explicit stream polylines). The behavior doesn't match what you thought you'd specified.

Cause: the two uses of the DEM look similar on the surface but produce very different physics. Level 1 is always active (DEM-as-drain via DrnL is a default surface drainage mechanism). Level 2 requires explicit stream polylines, which may auto-populate based on stream order or may be manually drawn — and those polylines use DEM elevation as their default stage.

Fix: explicitly check the stream representation for features of interest. Ch. 14 covers the three-level framework in full; §14.2 is the DEM-as-drain canonical treatment. See also §22.2.4 for a related trap specific to DEM-as-stage: when narrow streams in variable topography produce spurious mass-balance flux because DEM cells don't resolve river geometry, and §22.2.5 for the broader DEM resolution-and-artifacts pitfall.

How to check: the two mechanisms live in different places. DEM-as-drain (Level 1) is configured in the DrnL (surface-drainage leakance) dialog — look for non-zero leakance values that indicate the mechanism is active. DEM-as-stage (Level 2) is configured in the Stream Polyline dialog with the stream-active flag enabled and per-polyline stage values defaulting to DEM elevation. Both can coexist and affect different reaches. Use the model Info panel to see which SW level is active for each feature of interest.

22.5.2 MT3DMS vs Magnet-native transport engine ambiguity

Symptom: transport results differ between runs that you thought used the same settings. A feature (e.g., dual-domain mass transfer) works in one run and is ignored in another.

Cause: IGW-NET supports two transport engines — MT3DMS (external, capability-rich) and a Magnet-native transport module (built-in, lighter-weight). Engine selection is encoded in the gmSolverID field via a colon-split convention, and different transport features are supported by different engines. Changing transport options can silently switch engines.

Fix: verify the active transport engine explicitly — the solver configuration dialog and the parser output both show it. For features that are engine-specific (dual-domain, reactive chains, specific sorption forms), check that the selected engine supports them. See Ch. 12 and Ch. 24.

22.5.3 Numerical dispersion swamping physical dispersion — and why it hides peak concentrations

Symptom: a contaminant plume looks sensible at a glance — a smooth, approximately elliptical shape that spreads as it migrates, with concentration gradients that appear physically reasonable. But as you compare across grid resolutions, the plume shape changes noticeably — it becomes sharper and more intense at finer resolution, or blurrier at coarser resolution. The apparent dispersivity you'd back-calculate from the plume shape is several times larger than the physical dispersivity you specified. In transport calibration, you find that the specified dispersivity that “works” is suspiciously smaller than site-characterization values would suggest. A particularly consequential variant: the predicted peak concentration at a receptor is substantially lower than measured values at real plumes of similar origin and age.

Cause: numerical dispersion (also called numerical dilution) is artificial smearing of concentration gradients introduced by the averaging effect of the finite-difference / finite-volume discretization — any concentration that enters a cell is immediately averaged across that cell's volume, and this averaging mimics real dispersion mathematically but is an artifact of the discretization, not of the physics. In areas with large concentration gradients, the cell-averaging effect is largest, and so numerical dispersion is worst precisely where concentration contrast matters most — at the leading edge of a plume, at the peak of a concentration plume near a source, and at the transition between clean and contaminated water. This is why numerical dispersion cannot be fully eliminated through grid refinement alone (though refinement reduces it): it is intrinsic to the cell-averaging philosophy of the method, not a bug. The smaller the Peclet number (Pe = cell size × velocity / dispersion coefficient) is not, the worse the smearing — and for many realistic site-scale problems with modest physical dispersivity and moderate grid resolution, numerical dispersion dominates physical dispersion by a factor of 5 to 50. The simulation is then computing physics the user did not specify.

A critical directional asymmetry makes this especially insidious: numerical dispersion tends to be conservative for predictions about plume extent — the plume looks bigger and reaches farther than it really does, which overestimates the area of influence and therefore overstates risk in that dimension — but severely non-conservative for predictions about peak concentration at a receptor. Because numerical dilution smears concentration across the cell, the peak concentration at any given location is lower than the physical peak would actually be. A model used to answer “how far will the plume reach?” tends to err on the safe side; the same model used to answer “what concentration will arrive at the water-supply well?” tends to badly underestimate the true value. For regulatory, cleanup, or health-risk applications where peak concentration matters, this underestimation can understate the real risk by a large factor while appearing to produce a physically-reasonable plume picture. This asymmetry — conservative for one question, dangerously non-conservative for the other, using the same model run — is why every transport model used for peak-concentration predictions needs explicit numerical-dilution diagnosis.

Fix: four independent checks, used in combination. The first two diagnose whether you have a numerical-dilution problem; the second two reduce it.

Grid-refinement convergence test. Run the same transport simulation at progressively finer grids (halving cell size at each step). If plume shape and peak concentration change substantially between refinements, numerical dispersion is dominating and the coarser grid is underresolved. Refine until results stabilize, or use a submodel (Ch. 13) over the plume region.
Dual-simulation diagnostic with particles and concentration together — the most direct visualization of numerical dilution. IGW-NET supports simulating both a concentration source and a matching particle source in the same model run. Set up an identical release zone: at that zone, release continuous particles (say, 1000 particles per time step) and simultaneously inject a corresponding mass flux as a concentration source. Run both together. The particles trace where mass really goes — particle tracking has no numerical dispersion because particles follow the velocity field directly rather than a cell-averaged concentration; and IGW-NET's particle tracking includes a Brownian-motion component that captures the real physical diffusion/dispersion (Ch. 11). The concentration plume shows what the transport solver produces. Overlaying the two reveals the numerical-dilution error directly: anywhere the continuous plume extends beyond or has lower peak than the particle cloud suggests, the difference is numerical dilution. This is the cleanest available diagnostic for how much your transport solution suffers from numerical artifacts, and it costs one extra simulation step. Make this a routine check for any transport model used for peak-concentration predictions.
Reduce numerical dispersion with unstructured grids, submodels, and solver-scheme choice. Three practical approaches, usable alone or together. Unstructured grids (MODFLOW 6 DISV, Ch. 20) let you refine specifically where concentration gradients are large — near the plume front, near pumping wells, near the source zone — without refining the whole domain. Submodels (Ch. 13) serve the same purpose at a different scale, letting you run a fine-resolution local transport model with boundary conditions inherited from a coarser regional flow model. Higher-order solver schemes, specifically MT3DMS with the TVD (total variation diminishing) or MMOC (modified method of characteristics) scheme, substantially reduce numerical dispersion compared to the default upwinding scheme. IGW-NET's default transport scheme is upwinding specifically because it is the most numerically robust — convergence is reliable, the solver rarely fails, and the simulation produces results you can see. Higher-order schemes give you less numerical dilution but are less stable and can produce oscillations, negative concentrations, or convergence failures that mask the real solution. The recommended workflow is to run the default first, confirm the solution is sensible, then re-run with a higher-order scheme to sharpen peaks and reduce dilution. Visualizing the solution first is a prerequisite for being able to interpret what a higher-order scheme produces — if you can't see the robust solution, you can't tell whether the higher-order solution's sharper features are real or oscillation artifacts. See Ch. 12 and Ch. 24 for transport engine and scheme selection.
Cell-size heuristic for advection-dominated transport. Aim for grid Peclet number ≲ 2 within the plume region, which means cell size ≲ 2 × (physical dispersivity). If physical dispersivity is 5 m, cells should be ~10 m or smaller in the plume region. For realistic field-scale problems where dispersivity is tens of meters, this is usually achievable with a refined submodel.

The underlying principle: numerical dilution is intrinsic to cell-averaged transport methods and cannot be fully eliminated, only reduced. Its asymmetric effect — conservative for plume extent, non-conservative for peak concentration — means a transport model's interpretation depends on what question is being asked. Particle tracking with Brownian-motion dispersion serves as the numerical-dilution-free reference against which continuous-concentration results can be calibrated. For any question where peak concentration matters — health-risk assessment, regulatory compliance at specific receptors, cleanup target concentrations — pair your transport simulation with particle tracking from the same source, visualize both, and treat the difference as a bound on the real numerical-dilution error. When dilution needs to be reduced, use unstructured grids or submodels to target refinement where gradients are largest, and step up to higher-order solver schemes (TVD, MMOC) only after the robust default solution has been visualized and confirmed as sensible. See Ch. 11 for particle tracking with Brownian-motion dispersion, Ch. 12 for the transport engine framework, Ch. 13 for submodel-based refinement, and Ch. 20 for unstructured grids.

How to check: open the Transport Solver dialog and note the scheme in use (upwinding / TVD / MMOC / HMOC). Upwinding is the default — robust but maximum numerical dilution. Then open the Dispersivity dialog and note the longitudinal (αL) and transverse (αT) values; if these are small (sub-meter scale) and your cells are tens of meters or larger, numerical dilution likely dominates physical dispersion. To enable the dual-simulation diagnostic: add a Particle Source Zone at the same location as your Concentration Source Zone before the next simulation, then overlay the resulting particle cloud against the continuous concentration plume. Any spreading of the continuous plume beyond the particle cloud, or lower peak concentration than the particle density would predict, is numerical dilution visualized directly.

22.5.4 Free-surface trap in 3D transport

Symptom: a 3D transport model with pumping produces mass-balance errors that grow over time, or produces plume migration patterns that violate conservation near a partially-unsaturated cell.

Cause: in 3D unconfined flow, cells that become dry during transient pumping create a moving free surface. Transport computed on a grid whose cells are changing saturation state is easy to get wrong. IGW-NET's design (water-table-as-top, linearize-within/iterate-between) avoids this trap in the flow solver, but transport inherits flow and can still show free-surface artifacts if the flow solution has them.

Fix: see Ch. 12 for the "water flowing IN" rule and the free-surface-trap treatment. Often the practical fix is to stay in 2D for transport when the question doesn't specifically require 3D (see Ch. 12's "2D first, 3D when the question demands it" tip), or ensure the flow model has stabilized before initiating transport.

How to check: in the Layering dialog, confirm the water-table-as-top formulation is active — this is IGW-NET's protection against the classical free-surface trap (Ch. 7 §7.4). Confirm in the Transport configuration whether 3D transport is coupled with the flow solution. The combination of a moving free surface with 3D transport is where this pitfall surfaces; the guidance in Ch. 12 (“2D first, 3D when the question demands it”) is worth following before committing to 3D transport on a transient unconfined problem.

22.5.5 Simulate re-initializes everything

Symptom: you set up an intricate particle-tracking configuration, click Simulate, and the particles appear to reset — positions, counts, options all return to defaults.

Cause: the Simulate button re-initializes transient state across the model. Particles that have been advected during a previous simulation are reset; particle configurations are preserved but their run state is not.

Fix: configure particle options before clicking Simulate, and use the Analysis Tools to export/save particle results before re-running. See Ch. 11.

How to check: this is a workflow behavior of the SIMULATE button. After running, open the model Info panel or a grid-view to see the computational (ijk) grid that was produced — it may differ from the conceptual (xyz) grid that was visible before SIMULATE. The xyz→ijk remapping is described in Ch. 7.

22.5.6 DrawDomain erases — one active domain at a time

Symptom: you want to add a new submodel or modify a domain, and find that drawing a new domain erased your existing one.

Cause: IGW-NET maintains exactly one active domain. DrawDomain replaces; it doesn't add. To nest a finer model inside an existing one, use submodels (zones with "submodel" flag, see Ch. 13), not a second DrawDomain operation.

Fix: for nested workflows, use the zone-based submodel approach. For replacing a domain, save the current model first.

22.5.7 Forgetting that Data Center dataset names aren't stored

Symptom: you open an older model and can't tell which specific K raster or recharge raster it was built against. The values are in the file but the named dataset they came from is not.

Cause: the parser stores the values that came from a Data Center raster (in the appropriate fields) but not the dataset name or ID. This is by design — the model should be self-contained once built. But it means the provenance of the data is not encoded in the model file.

Fix: document the data provenance in the model's metadata (there are fields for this), or in an accompanying README. See Ch. 23 for the full Data Center treatment.

Related chapters

Ch. 25 — Error Messages & Diagnostics — for simulations that fail to run or return errors (this chapter is for simulations that run but produce wrong results).
Ch. 18 — Calibration — full treatment of structure-first, K-recharge non-uniqueness, and overfitting.
Ch. 24 — Solvers and Solver Options — when numerical adjustments don't work, the problem is structural.
Ch. 14 — Surface Water as Boundary Conditions — the 10×-rainfall pitfall canonical treatment.
Ch. 13 — Nested Modeling — 10-cell buffer rule in full context.
Platform Reference pillar — field-level catalog for anyone needing deeper detail on specific model-file fields.