Physics-constrained symbolic regression for discovering closed-form equations of multimodal water retention curves from experimental data

This paper introduces a physics-constrained machine learning framework using genetic programming to automatically discover interpretable, closed-form equations for multimodal water retention curves directly from experimental data, overcoming the limitations of traditional superposition-based models.

The Gist

Imagine you are trying to predict how much water a specific type of soil can hold as it dries out. This isn't just a simple "dry vs. wet" situation; soil is like a sponge made of different-sized holes (pores). Some holes are tiny and hold water tightly, while others are large and let water drain away quickly. When soil has a mix of these different hole sizes (a "multimodal" structure), the math to describe it gets incredibly messy.

Traditionally, scientists have tried to solve this in two ways:

1. The "Old School" Math: They use pre-made formulas (like the Van Genuchten model). It's like trying to fit a square peg into a round hole. It works okay for simple soils, but for complex, mixed soils, the formula just doesn't fit right.
2. The "Black Box" AI: They use powerful neural networks (deep learning). These are like super-smart students who can memorize the answer to every test question perfectly. But, they can't explain how they got the answer. If you ask them, "Why did you predict the soil is 110% full of water?" (which is physically impossible), they might just say, "Because the math said so." Engineers don't trust these "black boxes" because they can't be used in real-world safety calculations.

The New Solution: The "Physics-Constrained Detective"

This paper introduces a new method called Physics-Constrained Symbolic Regression (PCSR). Think of this not as a student memorizing answers, but as a detective trying to solve a mystery using a very specific set of rules.

Here is how it works, broken down with simple analogies:

1. The Genetic Programming (The "Evolutionary Chef")

Imagine a kitchen where you want to invent a new recipe (a mathematical equation) that tastes exactly like a specific dish (the experimental data).

• Instead of one chef, you have a whole population of chefs.
• They start with random ingredients (numbers, variables like "suction," and operations like "plus" or "multiply").
• They cook up thousands of random recipes.
• Selection: The ones that taste closest to the target dish are kept.
• Mutation & Crossover: The best chefs swap ingredients or tweak their recipes slightly to see if they can make it even better.
• Over many generations, the recipes evolve into a perfect mathematical equation. This is called Symbolic Regression.

2. The Problem: The "Wild Chef"

The problem with this evolutionary approach is that the chefs can get crazy.

• They might invent a recipe that says, "Add 500 cups of water to a cup of soil."
• Mathematically, it fits the data points you gave them perfectly.
• But physically, it's nonsense. Water can't be more than 100% of the soil, and it shouldn't magically appear out of nowhere.
• Without rules, the AI finds "solutions" that are mathematically correct but physically impossible.

3. The Fix: The "Physics Police"

This is where the paper's innovation shines. The authors act as Physics Police who stand over the chefs' shoulders during the cooking process. They enforce three strict rules (Constraints) that the recipe must follow:

• The "No Magic Water" Rule (Monotonicity): As the soil gets drier (suction increases), the water content must go down. It can never go up. If the recipe suggests water increases as the soil dries, the Police slap it with a fine (a penalty in the math).
• The "Full and Empty" Rule (Limiting):
  • When the soil is soaking wet, it must be 100% full.
  • When it is bone dry, it must stop at a specific "residual" amount (a tiny bit of water stuck in the tiny pores).
  • The recipe cannot predict 150% full or -10% full.
• The "Shape" Rule (Modes): If you know the soil has two distinct types of holes (big and small), the curve describing the water loss should have a specific "wavy" shape (two bumps). The Police ensure the recipe doesn't accidentally invent a third or fourth bump that doesn't exist.

The Result: The "Perfect Hybrid"

By combining the creative evolution of the chefs with the strict rules of the Physics Police, the system discovers a mathematical equation that:

1. Fits the data perfectly (it matches the experimental measurements).
2. Is physically possible (it obeys the laws of nature).
3. Is transparent (it's a clear, readable equation, not a black box).

Why This Matters

Imagine you are an engineer designing a dam or a building foundation. You need to know how the soil will behave when it rains or dries out.

• If you use the Black Box AI, you might get a prediction that looks good on a computer but fails in real life because it violated a law of physics.
• If you use the Old School Math, you might get a safe prediction, but it's too simple to handle complex soils.
• With PCSR, you get a custom-made, mathematically precise equation that is guaranteed to be physically safe. It's like giving the engineer a custom-tailored suit that fits perfectly and is made of fire-proof material.

In a nutshell: This paper teaches computers how to "learn" the laws of physics while they are learning from data, ensuring that the new formulas they discover are not just smart, but also sensible and safe for the real world.

Technical Summary

Here is a detailed technical summary of the paper "Physics-constrained symbolic regression for discovering closed-form equations of multimodal water retention curves from experimental data."

1. Problem Statement

Modeling the unsaturated behavior of porous materials with multimodal pore size distributions (e.g., soils with distinct sand and clay fractions) is a significant challenge in geotechnical engineering.

• Limitations of Standard Models: Traditional semi-empirical models (e.g., Van Genuchten, Brooks-Corey) assume a unimodal pore size distribution. They fail to accurately capture the complex, multi-scale retention characteristics of multimodal materials.
• Limitations of Current Workarounds: A common approach is to superpose multiple unimodal functions. However, this requires separate parameter identification for each mode, reducing interpretability and generalizability, especially when experimental data is sparse.
• Limitations of Pure Data-Driven Methods: While Deep Neural Networks (DNNs) can learn complex non-linear relationships, they act as "black boxes," lacking interpretability and physical consistency. They often fail to satisfy fundamental thermodynamic constraints (e.g., predicting saturation $>1$ or non-monotonic behavior).
• Limitations of Standard Symbolic Regression (SR): While SR offers interpretability by discovering closed-form equations, it suffers from scalability issues, sensitivity to noise, and a lack of physical constraints, often leading to overfitting or physically impossible solutions.

2. Methodology: Physics-Constrained Symbolic Regression (PCSR)

The authors propose a Physics-Constrained Symbolic Regression (PCSR) framework designed to automatically discover closed-form mathematical expressions for water retention curves directly from experimental data.

A. Problem Formulation

The task is framed as a multi-objective optimization problem. The goal is to find a function $\hat{S}_w(s)$ mapping suction pressure $s$ to water saturation ratio $S_w$.

• Data Preprocessing: Raw data $(s, S_w)$ is mapped to a normalized space $(s^*, S^*_w) \in [0, 1]$ using reference points (minimum suction $s_{min}$, residual suction $s_{res}$, max saturation $S_{w,max}$, residual saturation $S_{w,res}$). This stabilizes the learning process.
• Objective Function: The total loss function $L$ consists of three components:
  $$L = L_{data} + L_{phys} + L_{mode}$$
  1. $L_{data}$ (Data Loss): Mean Squared Error (MSE) between predictions and normalized experimental data.
  2. $L_{phys}$ (Physics Loss): Penalizes violations of physical constraints at specific collocation points.
  3. $L_{mode}$ (Mode Loss): Penalizes deviations from a pre-specified number of modes ($N_{mode}$) in the pore size distribution.

B. Physical Constraints

The framework enforces three critical thermodynamic constraints via the $L_{phys}$ term:

1. Monotonicity: Saturation must decrease as suction increases ($dS_w/ds \leq 0$).
2. Limiting Conditions:
   • At wet end ($s \leq s_{min}$): Saturation is constant at maximum ($S_w = S_{w,max}$) and the derivative is zero.
   • At dry end ($s \geq s_{res}$): Saturation approaches residual ($S_w = S_{w,res}$) and the derivative is zero.
3. Boundedness: Saturation must remain within $[0, 1]$.

C. Algorithm and Implementation

• Genetic Programming (GP): Mathematical expressions are represented as binary trees and evolved using selection, mutation, and crossover operations.
• Constraint Enforcement: Unlike differentiable neural networks, the constraints in PCSR are enforced via a custom penalization strategy within the loss function.
  • Monotonicity is checked using the ReLU function on the derivative at collocation points.
  • Mode counting is performed by analyzing sign changes in the second derivative (inflection points).
• Software: Implemented using the Julia package SymbolicRegression.jl.

3. Key Contributions

1. Novel Framework: Introduction of a PCSR framework specifically tailored for multimodal water retention curves, bridging the gap between data-driven discovery and physical laws.
2. Mode Constraint: A unique contribution is the inclusion of a mode-guiding loss term ($L_{mode}$). This allows the user to specify the number of pore-size modes (e.g., $N_{mode}=2$ for bimodal soils) as a structural prior, preventing the model from overfitting noise into spurious modes.
3. Closed-Form Interpretability: Unlike DNNs, the output is an explicit, analytical equation that can be directly integrated into existing hydro-mechanical simulation codes.
4. Robustness to Noise: The framework demonstrates superior performance in sparse and noisy data scenarios by suppressing non-physical fluctuations that standard SR or vanilla machine learning models would learn.
5. Open Source: The full implementation and datasets are made publicly available to ensure reproducibility.

4. Results

The framework was evaluated on both unimodal and multimodal datasets (experimental and synthetic).

• Unimodal Cases:

  • Compared against the Van Genuchten model and vanilla SR.
  • Vanilla SR achieved low MSE but exhibited overfitting and spurious fluctuations.
  • PCSR without mode constraint satisfied physics but sometimes failed to capture the correct unimodal shape (exhibiting extra inflection points).
  • PCSR with $L_{mode}$ successfully discovered concise, physically consistent equations with exactly one mode ($N_{mode}=1$), effectively preventing overfitting.

• Multimodal Cases ($N_{mode} = 2, 3, 4$):

  • Tested on experimental data (bimodal soils) and synthetic data (3 and 4 modes).
  • Vanilla SR and PCSR without $L_{mode}$ failed to match the target number of modes, often discovering 3–12 modes instead of the target 2, and exhibited noisy, non-physical curves.
  • PCSR with $L_{mode}$ consistently discovered equations with the exact target number of modes ($\hat{N}_{mode} = N_{mode}$). The resulting curves were smooth, physically consistent, and accurately captured complex multi-scale retention behaviors.

• Comparative Analysis: The study showed that while standard SR minimizes error, it often violates physical laws (e.g., negative derivatives, saturation > 1). The PCSR framework successfully suppressed these violations, even in the presence of Gaussian noise.

5. Significance and Future Outlook

• Engineering Impact: The ability to generate interpretable, closed-form equations for complex multimodal soils allows for seamless integration into large-scale hydro-mechanical simulation codes, which is currently difficult with black-box neural networks.
• Data Efficiency: The method is particularly valuable for data-sparse scenarios where traditional parameter calibration is unreliable. The physics constraints act as a regularizer, guiding the search toward physically plausible solutions.
• Future Work: The authors plan to extend the framework to account for hysteresis (wetting vs. drying paths) and explore applications in broader geotechnical and environmental problems. They also aim to investigate the physical significance of individual parameters within the discovered symbolic structures.

In summary, this paper presents a robust, physics-informed machine learning approach that overcomes the "interpretability-accuracy trade-off," delivering transparent, mathematically rigorous models for complex porous media behavior.