🔢 mathematics

Data-driven linear solver selection and performance tuning for multiphysics simulations in porous media

This paper presents a data-driven, adaptive machine learning algorithm that automatically selects and tunes preconditioned linear solvers for multiphysics porous media simulations by continuously updating its policy with runtime performance data, achieving efficiency and robustness comparable to expert-guided methods without requiring prior knowledge.

Original authors: Yury Zabegaev, Inga Berre, Eirik Keilegavlen

Published 2026-06-16

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: Yury Zabegaev, Inga Berre, Eirik Keilegavlen

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to solve a massive, incredibly complex puzzle. This isn't a jigsaw puzzle with pictures, but a "physics puzzle" that simulates how water, heat, and rock move and interact deep underground. Scientists call this "multiphysics simulation."

To solve this puzzle, the computer has to break it down into millions of tiny math problems (linear systems) that need to be solved one after another. The problem is, there are thousands of different "tools" (called solvers) the computer could use to solve each tiny math problem. Some tools are fast but might break the puzzle; others are slow but very safe.

Choosing the right tool for the right moment is like trying to pick the perfect wrench from a toolbox of 10,000 different wrenches, where the shape of the bolt changes every second. If you pick the wrong one, the computer wastes hours or even days trying to force it to work, or it gives up entirely.

The Problem: Too Many Choices, Not Enough Time

Traditionally, experts would sit down and guess which tool to use based on experience. But because the underground conditions change constantly (like the temperature or pressure shifting), a tool that worked yesterday might fail today. Collecting data on every single tool for every possible scenario would take forever—longer than the simulation itself.

The Solution: A Smart, Learning Assistant

The authors of this paper created a smart assistant (an algorithm) that acts like a seasoned mechanic who learns on the job. Here is how it works, using a simple analogy:

1. The "Trial Run" (Initial Exploration)
When the simulation starts, the assistant doesn't know the terrain yet. So, it randomly picks a few tools to try out. It's like a chef tasting a few different spices to see what the dish needs. It tries a tool, sees if it works, and records how fast it was.

2. The "Two-Step Brain" (Machine Learning Pipeline)
Once it has a little bit of data, the assistant uses a two-step thinking process for every new math problem:

Step 1: The Filter (The Classifier): Imagine a bouncer at a club. This part of the brain looks at the available tools and immediately says, "No, that one is going to fail," or "No, that one is too slow." It quickly throws out the bad options so the computer doesn't waste time on them.
Step 2: The Predictor (The Regressor): From the remaining "good" tools, this part predicts exactly which one will be the fastest for this specific moment. It's like a GPS that says, "Based on the traffic right now, take Route B, not Route A."

3. Learning as You Go (Incremental Updates)
The best part is that the assistant never stops learning. Every time it solves a math problem, it updates its memory. If it picks a tool and it turns out to be slower than expected, it learns from that mistake. If a tool works perfectly, it remembers that too. It gets better and better the longer the simulation runs, refining its choices without needing a massive library of past data.

The Experiments: Testing the Assistant

The researchers tested this assistant on two very difficult "puzzles":

Heat and Water: Simulating how hot water moves through porous rock (like a sponge).
Rock, Water, and Fractures: Simulating how rock cracks and shifts when water and heat are involved, which is crucial for things like geothermal energy or oil extraction.

They ran these simulations with two different approaches:

The "Random" Approach: Picking tools blindly.
The "Expert" Approach: Using a computer that had seen every possible tool on every possible scenario beforehand (the "perfect" expert).
The "Learning" Approach: The new assistant described above.

The Results: Almost as Good as the Expert

The results were impressive:

Avoiding Disaster: The learning assistant was very good at avoiding tools that would fail completely. It kept the simulation running smoothly.
Speed: It found tools that were in the top 3% of the best possible options.
Efficiency: It didn't need to try thousands of tools to figure it out. It only needed to try a tiny handful (about 64 random attempts) at the very beginning to get started.
Comparison: Surprisingly, the "Learning" assistant performed almost exactly as well as the "Expert" who had all the data in the world. The "Expert" was only slightly faster, but the cost of gathering all that data beforehand would have been impossible in real life.

Why This Matters

In the real world, simulation engineers often don't have a team of math experts standing by to tweak the settings for every new project. This "learning assistant" allows the computer to tune itself automatically. It saves massive amounts of time and computing power, making it possible to run complex simulations that were previously too difficult or expensive to solve.

In short, the paper shows that you don't need to know everything before you start. You just need a smart system that learns quickly from its mistakes and knows which tools to avoid.

Technical Summary: Data-Driven Linear Solver Selection and Performance Tuning for Multiphysics Simulations in Porous Media

1. Problem Statement

Simulating coupled multiphysics processes in porous media (e.g., fluid flow, heat transfer, and poromechanics) requires solving large, sparse, non-symmetric linear systems at each time step or Newton iteration. These systems are typically solved using preconditioned iterative methods, such as GMRES. The efficiency and robustness of these solvers depend heavily on the choice of preconditioner and its numerical parameters.

The core challenges identified in this work are:

Combinatorial Explosion: Preconditioners for multiphysics problems are often constructed by combining sub-algorithms for individual physics components (e.g., Algebraic Multigrid (AMG) for mechanics, Constrained Pressure Residual (CPR) for flow). This results in a vast configuration space, ranging from thousands to tens of thousands of possible solver configurations (e.g., $10^4$ in the studied cases).
Context Dependence: Solver performance is highly sensitive to the simulation state (e.g., pressure, temperature, fracture contact status). A configuration optimal for one time step or one set of physical parameters may be suboptimal or fail entirely for another.
Data Scarcity and Relevance: Traditional optimization relies on large datasets collected from prior simulations. However, as simulation setups change (e.g., different permeability fields or well locations), historical data becomes less representative. Collecting a comprehensive dataset for every new scenario is computationally prohibitive, as testing all configurations can take days.
Expertise Gap: Manually tuning these solvers requires significant expert knowledge, which is not always available to simulation engineers.

The research question posed is: How can we design a selection algorithm for preconditioned linear solvers that ensures strong simulation performance while minimizing the time spent exploring suboptimal solvers, particularly when historical data is insufficient or irrelevant?

2. Methodology

The authors propose a data-driven, adaptive solver selection algorithm that operates during the simulation runtime. The approach utilizes a two-stage machine learning (ML) pipeline that continuously updates its policy based on feedback from the current simulation.

2.1 Solver Configuration Space

The algorithm treats the selection of a solver as an optimization problem over a configuration space $\mathcal{A}$ . Each configuration $a \in \mathcal{A}$ is encoded as a vector containing:

Categorical decisions: Algorithm choices (e.g., GMRES vs. Direct, CPR vs. System-AMG, specific AMG variants like BoomerAMG or GAMG).
Numerical parameters: Tunable hyperparameters (e.g., AMG strong thresholds, relaxation types, GMRES restart values).
These are encoded into binary and real-valued vectors to serve as ML inputs.

2.2 Simulation Context

To address the context-dependence of solver performance, the algorithm uses a simulation context vector $c$ . This vector compresses the full simulation state $\chi$ into a fixed number of features, such as:

Time step size ( $\Delta t$ ).
CFL number.
Temperature and pressure bounds.
Flux magnitudes (enthalpy, diffusion).
Mesh quality metrics.
Permeability estimates (for heterogeneous cases).

2.3 The Machine Learning Pipeline

The selection process employs a sequential pipeline of two ML models:

Classifier: Predicts the probability of success ( $S$ ) for a given configuration $a$ and context $c$ . Configurations predicted to fail are discarded.
Regressor: Estimates the expected reward ( $R$ ), defined as the negative logarithm of the solve time ( $-\log T$ ), for the remaining configurations.

The final selection $a_{pred}$ is the configuration that maximizes the predicted reward among those predicted to succeed:
$a_{pred} = \arg\max_{a \in \mathcal{A}} r(a, c) \quad \text{s.t.} \quad s(a, c) = \text{true}$

2.4 Learning Strategy

The algorithm operates in two phases to balance exploration and exploitation:

Initial Exploration: At the start of a simulation (or when no relevant historical data exists), the algorithm randomly selects configurations for a fixed number of linear systems ( $num_{initial} = 64$ in the experiments). This gathers initial performance data specific to the current simulation setup.
Incremental Learning: The ML models (Gradient Boosting classifiers and regressors) are updated in batches. After the initial phase, the algorithm selects the best predicted configuration for each new linear system. If the solver fails, it retries with a new configuration. The performance data (context, configuration, success/failure, time) is accumulated and used to retrain the models periodically.

3. Key Contributions

Adaptive Selection Framework: The paper presents a novel framework that collects performance data on-the-fly during the target simulation, rather than relying on pre-collected datasets from different setups. This addresses the issue of data relevance in changing simulation environments.
Handling Large Configuration Spaces: The method is designed to operate within a space of $\approx 10^4$ configurations, utilizing a classifier to filter out likely failures before applying a regressor for performance estimation.
Context-Aware Optimization: By explicitly incorporating simulation context features, the algorithm adapts to the changing physics of the problem (e.g., transitions from diffusive to advective heat transfer regimes).
Low Overhead Implementation: The authors demonstrate that the computational cost of the ML inference and retraining is negligible compared to the time required to solve the linear systems.

4. Numerical Experiments and Results

The algorithm was evaluated on two time-dependent nonlinear model problems:

Sequence A: Coupled fluid flow and heat transfer in a heterogeneous porous medium (based on SPE10 benchmarks). Configuration space: ~14,415 options.
Sequence B: Thermo-poromechanics with fractures and contact mechanics. Configuration space: ~32,000 options.

Key Findings:

Robustness and Efficiency: The solver selection algorithm successfully avoided bad configurations (those that fail to converge or take excessive time).
- In Sequence A, the success rate after initial exploration reached 88.88%, compared to 64.5% in random selection.
- In Sequence B, the success rate reached 99.92%, compared to 91.3% in random selection.
Performance Proximity to Optimal: The algorithm consistently selected configurations within the top 1.9% to 2.7% of all possible configurations in terms of solve time.
Comparison with "Expert" Baseline: The authors compared their adaptive approach against an "expert" baseline that had access to all performance data collected from previous runs (simulating a scenario with perfect prior knowledge).
- The adaptive algorithm performed comparably to the expert baseline.
- In Sequence B, the adaptive approach sometimes outperformed the expert baseline because the expert's need to retrain on a massive dataset introduced slightly higher overhead, whereas the adaptive method used a smaller, more relevant dataset.
Overhead: The machine learning overhead was negligible (approx. 0.02–0.04 seconds per linear system), while the linear solve times were on the order of seconds to minutes.

5. Significance and Claims

The paper claims that the proposed approach effectively addresses the challenge of solver selection and tuning in multiphysics simulations where:

The configuration space is too large for exhaustive search.
Historical performance data is insufficient or irrelevant due to changes in simulation setup.
Expert knowledge for manual tuning is unavailable.

The authors emphasize that their method provides a practical solution for simulation engineers, enabling robust and efficient simulations with minimal manual intervention. The results suggest that "inexpensive tuning on the fly" can yield performance nearly identical to that of an expert with full prior knowledge, without the prohibitive cost of collecting massive historical datasets. The approach is particularly valuable for industrial-scale simulations where reliability and automation are critical.