Accelerating 3D Magnetotelluric Forward Modelling with Domain Decomposition and Order-Reduction Methods

Tao, Luis; Muixí, Alba; Zlotnik, Sergio; Zyserman, Fabio Ivan; Afonso, Juan Carlos; Diez, Pedro

doi:10.5194/egusphere-2026-700

Preprints

https://doi.org/10.5194/egusphere-2026-700

Preprints

20 Apr 2026

| 20 Apr 2026

Status: this preprint is open for discussion and under review for Geoscientific Model Development (GMD).

Accelerating 3D Magnetotelluric Forward Modelling with Domain Decomposition and Order-Reduction Methods

Luis Tao, Alba Muixí, Sergio Zlotnik, Fabio Ivan Zyserman, Juan Carlos Afonso, and Pedro Diez

Abstract. Three-dimensional (3D) magnetotelluric (MT) forward modelling is computationally demanding, limiting its use in global uncertainty quantification and sampling-based probabilistic inversion. Here, we introduce a novel forward-modelling framework that combines an iterative domain decomposition (DD) formulation with proper orthogonal decomposition (POD) reduced-order modelling to enable scalable and efficient 3D MT simulations. The DD component partitions the computational domain into subdomains, avoiding the factorization of a single global system, accelerating simulations by over 60 % compared to global solvers, and alleviating memory bottlenecks in large problems. The POD component leverages the local DD solutions to construct a reduced-order version of the problem that can deliver accurate and fast solutions to the 3D forward problem during subsequent evaluations. Using the DTM1 benchmark and a real-world conductivity model, we quantify runtime, memory, and accuracy in terms of MT quantities of interest (apparent resistivity and phase). DD–POD achieves speed-ups exceeding 90 % relative to full-order solvers and up to 70 % relative to existing ROM techniques, while maintaining acceptable accuracy. These results suggest that DD–POD can make higher-resolution 3D MT forward modelling practical within sampling-based workflows by substantially reducing both runtime and memory demands.

Received: 05 Feb 2026 – Discussion started: 20 Apr 2026

Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors. Views expressed in the text are those of the authors and do not necessarily reflect the views of the publisher.

Download & links

Luis Tao, Alba Muixí, Sergio Zlotnik, Fabio Ivan Zyserman, Juan Carlos Afonso, and Pedro Diez

Status: open (until 08 Jul 2026)

Post a comment Subscribe to comment alert

RC1:
'Comment on egusphere-2026-700', Rob L. Evans, 28 May 2026 reply
This paper presents a new approach to modeling three-dimensional electrical resistivity models with the application to the magnetotelluric method.
By dividing the model into subdomains, and connecting fields across the horizontal boundaries of these layered domains, the code is able to parallelize the computation working on the smaller domains individually.
The paper demonstrates how this approach achieves reliable results in much faster computation times, something that is essential for the kinds of Bayesian approaches that develop an ensemble of models that satisfy a data set, important to providing some estimates of parameter uncertainty.
The paper is well written and convincing and should be published more or less as is. Rather than specific critiques, I do have a couple of questions that would be worth diving into a bit more, should the authors so choose.
Topography: One of the issues in modeling data over complex topography – I’ll admit a bias to seafloor measurements, but similar issues apply to volcanic geothermal settings – is the need to finely mesh near the surface, potentially making the total model larger than perhaps it needs to if those elements extend to depth in a finite difference model. The current set up shown in the paper has the same (as far as I can tell) Nx and Ny mesh elements through all the subdomains of the model. Is it possible to have variable Nx and Ny, such that there is fine meshing in a surface (topographic) subdomain, and fewer elements at depth? Would the boundary field matching across subdomains work as well in this kind of approach?

As it stands, the errors in the data values, at least for the example shown, seem to be around 5% - that is on the order of typical error floors set for field data, so I’d be concerned trusting an inversion where the model errors were of the same order. I think this merits a bit more discussion in the paper, if the goal here is to move towards an ensemble-type inversion.

One of the impediments my group has found to the utilization of 3D MT codes are the tools (or lack thereof) for generating meshes, editing models, and inputting data into inversions. Could the authors speak a little as to how they generate the meshes and what tools, if any, are available for this?

Reply
Citation: https://doi.org/10.5194/egusphere-2026-700-RC1
- AC1:
  'Reply on RC1', Luis Tao, 02 Jun 2026 reply
  We appreciate the time taken to read the MS and for the constructive feedback. We also agree that diving into this questions would benefit the MS.
  Aside from the changes made in the MS (in the supplement) we also post here the answers for the questions:
  1) Thank you for this question. In our current implementation, which is based on structured meshes, the horizontal resolution remains constant with depth. We can, however, control the height of each layer. This allows us to increase the vertical resolution near the surface to better account for topography, while using larger vertical elements at greater depths. While we acknowledge this does not entirely resolve the issue you raised, it is important to note that the Domain Decomposition formulation we employ is highly general. It can accommodate the use of subdomains with varying Nx and Ny dimensions, including subdomains with less elements at depth, which could potentially lead to an even more efficient solver.
  
  We have updated the manuscript in Section 4 (Implementation Aspects), Subsection 4.1 (Domain Discretization), to include a comment noting that the formulation allows for more flexible partitions.
  2) This is a very fair point. The error is essentially controlled by four factors: the convergence threshold in the domain decomposition (δ), the penalty factor (β) applied to the electric field mismatch at the interfaces, the snapshots generated to create the basis, and the number of modes retained in that basis (nPOD). As shown in the manuscript, the error stabilizes at around 5% with increasing nPOD. If we use a tighter δ and a larger β penalty, we believe this approach will further decrease the error, though it introduces a slight trade-off in execution time. Additionally, the snapshots used to create the basis were generated via a grid search without optimizing the information they provided. We could likely improve this by employing more adaptive basis-enrichment strategies, a path we might have have to explore as we progress toward the inversion phase. Therefore, for the scope of this methodological paper, we considered our current error margin acceptable. For future inversions, we can enforce stricter hyperparameters to further decrease the error.
  
  This discussion can be found in the manuscript in Section 6 (Results)-Subsection 6.2 (Real-world model results)-SubsubSection 6.2.2 "Validation and accuracy", and in Subsection 6.3 (Implications for probabilistic inversion).
  3) Regarding mesh generation and model editing, our solver currently operates with parallelepipedal structured meshes. Consequently, the inputs for the mesh consist primarily of the dimensions of each element, and conductivity is assigned element-by-element, making it a largely manual process at this stage. Because the primary focus of this work is on advancing the numerical methodology rather than optimizing meshing capabilities, we have maintained this straightforward approach.
  
  We have added a comment to the manuscript in Section 4 (Implementation Aspects), Subsection 4.1 (Domain Discretization), to clarify this.
  
  Reply
  
  Citation: https://doi.org/10.5194/egusphere-2026-700-AC1

Luis Tao, Alba Muixí, Sergio Zlotnik, Fabio Ivan Zyserman, Juan Carlos Afonso, and Pedro Diez

Data sets

DTM1 input files for 20×20×18 (includes snapshots), 40×40×36, and 80×80×72 mesh Luis Tao https://zenodo.org/records/19567027

Model code and software

Code for (DD) and DD-POD + Code Guide Luis Tao and Fabio Zyserman https://zenodo.org/records/19567027

Luis Tao, Alba Muixí, Sergio Zlotnik, Fabio Ivan Zyserman, Juan Carlos Afonso, and Pedro Diez

Viewed

Total article views: 403 (including HTML, PDF, and XML)

HTML	PDF	XML	Total	BibTeX	EndNote
261	118	24	403	18	17

HTML: 261
PDF: 118
XML: 24
Total: 403
BibTeX: 18
EndNote: 17

Views and downloads (calculated since 20 Apr 2026)

Month	HTML	PDF	XML	Total
Apr 2026	115	45	14	174
May 2026	131	63	4	198
Jun 2026	15	10	6	31

Cumulative views and downloads (calculated since 20 Apr 2026)

Month	HTML	PDF	XML	Total
Apr 2026	115	45	14	174
May 2026	131	63	4	198
Jun 2026	15	10	6	31

Viewed (geographical distribution)

Total article views: 399 (including HTML, PDF, and XML) Thereof 399 with geography defined and 0 with unknown origin.

Country	#	Views	%

Latest update: 08 Jun 2026

Short summary

We present a new approach for performing 3D magnetotelluric forward simulations more efficiently. Conventional methods become increasingly demanding as model resolution increases. Our approach combines numerical techniques that reduce problem size and computational cost. Tests on benchmark examples and a real-world case demonstrate speed-ups of over 90% with acceptable loss of accuracy, enabling high-resolution simulations within practical time frames.


Total:	0
HTML:	0
PDF:	0
XML:	0