Autoregressive Physics-Informed Neural Operator Surrogates and Joint 4-D Seismic-Production Bayesian History Matching
Abstract
We present a surrogate-accelerated Bayesian history-matching framework, demonstrated on a realistic synthetic model (N=100,188 cells). The forward map 𝒢 is replaced by a Fourier neural operator 𝒢θ trained in a closed loop, x̂n+1=𝒢θ(u,x̂n), with the discrete black-oil residual R(xn+1,xn) as a physics-informed neural operator (PINO) penalty, truncated backpropagation and noise injection. We assimilate production data and 4-D seismic using adaptive regularised ensemble Kalman inversion (aREKI), and synthesise ΔVp inside the rollout via a closed-form Gassmann transform of the surrogate, under a two-stage schedule (aREKI on production, then an ensemble smoother with multiple data assimilation, ES-MDA, on both). To temper the ill-posedness, we invert (log K,ϕ) not in ℝ2N, where Ne≪ 2N leaves a rank-deficient covariance and a null space the data never inform, but on a chart u=Ψ(m̃), m̃∈ℝd, d≪ 2N, given by a Parseval-exact truncated cosine transform or a variational auto-encoder decoder. The ensemble then spans its own search space; the directions the data cannot constrain are cut to those the chart can express, and every member decodes to a geologically plausible field. For the rollout we prove an exposure bound W1≤∑j L n−1−j𝔼ej, a physics-penalty contraction |Dx𝒢θ|≤α<1 under a small residual defect, and an 𝒪(ρK) truncation bias. A Loewner-order inequality ℐp⪯ℐp+ℐs explains why joint assimilation is no less informative. On held-out data, 𝒢θ attains R2>0.999 on pressure and saturation and rolls out in $0.579$ s, a ∼2500× speed-up, bringing an Ne=190-member, Na=10-iteration joint inversion within one GPU-hour.
Related
- Bayesian spatiotemporal conditional autoregressive model for local temporal variations
- Operator Surrogates
- Bayesian Generalized Network Autoregressive Model with Structured Shrinkage and Persistence Priors
- Structure-Preserving Dimension Reduction and Surrogates for Partially Linear Bayesian Inverse Problems
- Pareto-Optimal Hypersonic Vehicle Design with Canonical-Coordinate Neural Operator Surrogates Under Ablation-Driven Shape Change and Incidence Variation
BibTeX
@misc{etienam2026autoregressive,
title = {{Autoregressive Physics-Informed Neural Operator Surrogates and Joint 4-D Seismic-Production Bayesian History Matching}},
author = {Clement Etienam and David Blom and Hari Pattathil and Tsubasa Onishi and Harpreet Sethi and Lalson Vincent and Suraj Pawar and Kishan Ramesh Kumar and Pandu Devarakota and Faruk O. Alpak},
howpublished = {arXiv},
year = {2026},
month = oct,
doi = {10.48550/arxiv.2610.31677},
eprint = {2610.31677},
archivePrefix = {arXiv},
url = {https://arxiv.org/abs/2610.31677},
note = {CC BY 4.0},
}
RIS
TY - UNPB
TI - Autoregressive Physics-Informed Neural Operator Surrogates and Joint 4-D Seismic-Production Bayesian History Matching
AU - Clement Etienam
AU - David Blom
AU - Hari Pattathil
AU - Tsubasa Onishi
AU - Harpreet Sethi
AU - Lalson Vincent
AU - Suraj Pawar
AU - Kishan Ramesh Kumar
AU - Pandu Devarakota
AU - Faruk O. Alpak
PY - 2026
DA - 2026/10/01
PB - arXiv
DO - 10.48550/arxiv.2610.31677
UR - https://arxiv.org/abs/2610.31677
AB - We present a surrogate-accelerated Bayesian history-matching framework, demonstrated on a realistic synthetic model (N=100,188 cells). The forward map 𝒢 is replaced by a Fourier neural operator 𝒢_θ trained in a closed loop, x̂_n+1=𝒢_θ(u,x̂_n), with the discrete black-oil residual R(x_n+1,x_n) as a physics-informed neural operator (PINO) penalty, truncated backpropagation and noise injection. We assimilate production data and 4-D seismic using adaptive regularised ensemble Kalman inversion (aREKI), and synthesise ΔV_p inside the rollout via a closed-form Gassmann transform of the surrogate, under a two-stage schedule (aREKI on production, then an ensemble smoother with multiple data assimilation, ES-MDA, on both). To temper the ill-posedness, we invert (log K,ϕ) not in ℝ^2N, where N_e≪ 2N leaves a rank-deficient covariance and a null space the data never inform, but on a chart u=Ψ(m̃), m̃∈ℝ^d, d≪ 2N, given by a Parseval-exact truncated cosine transform or a variational auto-encoder decoder. The ensemble then spans its own search space; the directions the data cannot constrain are cut to those the chart can express, and every member decodes to a geologically plausible field. For the rollout we prove an exposure bound W_1≤∑_j L^ n−1−j𝔼e_j, a physics-penalty contraction |D_x𝒢_θ|≤α<1 under a small residual defect, and an 𝒪(ρ^K) truncation bias. A Loewner-order inequality ℐ_p⪯ℐ_p+ℐ_s explains why joint assimilation is no less informative. On held-out data, 𝒢_θ attains R^2>0.999 on pressure and saturation and rolls out in $0.579$ s, a ∼2500× speed-up, bringing an N_e=190-member, N_a=10-iteration joint inversion within one GPU-hour.
ER -