Multiple physical constraints improve groundwater contaminant model predictions

A Bayesian data-driven framework that couples several physical constraints into the correction of groundwater contaminant transport models improved predictive accuracy and physical consistency in two tests: a synthetic three-dimensional reactive-transport simulation and a laboratory sand-column experiment. The method combines observations, a Gaussian process representation of structural error, and separate likelihood functions for physical constraints, then adaptively balances those likelihood components during parameter identification. In the synthetic case, coupling carbon and chlorine mass conservation reduced mean mass-balance error rates to 23.63% and 16.49%, respectively.

In the sand-column experiment, combining cadmium mass conservation with a concentration-bound constraint lowered the mean cadmium error rate to 10.36% and raised the validation-period Nash–Sutcliffe efficiency for predicted concentration to 0.9419. Groundwater transport models often reproduce the general behavior of contamination while systematically departing from observations because their representation of a real aquifer is incomplete. Small preferential pathways may be omitted, chemical reactions may be simplified, or processes such as phosphate–cadmium cotransport may be ignored.

These structural deficiencies differ from parameter uncertainty: even when parameter values are estimated carefully, a misspecified model structure can generate persistent prediction errors. The study addresses that problem by retaining the physics-based model and adding a Gaussian process regression model for its structural error, while using physical constraints to prevent the corrected predictions from violating mechanisms that should govern the system. The Bayesian correction model treats the observed data as the sum of a groundwater simulation, a systematic structural error represented by Gaussian process regression, and random measurement error.

Physical-model parameters and Gaussian process hyperparameters are estimated jointly through Markov chain Monte Carlo sampling with the DREAMzs algorithm. Equality constraints, such as conservation of carbon, chlorine, or cadmium, are represented as zero-mean Gaussian sub-likelihoods. Inequality constraints, including a bound on cadmium concentration, are converted into equivalent equality forms.

A combined likelihood multiplies the observational likelihood by these constraint likelihoods, allowing observations and physical mechanisms to influence the posterior distribution together rather than leaving the error-correction model as an unconstrained statistical component. A central technical problem is that constraint likelihoods associated with different mechanisms can differ greatly in magnitude. In the synthetic experiment, their values differed by as much as six orders of magnitude, meaning that a large likelihood could overwhelm a smaller one; in the column experiment, the difference was approximately one order of magnitude.

The proposed inverse-magnitude weighting scheme calculates a mean sub-likelihood, assigns each component a provisional weight inversely related to its magnitude, and normalizes those weights. Weights are updated separately in each Markov chain during burn-in, after which a common fixed set is used for subsequent sampling. A sliding-window variance test and a chi-squared criterion terminate the adaptation when weight changes are judged statistically insignificant at more than 99% probability.

The first test case simulated reactive transport of tetrachloroethylene, or PCE, through a 500 m by 300 m by 20 m synthetic aquifer with relatively low permeability and dendritic high-permeability pathways. A pumping well near the opposite boundary continuously monitored PCE and its sequential dechlorination products, trichloroethylene, dichloroethylene, and vinyl chloride. The working model retained only the main preferential pathways, deliberately omitting smaller branches present in the reference system.

After Gaussian noise was added, the experiment generated 50 concentration observations, with 30 used for identification and 20 reserved for validation. Scenarios compared an unconstrained correction, one coupled to carbon conservation, and one coupled to both carbon and chlorine conservation. In that simulation, adding the carbon constraint reduced the mean carbon mass-balance error rate from 48.62% to 35.03%, while the chlorine rate declined only slightly from 30.65% to 28.27%.

Coupling both elemental balances produced a larger change: mean carbon and chlorine error rates fell to 23.63% and 16.49%, corresponding to average mass errors of 0.86 and 3.52 kg. The multiple-constraint scenario also narrowed posterior distributions and produced more stable structural-error estimates. During validation, the Nash–Sutcliffe efficiencies for PCE, dichloroethylene, and vinyl chloride rose from negative values to about 0.8, and their mean absolute errors and root mean square errors decreased by approximately 95%, 81%, and 75%, respectively.

The second test used a 5 cm-long polypropylene sand column to examine cadmium transport in the presence of phosphate. A solution containing 0.1 mM cadmium and 1 mM phosphate was injected at 0.2 mL per minute, producing 58 concentration observations; the first 32 were used for identification and the remaining 26 for validation. The physical model represented cadmium sorption but omitted cotransport between cadmium and phosphate, creating a known source of structural error.

The tested scenarios compared no physical constraint, cadmium mass conservation alone, and combined mass conservation with a bound requiring normalized cadmium concentration to remain between zero and 0.75. The upper bound was defined from the experimental breakthrough curve and used as a physical restriction in the model, rather than treated as a field observation. The single cadmium conservation constraint lowered the mean mass-balance error rate from 30.14% to 13.64% and the peak error from 57.14% to 26.07%.

Adding the concentration bound reduced those values further to 10.36% and 19.33%. The unconstrained model reproduced the overall rise and fall of the concentration curve but overestimated the peak and underestimated later observations. The mass-constrained scenario improved the fit, yet some posterior predictions still fell below zero or above the maximum concentration allowed by the experimental constraint.

With both constraints, all posterior predictions remained within the admissible range, the 95% prediction intervals became narrower, and the validation-period Nash–Sutcliffe efficiency increased from 0.4874 without constraints to 0.6241 with mass conservation and 0.9419 with both constraints. Relative to the single-constraint case, the final scenario reduced mean absolute error by 56.52% and root mean square error by 60.61%. These results suggest that appropriately selected physical constraints can act like additional observations in a data-limited Bayesian analysis.

In the synthetic case, only 50 concentration measurements were available; in the column experiment, there were 58. Constraints reduced equifinality, restricted excessive adjustment of the Gaussian process hyperparameters, and narrowed posterior distributions without requiring changes to the architecture of the physical model. They also simplified the inferred structural error, making it less variable and easier for Gaussian process regression to characterize.

The effect was not simply a matter of adding more constraints: the study emphasizes that each constraint should be mechanistically relevant, quantitatively defined, and complementary to the observations and to the other constraints. The evidence comes from one synthetic field-scale analogue and one controlled column experiment, however, and the authors caution that real contaminated sites and watershed-scale applications would require representative observations, reliable diagnosis of structural error, and constraints directly linked to neglected physical processes. The adaptive weighting procedure also adds computational cost, and adding more constraints increases the complexity of the likelihood and may affect convergence.

The method therefore presents a generalizable framework tested in two settings, not a guarantee of equivalent performance at every field site. Its central advance is to place physical limits around the full corrected prediction—physical model, structural error, and measurement error—rather than imposing constraints only on the machine-learning error model alone.

Subject of Research: Physics-constrained Bayesian correction of structural errors in groundwater contaminant transport models

Article Title: A data-driven method coupling multiple physical constraints for correcting structural errors in groundwater contaminant transport models

Article References: Tian, J., Zeng, X., Pan, Y., Wang, D., & Wu, J. (2026). A data-driven method coupling multiple physical constraints for correcting structural errors in groundwater contaminant transport models. Hydrology and Earth System Sciences, 30(20), 6459-6479. https://doi.org/10.5194/hess-30-6459-2026

Image Credits: AI Generated

DOI: 10.5194/hess-30-6459-2026

Keywords: groundwater contamination, contaminant transport, Bayesian inference, model structural error, physical constraints, Gaussian process regression, mass conservation, predictive uncertainty, tetrachloroethylene, cadmium cotransport, data-driven, method