# Accelerating Phase Field Simulations Through a Hybrid Adaptive Fourier
Neural Operator with U-Net Backbone

Christophe Bonneville Sandia National Laboratories, Livermore, CA 94550;
<EMAIL_ADDRESS>Nathan Bieberdorf Materials Sciences
Division, Lawrence Berkeley National Laboratory and Department of Materials
Science and Engineering, University of California, Berkeley, CA 94720, USA;
<EMAIL_ADDRESS>Arun Hegde11footnotemark: 1 Mark
Asta22footnotemark: 2 Habib N. Najm11footnotemark: 1 Laurent Capolungo Los
Alamos National Laboratory, Los Alamos, NM 87544<EMAIL_ADDRESS>Cosmin
Safta11footnotemark: 1

###### Abstract

Prolonged contact between a corrosive liquid and metal alloys can cause
progressive dealloying. For one such process as liquid-metal dealloying (LMD),
phase field models have been developed to understand the mechanisms leading to
complex morphologies. However, the LMD governing equations in these models
often involve coupled non-linear partial differential equations (PDE), which
are challenging to solve numerically. In particular, numerical stiffness in
the PDEs requires an extremely refined time step size (on the order of
$10^{-12}s$ or smaller). This computational bottleneck is especially
problematic when running LMD simulation until a late time horizon is required.
This motivates the development of surrogate models capable of leaping forward
in time, by skipping several consecutive time steps at-once. In this paper, we
propose U-Shaped Adaptive Fourier Neural Operators (U-AFNO), a machine
learning (ML) based model inspired by recent advances in neural operator
learning. U-AFNO employs U-Nets for extracting and reconstructing local
features within the physical fields, and passes the latent space through a
vision transformer (ViT) implemented in the Fourier space (AFNO). We use
U-AFNOs to learn the dynamics mapping the field at a current time step into a
later time step. We also identify global quantities of interest (QoI)
describing the corrosion process (e.g. the deformation of the liquid-metal
interface, lost metal, etc.) and show that our proposed U-AFNO model is able
to accurately predict the field dynamics, in-spite of the chaotic nature of
LMD. In particular, our model reproduces the key micro-structure statistics
and QoIs with a level of accuracy on-par with the high-fidelity numerical
solver. Finally, we also investigate the opportunity of using hybrid
simulations, in which we alternate forward leap in time using the U-AFNO with
high-fidelity time stepping. We demonstrate that while advantageous for some
surrogate model design choices, our proposed U-AFNO model in fully auto-
regressive settings consistently outperforms hybrid schemes.

_K_ eywords Phase Field $\cdot$ Liquid-Metal-Dealloying $\cdot$ Corrosion
$\cdot$ Surrogate Model $\cdot$ Vision Transformers $\cdot$ U-Nets $\cdot$
Neural Operators

##  1 Introduction

Phase field modeling is a powerful tool for simulating microstructure and
morphological evolution of materials systems subjected to various driving
forces (e.g. chemical, mechanical, electrical, etc.). The phase field method
relies on a “diffuse-interface” approach to tracking phases, interfaces, as
well as conserved quantities (e.g. mole-fraction of an element). This allows
for temporal evolution of the system to be defined everywhere by a set of
partial differential equations (PDEs) representing well-known continuum and
interfacial kinetic processes [1, 2, 3]. While this avoids the numerical
complexity associated with explicitly tracking and updating a sharp interface
in addition to providing a route to describing interface reactions, phase
field modeling can still require a high-resolution spatial discretization to
accurately resolve nano-scale interfacial phenomena. This can render large
length- and time-scale phase field simulations prohibitively computationally-
expensive, thereby limiting the use of phase field modeling to coarse grain
microstructure evolution.

Several approaches have been developed to circumvent the computational expense
of phase field simulations. One approach is to use adaptive meshing to set the
spatial discretization to be very fine near the interface and much coarser in
the bulk phases where longer-range transport processes take place [4].
However, this is not guaranteed to be advantageous in situations where the
interface morphology is rapidly changing, since frequent re-initialization of
the adaptive mesh would be required. Another approach is to artificially
broaden interfaces – from the nanometer-scale to, for example, the micron-
scale – and thereby resolve the system with a relatively coarse mesh at
interfaces as well as within bulk phases. In this approach, transport
equations and effective kinetic parameters are rederived to ensure the phase
field simulations evolve consistently with a realistic thin-interface system
[5, 6, 7]. But any approach using a micron-scale interface is inherently
limited in its ability to resolve morphological features that initiate at the
nano-scale. And so, while these approaches are valuable for a wide variety of
phase field simulations, they are not necessarily suitable for situations with
fine-scale, dynamic morphologies.

One such example is the dealloying corrosion of metals, where the infiltration
of a corrosive agent into the alloy can rapidly lead to morphologically-
complex metal structures [8, 9, 10, 11, 12, 13]. Understanding corrosion-
induced dealloying is particularly critical in energetic applications, such as
next-generation molten-salt nuclear reactor design [14, 15]. In these
aforementioned metal structures, preferential dissolution of the alloy
precedes an instability that initiates nano-scale, compositionally-
heterogenous features along the interface [9]. In a coupled manner, these
features grow while the corrosive agent infiltrates as individual channels
into the alloy. Reemergence of the nano-scale interfacial instability can
cause these channels to bifurcate and advance tortuously into the alloy,
increasing the alloy’s topological complexity [16]. Simultaneously, capillary
forces act to coarsen the dealloyed structure and reduce local curvatures and
overall topological complexity [17]. Competition among these processes sets up
a dynamic morphology that underlies the rate of corrosion into the base metal.
Phase field simulations for dealloying corrosion, specifically for liquid
metal dealloying, have successfully reproduced the experimentally-observed
morphologies and uncovered the key mechanisms by which these structures evolve
[16, 18, 19, 20]. To simulate representative domains and time-scales with high
resolution, these phase field simulations have generally relied on using high
performance computing resources and advanced time-integration schemes [16,
20]. However, the length- and time-scales of such simulations are several
orders of magnitude below experiments, so that validation requires
extrapolation using scaling laws.

The computational cost of numerical simulations such as the ones involved in
liquid-metal dealloying is a major bottleneck. This is particularly true when
a large number of forward simulations are required, like in uncertainty
quantification [21, 22], inverse problems [21, 23, 24, 25], and design
optimization [26, 27, 28] applications. This computational challenge has
historically motivated the development of surrogate models that are faster
than the high-fidelity simulations, albeit being less accurate. This drop in
accuracy can be an acceptable compromise when the time gain is significant. A
particular class of surrogate models for PDEs, such as reduced-order-models
(ROMs) aim to decrease computational cost by reducing the dimensionality of
the problem. This can be done by projecting snapshots of high fidelity
simulations into a lower-dimensional latent space. Linear projection methods,
such as proper-orthogonal decomposition (POD), are widely used for building
accurate ROMs [29, 30, 31, 32, 33], and have been successfully applied to a
variety of problems, e.g. in fluid mechanics [34, 35, 36, 37, 38, 39],
structural dynamics [40, 41, 42], and control systems [43, 44, 45]. More
recently, employing machine learning (ML) approaches to compress high fidelity
data snapshots has gained significant interest. Indeed, non-linear mappings,
such as neural networks, can yield better performance especially in advection-
dominated problems [46, 47, 48, 49].

With recent advances in deep learning and computer vision, a second class of
surrogate models for PDEs has emerged, namely neural operators [50, 51, 52].
Neural operators learn mappings between function spaces through a resolution-
invariant tensor-to-tensor regression. Popular neural operator architectures
include Fourier neural operator (FNO) [51], DeepONets [52], and more recently
Laplace Neural Operator (LNO) [53] and Convolutional Neural Operator (CNO)
[54] constructions. In time-dependent problems, neural operators can be
employed to predict the physical state at a future time step, given the
current state. Often, this strategy is used auto-regressively, i.e. the
operator output becomes its own input, allowing for predicting consecutive
time steps. This strategy (sometimes referred to as auto-regressive roll-out)
is analogous to conventional explicit time integrators. Auto-regressive neural
operators roll-outs have the convenient advantage of being naturally well
suited for non-parametric initial conditions (e.g. random initial field),
unlike more traditional ROM approaches. As such, neural operators have been
applied to a variety of problems, including fluid dynamics [51, 55, 56],
electromagnetics and electro-convection [57, 58], fracture mechanics [59],
weather prediction [60, 61], phase field modeling [62, 63] and plasma physics
[64].

While neural operators have shown great accuracy in the aforementioned papers,
recent work suggests that this accuracy can be increased further by combining
key neural operator building blocks with U-Nets [65] and Vision transformers
(ViTs) [66]. Recent papers proposed to employ the encoder-decoder architecture
of a U-Net as a feature extraction mechanism, in combination with a neural
operator bottleneck [63, 67, 68]. Building on this work, ViTO [69] and DiTTO
[70] proposed a neural operator architecture based on ViTs to capture
interdependence relations within the input field compressed features, where
the compression is obtain by a U-Net. In this paper, we propose to build on
these ideas by strategically placing an Adaptive Fourier Neural Operator
(AFNO) [71] at the intermediate bottleneck of a U-Net. An AFNO is a special
implementation of a ViT, where the attention mechanism specific to transformer
architectures is computed in the Fourier space. This idea, initially
introduced with FNO [51] allows for significant computational efficiency gains
during training. We call our proposed model U-AFNO, and show that it is able
to capture the phase field dynamics occurring in LMD with excellent accuracy,
in-spite of the highly chaotic underlying physics.

LMD phase fields are typically used to to extract local and global quantities
of interest (QoIs). Such QoIs, describing, for example, the topology of
corroded/dealloyed metal interfaces or the chemical composition within the
alloy may be crucial for decision making, optimization, or uncertainty
quantification purposes. In the present study, we employ several physically
relevant QoIs to describe the dealloying process at both local and global
scales. Furthermore, we show that our proposed U-AFNO model is accurate, not
only according to state-of-the-art error metrics [51, 52, 50, 72, 73], but
also in correctly reproducing the aforementioned QoIs, unlike other commonly
used neural operator models. Lastly, we also investigate the performance of
blending our proposed surrogate model with high-fidelity simulations. In a
recent ML-based phase field surrogate study [73], Oommen et al. proposed to
augment a U-Net-based model with alternate high-fidelity solver time stepping.
After each leap in time through a U-Net forward pass, several of the following
time steps are computed with the high-fidelity solver. This hybrid approach
still partially relies on the high-fidelity solver, and is consequently
significantly slower than fully auto-regressive roll-outs, but allows for
greater (and tunable) accuracy. In this paper, we investigate the performance
of U-AFNOs with and without hybrid high-fidelity time stepping. In particular,
we show that even without hybrid time stepping augmentation, our proposed
U-AFNO model can outperform other augmented models, such as U-Nets.

The contributions of this paper are therefore two-fold: (1) a novel surrogate
model for LMD phase field simulation to accurately reproduce the micro-
structure dynamics and the QoIs describing dealloying processes and (2) a
thorough performance analysis of hybrid ML surrogates blended with high-
fidelity solvers. Note that we employ the term high-fidelity to refer to
simulations obtained from traditional PDE solvers (sometimes referred to as
full-order-models in the literature), as opposed to surrogate models. In
section 2, we first discuss the mathematical model for LMD processes and
provide brief details of our high-fidelity solver implementation. Then, in
section 3, we introduce our proposed U-AFNO surrogate model, and, in section
4, we discuss the various local and global QoIs employed in the paper.
Finally, in section 5, we present and discuss the performance of our proposed
model and summarize this work in section 6.

##  2 High-Fidelity Phase Field Model for Liquid-Metal Dealloying

Consider a model dealloying system, where a binary alloy composed of two
species, denoted as ‘A’ and ‘B’, is in contact with a liquid dealloying agent
of pure species ‘C’ (e.g. molten salt). The model is parameterized to
represent a liquid metal dealloying system, following [16, 20]. The solid and
liquid phases are tracked using a non-conserved phase field variable $\phi$,
where $\phi=1$ in the solid phase, $\phi=0$ in the liquid phase, and
$0<\phi<1$ represents their diffuse interface. Conserved phase field variables
are used to track each species, with $c_{\text{A}}$, $c_{\text{B}}$, and
$c_{\text{C}}$ representing the mole-fraction of species A, B, and C,
respectively. Everywhere in the system,
$c_{\text{A}}+c_{\text{B}}+c_{\text{C}}=1$.

While the system in principle exhibits dynamics in 3-dimensional (3D) space,
we focus here on a 2D context for computational convenience, presuming
uniformity in the third dimension. The initial condition for the system is
shown in figure 1. The spatial domain $\Omega$ is defined as a square of width
102.4 nm ($\Omega=[0,102.4]\times[0,102.4]$). The binary alloy is initialized
in the lower part of the domain, while the liquid phase is on top. The alloy
initial composition is nominally $(c_{\text{A}},c_{\text{B}})=(0.3,0.7)$
everywhere, with superimposed white noise sampled from a uniform distribution
in the interval $[-0.025,0.025]$ at each grid-point. The liquid is initialized
as pure $c_{\text{C}}=1$. Dirichlet boundary conditions of $\phi=0$ and
$c_{\text{C}}=1$ are enforced along the top edge of the simulation. Neumann
boundary conditions (zero normal gradient) are enforced along the bottom edge
of the simulation. Periodic boundary conditions are enforced on the left and
right edges of the domain. We discretize the system using a regular
512$\times$512 grid, with a grid spacing of $\Delta x=0.2$ nm.

Figure 1: Liquid metal dealloying with two sets of initial conditions. The
initial species fields ($c_{A}$ and $c_{B}$ at $t=0$) are contaminated with
low-amplitude random white noise (left fields). Given the chaotic nature of
the dealloying process, the initial noise perturbation eventually leads to
widely different solid phase fields at late time (e.g. right fields, at
$t=6\,\mu$s, after 6 million time steps).

The evolution of the phase field is defined using

$\frac{\partial\phi}{\partial
t}=-\tilde{M}_{\phi}\frac{\pi^{2}}{8\eta}\frac{\delta F}{\delta\phi}$ (1)

where $\tilde{M}_{\phi}$ is the mobility of the solid-liquid interface, $\eta$
is the diffuse interface width, and ${\delta F}/{\delta\phi}$ is the
functional derivative of the free energy functional, $F$, with respect to the
solid phase [3, 74]. Conserved-species mole fraction evolution is governed by
the continuity equation

$\frac{\partial c_{i}}{\partial
t}=\nabla\cdot\sum_{j=\text{A,B}}{M_{ij}(\phi)\nabla\bigg{(}\frac{\delta
F}{\delta c_{j}}\bigg{)}}$ (2)

with a phase-dependent solute mobility $M_{ij}$ and the functional derivative
${\delta F}/{\delta c_{j}}$. Eq. (2) is solved for species A and B evolution,
while species C can always be determined by the condition
$c_{\text{C}}=1-c_{\text{A}}-c_{\text{B}}$. The summation in Eq. (2) is only
accomplished for species A and B since ${\delta F}/{\delta c_{j}}$ represents
the diffusion potential of species $j$ with respect to species C. The solute
mobility is given by

$M_{ij}(\phi)=D_{\text{liq}}(1-\phi)\frac{V_{a}}{kT}c_{i}(\delta_{ij}-c_{j})$
(3)

where the diffusivity is assumed to be negligibly small in the solid phase and
smoothly increases to $D_{\text{liq}}$ in the liquid phase. $V_{\text{a}}$ is
an atomic volume, $k$ is the Boltzmann constant, and $T$ is temperature.

The free energy functional is defined as

$F=\int_{\Omega}(f_{\text{phase}}(\phi)+f_{\text{chem}}(\phi,c_{\text{A}},c_{\text{B}}))dV$
(4)

with the free energy density $f_{\text{phase}}$ that sets the diffuse
interface, defined as

$f_{\text{phase}}=\frac{4\sigma_{\text{sl}}}{\eta}\bigg{[}\frac{\eta^{2}}{\pi^{2}}(\nabla\phi)^{2}+\phi(1-\phi)\bigg{]}$
(5)

where $\sigma_{\text{sl}}$ is the interfacial energy of a pure solid-liquid
interface. The chemical energy density $f_{\text{chem}}$ is given by a regular
solution model, which has an ideal entropy of mixing and non-ideal enthalpy of
mixing.

$f_{\text{chem}}=\sum_{i=\text{A,B,C}}\bigg{[}\frac{kT}{V_{\text{a}}}c_{i}\ln(c_{i})\bigg{]}+\Omega_{\text{AC}}c_{\text{A}}c_{\text{C}}+h(\phi)\Delta
g_{\text{sl}}+\frac{1}{2}\sum_{i=\text{A,B,C}}\big{[}\kappa(\nabla
c_{i})^{2}\big{]}$ (6)

where $\Omega_{\text{AC}}$ is the excess enthalpy of mixing between species A
and C, $h(\phi)$ is an interpolation function, $\Delta g_{\text{sl}}$ is the
difference between solid and liquid reference chemical energies, and $\kappa$
is the energy penalty for spatial gradients of each species. The interpolation
function, defined to smoothly vary between $h(\phi=1)=1$ in the solid and
$h(\phi=0)=0$ in the liquid, is given by

$h(\phi)=\frac{1}{2}+\frac{2}{\pi}\bigg{[}(2\phi-1)\sqrt{\phi(1-\phi)}+\frac{1}{2}\arcsin{(2\phi-1)}\bigg{]}$
(7)

while $\Delta g_{\text{sl}}$ is defined as

$\Delta
g_{\text{sl}}=\sum_{i=\text{A,B,C}}\bigg{[}c_{i}L_{i}\bigg{(}\frac{T-T_{i}}{T_{i}}\bigg{)}\bigg{]}$
(8)

where $L_{i}$ is the latent heat of melting and $T_{i}$ is the melting
temperature for each species $i$.

The parameters for the model are given in Table 1. In this work, species A is
assumed to have a very small solubility with the liquid C bath, which is set
by the large, positive excess enthalpy of mixing $\Omega_{\text{AC}}\gg 0$ and
the large melting temperature of species A, $T_{\text{A}}\gg T$. This leads to
the dealloying condition where species B is selectively dissolved from the
alloy into the liquid C bath, while A rearranges laterally along the interface
to form a topologically-complex solid structure. Importantly, diffusion along
the interface (where $0<\phi<1$) is significantly fast, according to Eq. (3),
to enable reorganization of species A. We assume the dealloying is rate-
limited by diffusion in the liquid phase (based on [11, 16]), and therefore
set the interface mobility $\tilde{M}_{\phi}$ to be very large so that
dissolution is fast compared to diffusion in the liquid.

|  |  |  |
---|---|---|---|---
Parameter | Value |  | Parameter | Value
$T$ | 1775 K |  | $\sigma_{\text{sl}}$ | 0.2 J$\cdot$m-2
$\eta$ | 4 $\times 10^{-9}\text{m}$ |  | $\kappa$ | 2.4 $\times 10^{-9}$ J$\cdot$m-1
$L_{\text{A}}$ | 2.82 $\times 10^{9}$ J$\cdot$m-3 |  | $L_{\text{B}}$ | 1.89 $\times 10^{9}$ J$\cdot$m-3
$L_{\text{C}}$ | 1.84 $\times 10^{9}$ J$\cdot$m-3 |  | $T_{\text{A}}$ | 3290 K
$T_{\text{B}}$ | 1941 K |  | $T_{\text{C}}$ | 1358 K
$V_{\text{a}}$ | 0.01 $\times 10^{-27}\text{m}^{3}$ |  | $\Omega_{\text{AC}}$ | 1.44 $\times 10^{10}$ J$\cdot$m-3
$\tilde{M}_{\phi}$ | 12.0 m$\cdot$s-1$\cdot$GPa-1 |  | $D_{\text{liq}}$ |  $7\times 10^{-9}$ m2$\cdot$s-1
|  |  |  |
Table 1: Thermodynamic and kinetic parameters used in the high fidelity phase
field simulations [16, 20]

In the high-fidelity phase field simulations, non-conserved phases are updated
using a simple forward Euler time-integration of Eq. (1). However, for
conserved species, Eq. (2) is numerically very stiff, and using forward Euler
would require prohibitively small time-steps to maintain numerical stability.
This is, in part due to fine spatial-resolution of the numerical mesh and also
the fourth-derivative on each species field in Eq. (2) that arises from the
Laplacian in Eq. (6). Instead of using forward Euler, conserved species are
updated using a semi-implicit time-integration of Eq. (2), which is based on
the approach of [75] and uses spectral methods. The reader is referred to [20]
for details of this semi-implicit time-integration scheme and implementation.
These high fidelity phase field simulations use a time-step size of $\Delta
t=10^{-12}$ s.

##  3 Leaping forward in time with U-AFNO

In this section, we introduce our proposed U-AFNO model. We first cover
attention mechanisms and Fourier neural operators (section 3.1), and further
motivate the need for a U-Net backbone augmentation (section 3.2). Finally, we
also briefly introduce hybrid high-fidelity and surrogate schemes (section
3.3).

###  3.1 Adaptive Fourier Neural Operator

Vision transformers (ViT) [66] have become widely popular in computer vision
[76, 77], and recent attempts to use them for operator learning have shown
promising early results [78, 79, 80, 69]. Originally introduced for natural
language processing applications [81], the key idea of transformers is to
compute a similarity metric (attention) across sub-samples of the training
data (tokenisation). In ViTs, this is done by splitting the input image
(denoted $X$) into $N$ smaller patches $x\in\mathbb{R}^{d}$ (each patch is
reformatted into a vector of dimension $d$). The input image can thus be
represented by a tensor $X\in\mathbb{R}^{N\times d}$. The attention can be
computed with the following expression:

$Att:\mathbb{R}^{N\times d}\mapsto\mathbb{R}^{N\times d}\hskip
36.135ptAtt(X)=\sigma\bigg{(}\frac{XW_{q}(XW_{k})^{\top}}{\sqrt{d}}\bigg{)}XW_{v}$
(9)

Where $\sigma$ is an activation function, generally taken as softmax. The
query, key and value matrices, respectively $W_{q}$, $W_{k}$ and $W_{v}$ are
learnable weights. Since ViT are designed for 2D image inputs, Eq. (9)
complexity is quadratic in nature. Thus, the self-attention mechanism of ViTs
scales poorly with the patch size and the resolution of the input image (or
equivalently in this paper, the input physical field). This is particularly
problematic for two reasons. First, high resolution images of physical fields
are sometimes necessary to capture small-scale local features. Second, a small
patch size is desirable as it generally yields higher prediction accuracy
[82]. As a result, training ViTs can be computationally challenging. To
alleviate this issue, Adaptive Fourier Neural Operators (AFNOs) [71] have been
proposed. The main idea of AFNOs is to formalize Eq. (9) as a kernel
summation. Each row of $Att(X)$ can be re-written as:

$Att(X)(s,:)=\sum_{t=1}^{N}X(t,:)\cdot\underbrace{K(s,t)W_{v}}_{\kappa(s,t)}\hskip
36.135ptK=\sigma\bigg{(}\frac{XW_{q}(XW_{k})^{\top}}{\sqrt{d}}\bigg{)}$ (10)

Where $t\in[\\![1,N]\\!]$ represents a dummy row index (in this section only).
By assuming that the kernel is a Green’s kernel, i.e.,
$\kappa(s,t)=\kappa(s-t)$, self attention can be seen as a global convolution
operation. Building on an idea originally introduced with Fourier Neural
Operators (FNOs) [51], the attention operation is performed in the Fourier
space, which allows reducing the convolution to a simpler and computationally
faster element-wise matrix product:

$Att(X)(s,:)\Longleftrightarrow\mathcal{F}^{-1}(\mathcal{F}(\kappa)\cdot\mathcal{F}(X))(s,:)$
(11)

Where $\mathcal{F}$ refers to the discrete FFT operation. Additional details
on AFNOs can be found in [71]. It should be noted that unlike other popular
neural operator models, such as FNOs [51] and DeepONets [52], AFNOs are not
neural operators in the true sense. Indeed, neural operators are generally
considered to be resolution independent. While the operation described in Eq.
(11) is resolution independent by construction, in AFNOs the input image is
processed through an embedding layer before computing attention. This
embedding layer is generally fully connected, with a fixed dimension,
resulting in the loss of resolution invariance.

###  3.2 U-AFNO: Adaptive Fourier Neural Operator with a U-Net Backbone

AFNOs have been successfully employed to build surrogates in a variety of
physical simulation problems, such, as and most notably, weather forecasting
through auto-regressive roll-out [60, 61]. However, AFNOs are prone to exhibit
squared-shape artifacts and mismatch between patch interfaces in the predicted
field. In our experiments, we have found this to be particularly true with the
binary phase fields describing LMD simulations. To address this issue and
ensure both meaningful feature extraction and smooth outputs, we propose to
combine a simple AFNO with a U-Net backbone (U-AFNO). U-Nets [65] employ an
encoder-decoder architecture, with U-shaped successive convolutions,
downsampling, and upsampling layers, and skipped connections between the
encoder and decoder. U-Nets have proven to be highly successful in image
regression tasks, and combining the strengths of U-Nets with neural operators
has been recently found to be useful [63, 67, 68, 73]. In particular, ViTO
[69] and DiTTO [70] studies have proposed to use vision transformer blocks
within the standard U-Net bottleneck, along with a trunk network for enforcing
resolution independence (as originally introduced with DeepONets [52]). In the
present study, we adopt a similar, yet simpler, approach by replacing the
U-Net bottleneck with an AFNO. Note that downsampling the input image through
the U-Net encoder has the convenient advantage of allowing for patch sizes of
$1\times 1$ (without downsampling, such small patch size would be
computationally intractable).

Figure 2: U-AFNO Architecture - The model takes as input the field at time $t$
and outputs the field at a future time step $t+\Delta\tau$. The U-Net color
scheme is analogous to the one employed in the original U-Net paper [65]. The
abbreviations in the AFNO layers are the following: FFT: Fast Fourier
Transform, IFFT: Inverse Fast Fourier Transform, MHA: Mutli-Head Attention,
FC: Fully Connected.

The architecture of the U-AFNO is represented in Figure 2. The input field has
dimensions $3\times 512\times 512$ (the 3 channels correspond to $\phi$,
$c_{A}$ and $c_{B}$ respectively – note that $c_{C}$ is not needed here since
its computation knowing $c_{A}$ and $c_{B}$ is trivial), and is downsampled to
dimensions $256\times 64\times 64$ through the U-Net encoder. The AFNO
bottleneck is composed of 12 back-to-back Attention blocks, each of which
employs 16 attention heads, a patch size of $1\times 1$, a skipped connection
and two fully-connected layers with 3072 neurons each. The AFNO output is then
upsampled back into the original dimension through the U-Net decoder. Both the
U-Net and the AFNO layers use GeLU activation functions, and the final output
of the U-AFNO is passed through a sigmoid (logistic) function. This guarantees
that the U-AFNO predictions are bounded between 0 and 1, and hence consistent
with the phase field physics.

###  3.3 Hybrid High-Fidelity and Surrogate Model Coupled Solver

Neural operators and adjacent regression models have been extensively used
with an auto-regressive roll-out to predict PDE dynamics over time [51, 52,
50, 60]. While this approach is simple and fast, errors tend to grow quickly,
potentially unbounded. Similarly to traditional explicit numerical solvers,
auto-regressive roll-out may exhibit instabilities and/or un-physical behavior
over time. This is in part because the outputs of neural operators are
generally not guaranteed to lie within the correct manifold of physical fields
(especially with fully data-driven operators). In order to slow-down error
growth, a recently proposed approach [73, 62] alternates forward passes
through the neural operator (forward leap in time) with high-fidelity time
stepping (relaxation steps). This strategy allows for greater flexibility in
balancing the trade-off of speed-up vs. accuracy. More relaxation steps means
greater reliance on the high-fidelity solver and fewer neural operator forward
passes, which mitigates error build-up (as each new operator forward pass
inevitably introduces some error), but also takes more time. Conversely, fewer
relaxation steps means faster surrogate predictions but larger error
accumulation in the long term. In this paper, we investigate the performance
of U-AFNOs with and without these high-fidelity relaxation steps (in the
former case, we employ the term hybrid simulation),

##  4 Phase Field Quantities of Interest

An essential motivation for developing physical simulations (whether high-
fidelity or surrogate ones) is that it is often a necessary first step for
extracting informative quantities of interest (QoIs) from the physical fields.
These QoIs, typically computed locally and/or globally across the spatial
field and over time may be crucial for decision making, optimization or
uncertainty quantification purposes. In some applications, the extracted QoIs
may be more important than the simulated physical field itself. In dealloying
corrosion, some QoIs are particularly informative, such as the local curvature
[10, 12, 83] and the total length (perimeter) [84] of the solid-liquid
interface, the penetration depth of the corrosive liquid into the metal [19,
11], and the amount of remaining metal, species A and B, throughout the
corrosion process.

###  4.1 Curvature and Perimeter

The first collection of QoIs we consider involve the curvature of the liquid-
metal interface. We formally define the interface in this 2D system as any
spatial coordinate $(x,y)\in\Omega$ such that $\phi(t,x,y)=0.5$, $\forall t$.
The interface coordinates may then be described at any time $t$ as a
parameterized curve, $\gamma(s,t)=(x(s,t),y(s,t))$ with
$s\in[0,1]\subset\mathbb{R}$, and the signed curvature along the interface is:

$k(s,t)=\frac{x^{\prime}y^{\prime\prime}-y^{\prime}x^{\prime\prime}}{(x^{\prime
2}+y^{\prime 2})^{3/2}}\hskip 21.68121ptx^{\prime}\equiv\frac{\partial
x}{\partial s}\hskip 21.68121pty^{\prime}\equiv\frac{\partial y}{\partial s}$
(12)

The curvature is defined locally, but we can quantify it globally by
considering statistics such as the mean and standard deviation of the absolute
curvature over the entire interface, at any moment in time:

$\mu_{k}(t)=\int_{0}^{1}|k(s,t)|\text{d}s\hskip
36.135pt\sigma_{k}(t)=\bigg{[}\int_{0}^{1}(|k(s,t)|-\mu_{k}(t))^{2}\text{d}s)\bigg{]}^{1/2}$
(13)

Furthermore, we also define the interface perimeter at any time as:

$p(t)=\int_{0}^{1}\sqrt{(x^{\prime})^{2}+(y^{\prime})^{2}}\text{d}s$ (14)

Eqs. (12), (13) and (14) are all computed discretely using finite difference.
Figure 3 represents the general approach to compute the curvature locally at
each time step, and extract statistical quantities to represent more global
curvature dynamics over time. Note that in liquid metal dealloying, metal
ligaments may fully separate from the rest of the solid phase and form
islands. Similarly, some metal ligaments may rejoin and close-up, trapping
liquid bubbles inside the metal. In some cases, these features can shrink via
capillary forces until they are eventually eliminated. During such events, the
mean and standard deviation of the curvature is naturally skewed upward,
before sharply dropping once an island or bubble disappears. Such patterns can
be observed in Figure 3.

###  4.2 Penetration Depth

The second collection of QoIs considered describe corrosive liquid penetration
within the metal. As the species C-rich liquid progressively penetrates into
the metal, the latter reorganizes to form narrow ligaments that typically grow
as the liquid penetrates further into the bulk metal. We can quantify the
average ligament length over time, as well as the maximum penetration depth,
by considering the local extrema of the parametric curve representing the
interface (along the vertical axis). The set of local minima and maxima at
time $t$ are defined as:

$\mathcal{S}_{\text{min}}(t)=\\{y(s,t)\,\,|\,\,y^{\prime}=0\,\,\,\text{and}\,\,\,y^{\prime\prime}>0\\}\hskip
36.135pt\mathcal{S}_{\text{max}}(t)=\\{y(s,t)\,\,|\,\,y^{\prime}=0\,\,\,\text{and}\,\,\,y^{\prime\prime}<0\\}$
(15)

Each local maximum is representative of a ligament tip $y-$axis coordinate,
while each local minimum is representative of a ligament base $y-$axis
coordinate. Thus, we approximate the mean ligament height $\mu_{d}$ as the
difference between the mean of the set of local minima and local maxima, thus

$\mu_{d}(t)=\mathbb{E}[\mathcal{S}_{\text{max}}(t)]-\mathbb{E}[\mathcal{S}_{\text{min}}(t)].$
(16)

We also define the maximum penetration depth as

$\text{max}_{p}(t)=1-\min(\mathcal{S}_{\text{min}}(t)).$ (17)

Figure 3: Curvature QoIs - At each time $t$, the parametric curve
$\gamma(s,t)$ delineating the interface between the two phases is computed,
from which we compute the curvature $k(s,t)$. Histograms of the curvature
distribution at two different times ($t=2.5~{}\mu$s and $t=3.5~{}\mu$s) are
shown, from which means and standard deviations are estimated and computed
over time (blue and orange plots).

###  4.3 Mass Loss

The third and last set of QoIs considered describes the dynamics of lost
material throughout the corrosion process. We define the total “mass" of
metal, species A, and species B as the integrals over the spatial domain of
$\phi$, $c_{A}$ and $c_{B}$ respectively, thus

$m_{\phi}(t)=\int_{\Omega}\phi(t,x,y)\text{d}x\text{d}y\hskip
36.135ptm_{A}(t)=\int_{\Omega}c_{\text{A}}(t,x,y)\text{d}x\text{d}y\hskip
36.135ptm_{B}(t)=\int_{\Omega}c_{\text{B}}(t,x,y)\text{d}x\text{d}y$ (18)

The actual masses in proper mass-units are related to these quantities with
proportionality constants, these being the metal density in the case of
$m_{\phi}$, and the molar weights for $m_{\text{A}}$ and $m_{\text{B}}$.

##  5 Application

In this section, we present and discuss the performances of our U-AFNO model,
both for predicting the full micro-structure fields as well as the QoIs
described in the previous section. We also provide details about the model
training and review the different error metrics used in this study.

###  5.1 Training

The training set comprises 87 high-fidelity simulations, each containing 100
field snapshots spanning from time $t=1\mu$s to $t=6\mu$s. The U-AFNO is
trained to leap $50000$ (high-fidelity) time steps ahead, i.e.
$\Delta\tau=5\cdot 10^{4}\times\Delta t=0.05\mu$s every forward pass. The test
set comprises 20 high-fidelity simulations of similar form. The model is
trained for 20 epochs, using the Adam optimizer [85] with a learning rate of
$10^{-4}$. The code is implemented with Pytorch [86] and the AFNO bottleneck
is based on the FourCastNet implementation [60]. For baseline comparison, we
also train a standard U-Net, an AFNO (with no U-Net backbone), and an FNO. All
our models are trained on a single Nvidia A100 80GB GPU. The U-AFNO training
took 35 hours.

In this paper, we run test surrogate simulations for up to time step $6\times
10^{6}$ (i.e. $t=6\mu$s), all starting with quasi-identical initial
conditions. We initialize the surrogate simulations with time integration of
the high-fidelity model (in the following sections,
$n_{\text{init}}=$1000000$$ steps). This allows for the instabilities at the
interface between the phases to grow and be sufficiently distinguishable
across different initial conditions. We test both fully auto-regressive
simulations (i.e. high-fidelity relaxation steps are not employed, except
during initialization), and hybrid simulations with different numbers of high-
fidelity relaxation steps (as presented in section 3.3).

###  5.2 Error Metrics

We consider several error metrics to assess the accuracy of our model. While
the relative error between the ground truth (i.e. high-fidelity simulation)
and the predicted field is a common error metric in the neural operator
literature [51, 52, 50], it is not directly applicable here due to the chaotic
nature of LMD simulations. In chaotic systems, comparing invariant statistical
properties, such as field auto-correlations is more meaningful and informative
[72]. In phase field simulations, the field spatial auto-correlation is
equivalent to the probability that two points in the field will be part of the
same phase [87]. As such, the spatial auto-correlation relative error has been
shown to be more relevant for describing the micro-structure statistics than
other relative error metrics, and has been recently used in neural-operator-
based phase field models [73]. The spatial auto-correlation relative error is
defined as:

$e_{\text{AC}}(\mathbf{u},\mathbf{\hat{u}},t)=\frac{|\\!|S_{\mathbf{\hat{u}}\mathbf{\hat{u}}}(\boldsymbol{r},t)-S_{\mathbf{u}\mathbf{u}}(\boldsymbol{r},t)|\\!|_{2}}{|\\!|S_{\mathbf{u}\mathbf{u}}(\boldsymbol{r},t)|\\!|_{2}}$
(19)

where $\mathbf{u}$ can represent either a ground truth species mole fraction
or phase field variable at time $t$, and $\mathbf{\hat{u}}$ represents the
corresponding surrogate prediction. $S$ is the corresponding spatial auto-
correlation function and $\boldsymbol{r}$ is a vector connecting any two
points in the spatial field. Note that, for baseline comparison, we also
compute the average auto-correlation relative discrepancy between each pair of
ground truth simulation contained in the test set. In this later case,
equation 19 is used as well, but we employ the term discrepancy rather than
error since we only compare ground truth simulations together. This
discrepancy metric allows for meaningful judgement of the difference between
ground truth dynamics and our surrogate, compared to a high-fidelity
simulation with a different random initial condition. To assess the ability of
our model to predict the key QoIs described in section 4, we also consider a
second set of error metrics, the relative error of each QoI across all time
steps. For example, we define the mean curvature relative error as:

$e(\boldsymbol{\mu_{k}},\boldsymbol{\hat{\mu}_{k}})=\frac{|\\!|\boldsymbol{\hat{\mu}_{k}}-\boldsymbol{\mu_{k}}|\\!|_{2}}{|\\!|\boldsymbol{\mu_{k}}|\\!|_{2}}$
(20)

where $\boldsymbol{\mu_{k}}$ is the vector of ground truth mean curvature at
each time step:
$\boldsymbol{\mu_{k}}=(\mu_{k}(t_{0}),\dots,\mu_{k}(t_{\text{max}}))$ (and
$\boldsymbol{\hat{\mu}_{k}}$ is the equivalent vector for the surrogate
prediction). In a similar fashion, we define the relative error for the
curvature standard deviation, the interface perimeter, the maximum penetration
depth, the mean ligament height, and the remaining metal mass.

###  5.3 Results

We first consider the fully auto-regressive case. Figure 4 shows the solid
phase $\phi$ and species A and B predicted with a U-AFNO-B/1 surrogate (where
B/1 refers to a $1\times 1$ patch size), at different moments in time from
$t=0\,\mu$s to $t=6\,\mu$s. The total number of auto-regressive forward passes
through the U-AFNO is $100$. Due to the chaotic behavior of the LMD physics,
the surrogate model is not able to reproduce the exact ligament shapes
available from the ground truth data, but it captures patterns that are
visually consistent with them. Specifically, based on visual inspection, the
ligaments of the solid phase exhibit satisfactory thickness, height and size.
Further, the liquid phase penetrates into the metal in a manner consistent
with the ground truth results, making its way steadily towards the bottom
boundary. The dynamics of both species are similarly consistent; especially at
the solid-liquid interface where species A tends to segregate to the growing
ligaments while species B is dissolved into the liquid. This behavior,
consistent with dealloying phenomena, is correctly reproduced here. This
empirical observation can be confirmed by looking at the associated auto-
correlation maps, represented in Figure 5. The predicted auto-correlations
closely match the ground truth, indicating that the U-AFNO reproduces
invariant field statistics well.

The auto-correlation relative error metrics for each surrogate model are
represented in Figure 6. The auto-correlation error for the solid phase, and
species A and B are shown every $500000$ time steps. For each field, the
U-AFNO-B/1 clearly outperforms all the baseline models. At the end of the
simulation, the relative error for $\phi$ and $c_{\text{A}}$ does not exceed
$15\%$, while that for $c_{\text{B}}$ does not exceed $20\%$. Hence, even
after 100 forward pass through the U-AFNO model, the error remains remarkably
well controlled. When employing a larger patch size ($2\times 2$ and $4\times
4$), the errors slightly deteriorate, which is expected since larger
transformer patch sizes generally yield worse accuracy. Yet, each U-AFNO model
outperforms the other baseline models. In particular, they are more accurate
than both U-Nets and vanilla AFNOs, confirming that combining these two model
architectures can significantly improve prediction accuracy. The U-AFNO-B/1
relative error is on-par with the average discrepancy between distinct pairs
of ground truth simulation (shown in hatched black). This is significant: in
practice, it means that, on the auto-correlation relative error metric, and
considering a given high-fidelity simulation, our surrogate will reproduce the
dynamics equally well as any other high-fidelity run (on perturbatively
different initial condition).

The QoI relative errors are provided in Table 2, and the predicted dynamics of
each QoI are shown in Figure 7 and Figure 8 for the U-AFNO-B/1 and the U-Net,
respectively. Note that the QoIs and error metrics are all computed across 20
test simulations, therefore the values in Table 2 and curves in Figure 7 and 8
are mean values (with standard deviations). The U-AFNO clearly outperforms the
baseline models and predicts the evolution of each QoI with very reasonable
accuracy, as it rarely exceeds $5$ to $10\%$ relative error (except for the
mean ligament height, with a mean relative error of $36.5\%$). As shown in
Figure 7, our proposed model reproduces the mean curvature and mass envelopes
with remarkable accuracy, and is also able to capture the curvature standard-
deviation, interface perimeter and maximum penetration depth true patterns
very well. Conversely, in Figure 8, the U-Net underestimates the mean
curvature and mass envelope values, while overestimating the curvature
standard-deviation, interface perimeter and maximum penetration depth. Most
notably, with the U-Net, the liquid penetrates into the solid exceedingly fast
and reaches the bottom boundary much earlier then it should (the maximum
penetration depth steeply rises after 2 $\mu$s). This is due to the fact that
with the U-Net, the liquid-solid interface tends to collapse at both the right
and left edges of the domain. It is not clear why this phenomena occurs, but a
possible explanation is that the U-Net fails to accurately capture the effects
of the boundary conditions on the field dynamics.

Figure 4: U-AFNO-B/1 field predictions in the fully auto-regressive case. The
predicted and ground truth fields are shown at $t=(0,1,2,3,4,5,6)~{}\mu$s.
Since the initial $10^{6}$ time steps are computed with the high-fidelity
solver, time $t=0-1~{}\mu$s for the U-AFNO and the ground truth are identical.
Figure 5: U-AFNO field auto-correlation maps in the fully auto-regressive
case. The auto-correlation maps correspond to the fields shown in figure 4.
The error maps are the absolute error between the predicted and ground truth
auto-correlation. Figure 6: Auto-correlation relative errors for fully auto-
regressive simulations. Each model is trained to skip $50000$ time steps at-a-
time, but for simplicity we only show the error values every $500000$ steps.
The hatched grey regions represent the initial $10^{6}$ time steps computed
with the high fidelity solver. The cross-hatched black lines in each plot
represent the average auto-correlation discrepancy between each distinct pair
of ground truth simulations in the test set.

|  |  |  |  |  |  |  |
---|---|---|---|---|---|---|---|---
Model |  Mean Curvature |  Curvature Standard Deviation |  Interface Perimeter | Total Mass |  $c_{A}$ Mass |  $c_{B}$ Mass |  Maximum Penetration Depth |  Mean Ligament Height
U-AFNO-B/1 (Ours) | $\boldsymbol{0.057\pm 0.019}$ | $\boldsymbol{0.117\pm 0.032}$ | $\boldsymbol{0.068\pm 0.029}$ | $\boldsymbol{0.011\pm 0.006}$ | $\boldsymbol{0.023\pm 0.011}$ | $\boldsymbol{0.007\pm 0.003}$ | $\boldsymbol{0.040\pm 0.017}$ | $0.358\pm 0.238$
U-AFNO-B/2 (Ours) | $0.082\pm 0.030$ | $0.128\pm 0.044$ | $0.182\pm 0.038$ | $0.032\pm 0.009$ | $0.057\pm 0.017$ | $0.014\pm 0.004$ | $0.055\pm 0.023$ | $0.345\pm 0.096$
U-AFNO-B/4 (Ours) | $0.058\pm 0.017$ | $0.122\pm 0.030$ | $0.073\pm 0.033$ | $0.036\pm 0.012$ | $0.055\pm 0.019$ | $0.029\pm 0.004$ | $0.052\pm 0.022$ | $\boldsymbol{0.328\pm 0.150}$
U-Net | $0.069\pm 0.022$ | $0.144\pm 0.043$ | $0.202\pm 0.076$ | $0.096\pm 0.015$ | $0.119\pm 0.026$ | $0.033\pm 0.006$ | $0.364\pm 0.037$ | $0.368\pm 0.165$
FNO | $1.862\pm 0.290$ | $3.798\pm 0.310$ | $0.282\pm 0.073$ | $0.107\pm 0.021$ | $0.164\pm 0.031$ | $0.054\pm 0.019$ | $0.326\pm 0.107$ | $0.841\pm 0.047$
AFNO-B/16 | $0.790\pm 0.092$ | $1.895\pm 0.173$ | $0.132\pm 0.053$ | $0.123\pm 0.026$ | $0.106\pm 0.022$ | $0.133\pm 0.036$ | $0.040\pm 0.016$ | $0.366\pm 0.146$
AFNO-B/32 | $0.777\pm 0.081$ | $2.110\pm 0.269$ | $0.286\pm 0.085$ | $0.121\pm 0.035$ | $0.211\pm 0.049$ | $0.045\pm 0.032$ | $0.123\pm 0.059$ | $0.581\pm 0.109$
|  |  |  |  |  |  |  |

Table 2: QoIs relative errors with fully auto-regressive roll-out. For U-AFNOs
and vanilla AFNOs, B/1, B/2 and B/4 indicate a patch size of $1\times 1$,
$2\times 2$ and $4\times 4$, respectively. U-AFNOs can employ significantly
smaller patch sizes than AFNOs, since the U-Net encoder architecture
downsamples the AFNO input. Note that for fair comparison, every model output
is wrapped into the interval $[0,1]$, as described in section 3.2 Figure 7:
QoI dynamics predicted with the U-AFNO-B/1 (fully auto-regressive case). For
each QoI, the solid blue line is the mean QoI value over time, taken across
all the ground truth simulations in the test set (and the dotted lines are the
corresponding standard deviations). The solid red lines represent the QoI
dynamics predicted by the U-AFNO (averaged over all test surrogate
simulations). Figure 8: QoI dynamics predicted with the U-Net (fully auto-
regressive case). For each QoI, the solid blue line is the mean QoI value over
time, taken across all the ground truth simulations in the test set (and the
dotted lines are the corresponding standard deviations). The solid red lines
represent the QoI dynamics predicted by the U-Net (averaged over all test
surrogate simulations).

We now investigate the effects of adding hybrid time stepping, as briefly
described in section 3.3. We use the same trained models as in the fully auto-
regressive case (leaping $50000$ time steps at-a-time), but after each forward
pass through the surrogate model, we now employ the high-fidelity solver to
advance in time for a limited number of relaxation time steps. For both the
U-AFNO and the baseline models, we run test simulations for $1000$, $10000$
and $50000$ relaxation time steps. Figure 9 and Figure 10 show the auto-
correlation relative errors with and without relaxation steps, for the U-AFNO
and the U-Net, respectively. With the U-AFNO, the accuracy of the surrogate
simulation does not significantly improve when using high-fidelity relaxation
time stepping. While this may seem surprising, it is in fact reasonable since
the margin for improvement is limited. Indeed, as seen in the fully auto-
regressive case, the U-AFNO average error is already on par with the average
auto-correlation discrepancy between distinct ground truth simulations. Thus,
reducing the error further with a model that still relies on surrogate
modeling (whether moderately or not) is unlikely to succeed.

For the U-Net however, we observe a different trend. As seen in Figure 10, the
error significantly decreases when employing high-fidelity relaxation steps,
likely because there is more room for improvement (compared to the fully auto-
regressive case). With just 1,000 relaxation time steps, the errors sharply
decrease (especially for $\phi$ and $c_{B}$), indicating that even just a few
relaxation time steps are sufficient to correct the field physics. When
increasing the number of relaxation time steps, the errors eventually decrease
further, but with diminishing returns. This is expected since (as discussed in
the previous paragraph), it is unlikely for any surrogate model to achieve
significantly lower errors than the average auto-correlation discrepancy
between ground truth simulations. The observation that blending U-Net forward
passes with high-fidelity time stepping to improve accuracy is very consistent
with findings in Oommen et al. [73]. However, as seen with our proposed U-AFNO
model and in Figure 9, the gains from employing high-fidelity time stepping
may vary depending on the specific surrogate model architecture.

Tables 3, 4, and 5 show the QoI relative errors for $1000$, $10000$ and
$50000$ high-fidelity relaxation time steps, respectively. While the U-AFNO
still generally outperforms the baseline models (in particular with a patch
size of 1$\times$1), the improvements compared to the fully auto-regressive
case are marginal, and sometimes the accuracy even degrades (as already
observed with the auto-correlation errors in Figure 9). For example, the
relative error for the mass of species B increases from $0.7\%$ in the fully
auto-regressive case to $3.2\%$ in the hybrid (1000 relaxation steps) case.
With some of the baseline models, such as the FNO and the vanilla AFNOs, the
accuracy slightly improves between the fully auto-regressive case and the
hybrid case (with 1000 relaxation steps), but then significantly deteriorates
with $10000$ and $50000$ relaxation steps. In these latter two cases, the
metal tends to be consumed by the corrosive liquid much faster than it
normally should, and fully disappears early-on. As a result, the liquid-metal
interface is not properly defined anymore, and the QoIs depending on this
interface are unavailable (e.g. curvature QoIs, etc.). Intuitively, the
accuracy should improve when employing more and more relaxation time steps, as
a larger and larger portion of the simulation is computed through the high-
fidelity solver. However as discussed earlier, this assertion does not always
hold true. Except for the vanilla U-Net, the accuracy either remains stable,
or worsens. Additional investigation as to why and when accuracy may degrade
is left to future work, but we conjecture that this degradation is linked to
the robustness of the high-fidelity solver to un-physical initial conditions.
Since the neural operator is purely data-driven and no physical constraints
are enforced, the neural operator outputs may violate some underlying
assumptions of the high-fidelity solver, leading to instability or un-physical
behavior build-up by the high-fidelity solver itself. For instance, in the
high-fidelity LMD simulations, periodic boundary conditions are enforced
between the left and right edges of the domain. However, these boundary
conditions are not hard-enforced during the neural operator training. Thus,
the neural operator outputs, which serve as initial conditions for the high-
fidelity solver may not be strictly periodic, leading to unexpected behaviors.

The high-fidelity solver average wall clock run time is 0.026$s$ per time step
(on 128 CPU cores). That is, 1303$s$ per $50000$ steps. The average wall clock
run-time for a forward pass through the U-AFNO (with a $1\times 1$ patch size)
is 0.116$s$ on a single GPU. Thus, the speed-up to compute $50000$ time steps
is about $11200$$\times$. The U-AFNO run time is quasi-negligible compared to
the high-fidelity solver run time, so the speed-up of a full surrogate
simulation (in the fully auto-regressive case) comes down to the time taken to
initialize the simulation. In this paper, we have run 6$\mu$s worth of real-
time simulation, where the initial 1$\mu$s is run with the high-fidelity
solver (i.e. $10^{6}$ time steps), and the remaining 5$\mu$s are obtained with
U-AFNO auto-regressive roll-out. This is equivalent to a 6$\times$ speed-up.
Note that augmenting the surrogate simulation with high-fidelity relaxation
steps would considerably lower the speed up gains, and as discussed in the
previous paragraph, the U-AFNO already achieve excellent accuracy without
needing such hybrid time stepping.

|  |  |  |  |  |  |  |
---|---|---|---|---|---|---|---|---
Model |  Mean Curvature |  Curvature Standard Deviation |  Interface Perimeter | Total Mass |  $c_{A}$ Mass |  $c_{B}$ Mass |  Maximum Penetration Depth |  Mean Ligament Height
U-AFNO-B/1 (Ours) | $\boldsymbol{0.054\pm 0.017}$ | $\boldsymbol{0.114\pm 0.035}$ | $\boldsymbol{0.083\pm 0.045}$ | $\boldsymbol{0.017\pm 0.009}$ | $\boldsymbol{0.028\pm 0.012}$ | $0.032\pm 0.006$ | $\boldsymbol{0.031\pm 0.019}$ | $\boldsymbol{0.317\pm 0.221}$
U-AFNO-B/2 (Ours) | $0.056\pm 0.017$ | $0.138\pm 0.042$ | $0.139\pm 0.047$ | $0.036\pm 0.011$ | $0.060\pm 0.019$ | $0.015\pm 0.003$ | $0.059\pm 0.024$ | $0.355\pm 0.220$
U-AFNO-B/4 (Ours) | $0.060\pm 0.017$ | $0.117\pm 0.037$ | $0.105\pm 0.063$ | $0.019\pm 0.013$ | $0.046\pm 0.021$ | $0.030\pm 0.005$ | $0.036\pm 0.021$ | $0.350\pm 0.201$
U-Net | $0.066\pm 0.024$ | $0.145\pm 0.052$ | $0.129\pm 0.062$ | $0.040\pm 0.021$ | $0.067\pm 0.032$ | $\boldsymbol{0.012\pm 0.006}$ | $0.104\pm 0.098$ | $0.381\pm 0.202$
FNO | $0.735\pm 0.077$ | $1.961\pm 0.215$ | $0.342\pm 0.053$ | $0.342\pm 0.062$ | $0.375\pm 0.060$ | $0.320\pm 0.059$ | $0.266\pm 0.104$ | $0.759\pm 0.073$
AFNO-B/16 | $0.137\pm 0.039$ | $0.471\pm 0.105$ | $0.406\pm 0.051$ | $0.250\pm 0.021$ | $0.185\pm 0.033$ | $0.336\pm 0.019$ | $0.054\pm 0.038$ | $0.433\pm 0.167$
AFNO-B/32 | $0.202\pm 0.034$ | $0.651\pm 0.113$ | $0.357\pm 0.080$ | $0.268\pm 0.057$ | $0.229\pm 0.051$ | $0.329\pm 0.066$ | $0.096\pm 0.041$ | $0.400\pm 0.123$
|  |  |  |  |  |  |  |

Table 3: QoIs relative errors with a hybrid surrogate/high-fidelity roll-out.
Here, 1000 high-fidelity relaxation time steps are employed after each forward
pass through the surrogate model. In some cases, the surrogate simulation is
too inaccurate to exhibit a properly defined liquid-metal interface, in which
case any QoI depending on this interface is unavailable (listed as N/A).
Figure 9: Auto-correlation relative errors for the U-AFNO-B/1, using either
fully auto-regressive simulations, or hybrid simulations with 1000, $10000$,
and $50000$ relaxation time steps. The cross-hatched black lines in each plot
represent the average auto-correlation discrepancy between each distinct pair
of ground truth simulations in the test set. Note that when using either 1000
or $10000$ relaxation steps, the timeline of the available snapshot is
progressively shifting. For example, in the fully auto-regressive case (or the
$50000$ relaxation steps case), the $1050000$th step and the $1100000$th step
will be consecutively available. Conversely, with $10000$ relaxation steps,
the $1050000$th, $1060000$th and $1110000$th steps will be consecutively
available. This is the reason why the 1000 relaxation steps (green) and
$10000$ relaxation step (orange) markers are slightly shifted. Figure 10:
Auto-correlation relative errors for the U-Net, using either fully auto-
regressive simulations, or hybrid simulations with 1,000, 10,000, and 50,000
relaxation time steps. The cross-hatched black lines in each plot represent
the average auto-correlation discrepancy between each distinct pair of ground
truth simulations in the test set.

|  |  |  |  |  |  |  |
---|---|---|---|---|---|---|---|---
Model |  Mean Curvature |  Curvature Standard Deviation |  Interface Perimeter | Total Mass |  $c_{A}$ Mass |  $c_{B}$ Mass |  Maximum Penetration Depth |  Mean Ligament Height
U-AFNO-B/1 (Ours) | $\boldsymbol{0.056\pm 0.016}$ | $\boldsymbol{0.114\pm 0.038}$ | $\boldsymbol{0.084\pm 0.045}$ | $\boldsymbol{0.013\pm 0.007}$ | $0.035\pm 0.015$ | $0.028\pm 0.006$ | $0.042\pm 0.025$ | $\boldsymbol{0.261\pm 0.129}$
U-AFNO-B/2 (Ours) | $\boldsymbol{0.056\pm 0.018}$ | $0.116\pm 0.033$ | $0.089\pm 0.047$ | $0.023\pm 0.010$ | $0.045\pm 0.019$ | $\boldsymbol{0.018\pm 0.004}$ | $0.047\pm 0.021$ | $0.270\pm 0.194$
U-AFNO-B/4 (Ours) | $0.057\pm 0.015$ | $0.123\pm 0.051$ | $0.110\pm 0.054$ | $0.025\pm 0.014$ | $\boldsymbol{0.034\pm 0.018}$ | $0.029\pm 0.006$ | $\boldsymbol{0.040\pm 0.027}$ | $0.291\pm 0.166$
U-Net | $0.060\pm 0.018$ | $0.125\pm 0.037$ | $0.087\pm 0.033$ | $0.044\pm 0.030$ | $0.062\pm 0.035$ | $0.021\pm 0.015$ | $0.101\pm 0.109$ | $0.336\pm 0.165$
FNO | N/A | N/A | N/A | $0.540\pm 0.055$ | $0.540\pm 0.040$ | $0.534\pm 0.079$ | N/A | N/A
AFNO-B/16 | N/A | N/A | N/A | $0.580\pm 0.124$ | $0.522\pm 0.117$ | $0.652\pm 0.127$ | N/A | N/A
AFNO-B/32 | N/A | N/A | N/A | $0.791\pm 0.020$ | $0.703\pm 0.024$ | $0.843\pm 0.015$ | N/A | N/A
|  |  |  |  |  |  |  |

Table 4: QoIs relative errors with a hybrid surrogate/high-fidelity roll-out.
Here, 10,000 high-fidelity relaxation time steps are employed after each
forward pass through the surrogate model. In some cases, the surrogate
simulation is too inaccurate to exhibit a properly defined liquid-metal
interface, in which case any QoI depending on this interface is unavailable
(listed as N/A).

|  |  |  |  |  |  |  |
---|---|---|---|---|---|---|---|---
Model |  Mean Curvature |  Curvature Standard Deviation |  Interface Perimeter | Total Mass |  $c_{A}$ Mass |  $c_{B}$ Mass |  Maximum Penetration Depth |  Mean Ligament Height
U-AFNO-B/1 (Ours) | $0.061\pm 0.018$ | $\boldsymbol{0.107\pm 0.036}$ | $0.092\pm 0.042$ | $\boldsymbol{0.029\pm 0.013}$ | $0.053\pm 0.019$ | $\boldsymbol{0.014\pm 0.004}$ | $\boldsymbol{0.058\pm 0.028}$ | $0.272\pm 0.132$
U-AFNO-B/2 (Ours) | $\boldsymbol{0.051\pm 0.013}$ | $0.108\pm 0.037$ | $\boldsymbol{0.076\pm 0.032}$ | $0.034\pm 0.011$ | $0.041\pm 0.016$ | $0.021\pm 0.004$ | $0.059\pm 0.027$ | $0.271\pm 0.151$
U-AFNO-B/4 (Ours) | $0.055\pm 0.017$ | $\boldsymbol{0.107\pm 0.036}$ | $0.109\pm 0.051$ | $0.041\pm 0.012$ | $0.048\pm 0.018$ | $0.018\pm 0.005$ | $\boldsymbol{0.058\pm 0.024}$ | $\boldsymbol{0.263\pm 0.152}$
U-Net | $0.065\pm 0.027$ | $0.110\pm 0.039$ | $0.102\pm 0.064$ | $0.030\pm 0.019$ | $\boldsymbol{0.039\pm 0.023}$ | $0.017\pm 0.006$ | $0.083\pm 0.067$ | $0.266\pm 0.101$
FNO | N/A | N/A | N/A | $0.400\pm 0.076$ | $0.394\pm 0.073$ | $0.414\pm 0.075$ | N/A | N/A
AFNO-B/16 | N/A | N/A | N/A | $0.493\pm 0.108$ | $0.443\pm 0.124$ | $0.635\pm 0.228$ | N/A | N/A
AFNO-B/32 | N/A | N/A | N/A | $0.734\pm 0.027$ | $0.669\pm 0.030$ | $0.799\pm 0.020$ | N/A | N/A
|  |  |  |  |  |  |  |

Table 5: QoIs relative errors with a hybrid surrogate/high-fidelity roll-out.
Here, 50,000 high-fidelity relaxation time steps are employed after each
forward pass through the surrogate model. In some cases, the surrogate
simulation is too inaccurate to exhibit a properly defined liquid-metal
interface, in which case any QoI depending on this interface is unavailable
(listed as N/A).

##  6 Conclusion

In this paper, we have introduced U-AFNOs, a new ML-based surrogate model for
fast prediction of time-dependant PDEs. This architecture was applied to phase
field simulations, specifically the dealloying corrosion of metals, where the
infiltration of a corrosive liquid into the alloy can rapidly lead to
morphologically-complex metal structures. Employed within an auto-regressive
roll-out, our proposed model is able to reproduce the correct LMD patterns,
and accurately captures invariant statistics, even in high-chaotic
simulations. Our work supports the idea that combining U-Nets and vision
transformers may provide very promising architectures for operator learning.
Furthermore, we have identified meaningful QoIs to accurately describe the LMD
dynamics and showed that they can be accurately predicted by our model. We
have also investigated the relevance of augmenting auto-regressive surrogate
simulations with hybrid high-fidelity time stepping. We showed that the
potential gains in accuracy are not systematic, and may be dependent on the
architecture of the ML surrogate and/or the high-fidelity solver (i.e.
problem-dependent).

## Acknowledgements

This work was supported by the U.S. Department of Energy, Office of Nuclear
Energy, and Office of Science, Office of Advanced Scientific Computing
Research through the Scientific Discovery through Advanced Computing project
on Simulation of the Response of Structural Metals in Molten Salt Environment.
This article has been co-authored by employees of National Technology and
Engineering Solutions of Sandia, LLC under Contract No. DE-NA0003525 with the
U.S. Department of Energy (DOE). The employees co-own right, title and
interest in and to the article and are responsible for its contents. The
United States Government retains and the publisher, by accepting the article
for publication, acknowledges that the United States Government retains a non-
exclusive, paid-up, irrevocable, world-wide license to publish or reproduce
the published form of this article or allow others to do so, for United States
Government purposes. The DOE will provide public access to these results of
federally sponsored research in accordance with the DOE Public Access Plan
https://www.energy.gov/downloads/doe-public-access-plan. Sandia Release
Number: SAND2024-07895O. Los Alamos National Laboratory, United States, an
affirmative action/equal opportunity employer, is operated by Triad National
Security, LLC, for the National Nuclear Security Administration of the U.S.
Department of Energy under Contract No. 89233218CNA000001. Lawrence Berkeley
National Laboratory is supported by the DOE Office of Science under contract
no. DE-AC02-05CH11231. This study made use of computational resources of the
National Energy Research Scientific Computing Center (NERSC), which is also
supported by the Office of Basic Energy Sciences of the US Department of
Energy under the same contract number.

## References

  * [1] Long-Qing Chen.  Phase-field models for microstructure evolution.  Annual Review of Materials Research, 32(Volume 32, 2002):113–140, 2002.
  * [2] A. Karma.  Phase field methods.  In K.H. Jürgen Buschow, Robert W. Cahn, Merton C. Flemings, Bernhard Ilschner, Edward J. Kramer, Subhash Mahajan, and Patrick Veyssière, editors, Encyclopedia of Materials: Science and Technology, pages 6873–6886. Elsevier, Oxford, 2001.
  * [3] Ingo Steinbach.  Phase-field models in materials science.  Modelling and Simulation in Materials Science and Engineering, 17(7):073001, jul 2009.
  * [4] Nikolas Provatas, Nigel Goldenfeld, and Jonathan Dantzig.  Adaptive mesh refinement computation of solidification microstructures using dynamic data structures.  Journal of Computational Physics, 148(1):265–290, 1999.
  * [5] Alain Karma and Wouter-Jan Rappel.  Phase-field method for computationally efficient modeling of solidification with arbitrary interface kinetics.  Phys. Rev. E, 53:R3017–R3020, Apr 1996.
  * [6] Alain Karma and Wouter-Jan Rappel.  Quantitative phase-field modeling of dendritic growth in two and three dimensions.  Phys. Rev. E, 57:4323–4349, Apr 1998.
  * [7] Alain Karma.  Phase-field formulation for quantitative modeling of alloy solidification.  Phys. Rev. Lett., 87:115701, Aug 2001.
  * [8] J.D Harrison and C Wagner.  The attack of solid alloys by liquid metals and salt melts.  Acta Metall., 7(11):722–735, 1959.
  * [9] Jonah Erlebacher, Michael J Aziz, Alain Karma, Nikolay Dimitrov, and Karl Sieradzki.  Evolution of nanoporosity in dealloying.  Nature, 410(6827):450–453, 2001.
  * [10] Ian McCue, Alain Karma, and Jonah Erlebacher.  Pattern formation during electrochemical and liquid metal dealloying.  MRS Bulletin, 43(1):27–34, 2018.
  * [11] Ian McCue, Bernard Gaskey, Pierre-Antoine Geslin, Alain Karma, and Jonah Erlebacher.  Kinetics and morphological evolution of liquid metal dealloying.  Acta Materialia, 115:10–23, 2016.
  * [12] Takeshi Wada, Pierre-Antoine Geslin, Daixu Wei, and Hidemi Kato.  Partial liquid metal dealloying to synthesize nickel-containing porous and composite ferrous and high-entropy alloys.  Communications Materials, 4(1):43, December 2023.
  * [13] Xiaoyang Liu, Arthur Ronne, Lin-Chieh Yu, Yang Liu, Mingyuan Ge, Cheng-Hung Lin, Bobby Layne, Phillip Halstenberg, Dmitry S Maltsev, Alexander S Ivanov, et al.  Formation of three-dimensional bicontinuous structures via molten salt dealloying studied in real-time by in situ synchrotron x-ray nano-tomography.  Nat. Commun., 12(1):1–12, 2021.
  * [14] P. Calderoni and C. Cabet.  23 - corrosion issues in molten salt reactor (msr) systems.  In Damien Féron, editor, Nuclear Corrosion Science and Engineering, Woodhead Publishing Series in Energy, pages 842–865. Woodhead Publishing, 2012.
  * [15] Robin Roper, Megan Harkema, Piyush Sabharwall, Catherine Riddle, Brandon Chisholm, Brandon Day, and Paul Marotta.  Molten salt for advanced energy applications: A review.  Annals of Nuclear Energy, 169:108924, 2022.
  * [16] Pierre-Antoine Geslin, Ian Mccue, Bernard Gaskey, Jonah Erlebacher, and Alain Karma.  Topology-generating interfacial pattern formation during liquid metal dealloying.  Nature communications, 6:8887, 11 2015.
  * [17] Yu chen Karen Chen-Wiegart, Steve Wang, Yong S. Chu, Wenjun Liu, Ian McNulty, Peter W. Voorhees, and David C. Dunand.  Structural evolution of nanoporous gold during thermal coarsening.  Acta Materialia, 60(12):4972–4981, 2012.
  * [18] Longhai Lai, Bernard Gaskey, Alyssa Chuang, Jonah Erlebacher, and Alain Karma.  Topological control of liquid-metal-dealloyed structures.  Nature Communications, 13(1):2918, 2022.
  * [19] Longhai Lai, Pierre-Antoine Geslin, and Alain Karma.  Microstructural pattern formation during liquid metal dealloying: Phase-field simulations and theoretical analyses.  Phys. Rev. Mater., 6:093803, Sep 2022.
  * [20] Nathan Bieberdorf, Mark Asta, and Laurent Capolungo.  Grain boundary effects in high-temperature liquid-metal dealloying: a multi-phase field study.  npj Computational Materials, 9, 07 2023.
  * [21] Lorenz Biegler, George Biros, Omar Ghattas, Matthias Heinkenschloss, David Keyes, Bani Mallick, Youssef Marzouk, Luis Tenorio, Bart van Bloemen Waanders, and Karen Willcox.  Large-scale inverse problems and quantification of uncertainty, 2010.
  * [22] Ralph C. Smith.  Uncertainty quantification - theory, implementation, and applications.  In Computational science and engineering, 2013.
  * [23] D. Galbally, K. Fidkowski, K. Willcox, and O. Ghattas.  Non-linear model reduction for uncertainty quantification in large-scale inverse problems.  International Journal for Numerical Methods in Engineering, 81(12):1581–1608, 2010.
  * [24] Vasileios Fountoulakis and Christopher Earls.  Duct heights inferred from radar sea clutter using proper orthogonal bases.  Radio Science, 51(10):1614–1626, 2016.
  * [25] Benjamin Peherstorfer, Karen Willcox, and Max Gunzburger.  Survey of multifidelity methods in uncertainty propagation, inference, and optimization.  Siam Review, 60(3):550–591, 2018.
  * [26] Shun Wang, Eric de Sturler, and Glaucio H. Paulino.  Large-scale topology optimization using preconditioned krylov subspace methods with recycling.  International Journal for Numerical Methods in Engineering, 69(12):2441–2468, 2007.
  * [27] Daniel White, Youngsoo Choi, and Jun Kudo.  A dual mesh method with adaptivity for stress-constrained topology optimization.  Structural and Multidisciplinary Optimization, 61, 02 2020.
  * [28] Alexander Forrester, Andras Sobester, and Andy Keane.  Engineering design via surrogate modelling: a practical guide.  John Wiley & Sons, 2008.
  * [29] G Berkooz, P Holmes, and J L Lumley.  The proper orthogonal decomposition in the analysis of turbulent flows.  Annual Review of Fluid Mechanics, 25(1):539–575, 1993.
  * [30] Gianluigi Rozza, D. Huynh, and A. Patera.  Reduced basis approximation and a posteriori error estimation for affinely parametrized elliptic coercive partial differential equations.  Archives of Computational Methods in Engineering, 15:1–47, 09 2007.
  * [31] Peter Benner, Serkan Gugercin, and Karen Willcox.  A survey of projection-based model reduction methods for parametric dynamical systems.  SIAM review, 57(4):483–531, 2015.
  * [32] Model reduction and approximation: theory and algorithms, 2017.
  * [33] Michael G. Safonov and Richard Y. Chiang.  A schur method for balanced model reduction.  1988 American Control Conference, pages 1036–1040, 1988.
  * [34] Giovanni Stabile and Gianluigi Rozza.  Finite volume POD-galerkin stabilised reduced order methods for the parametrised incompressible navier–stokes equations.  Computers & Fluids, 173:273–284, sep 2018.
  * [35] Traian Iliescu and Zhu Wang.  Variational multiscale proper orthogonal decomposition: Navier-stokes equations.  Numerical Methods for Partial Differential Equations, 30(2):641–663, 2014.
  * [36] Dylan Matthew Copeland, Siu Wun Cheung, Kevin Huynh, and Youngsoo Choi.  Reduced order models for lagrangian hydrodynamics.  Computer Methods in Applied Mechanics and Engineering, 388:114259, jan 2022.
  * [37] Siu Wun Cheung, Youngsoo Choi, Dylan Matthew Copeland, and Kevin Huynh.  Local lagrangian reduced-order modeling for rayleigh-taylor instability by solution manifold decomposition, 2022.
  * [38] Benjamin McLaughlin, Janet Peterson, and Ming Ye.  Stabilized reduced order models for the advection–diffusion–reaction equation using operator splitting.  Computers & Mathematics with Applications, 71(11):2407–2420, 2016.  Proceedings of the conference on Advances in Scientific Computing and Applied Mathematics. A special issue in honor of Max Gunzburger’s 70th birthday.
  * [39] Philip Holmes.  Turbulence, coherent structures, dynamical systems and symmetry.  Cambridge university press, 2012.
  * [40] David Amsallem, Julien Cortial, Kevin Carlberg, and Charbel Farhat.  A method for interpolating on manifolds structural dynamics reduced-order models.  International journal for numerical methods in engineering, 80(9):1241–1258, 2009.
  * [41] Ioannis T Georgiou and Ira B Schwartz.  Dynamics of large scale coupled structural/mechanical systems: A singular perturbation/proper orthogonal decomposition approach.  SIAM Journal on Applied Mathematics, 59(4):1178–1207, 1999.
  * [42] Robert J Guyan.  Reduction of stiffness and mass matrices.  AIAA journal, 3(2):380–380, 1965.
  * [43] Athanasios C Antoulas.  Approximation of large-scale dynamical systems.  SIAM, 2005.
  * [44] Arthur J Krener.  Reduced order modeling of nonlinear control systems.  In Analysis and Design of Nonlinear Control Systems: In Honor of Alberto Isidori, pages 41–62. Springer, 2008.
  * [45] Maria Strazzullo, Francesco Ballarin, Renzo Mosetti, and Gianluigi Rozza.  Model reduction for parametrized optimal control problems in environmental marine sciences and engineering.  SIAM Journal on Scientific Computing, 40(4):B1055–B1079, 2018.
  * [46] Kookjin Lee and Kevin T. Carlberg.  Model reduction of dynamical systems on nonlinear manifolds using deep convolutional autoencoders.  Journal of Computational Physics, 404:108973, 2020.
  * [47] J. Nathan Kutz.  Deep learning in fluid dynamics.  Journal of Fluid Mechanics, 814:1–4, 2017.
  * [48] Christophe Bonneville, Youngsoo Choi, Debojyoti Ghosh, and Jonathan L. Belof.  Gplasdi: Gaussian process-based interpretable latent space dynamics identification through deep autoencoder.  Computer Methods in Applied Mechanics and Engineering, 418:116535, 2024.
  * [49] Christophe Bonneville, Xiaolong He, April Tran, Jun Sur Park, William Fries, Daniel A. Messenger, Siu Wun Cheung, Yeonjong Shin, David M. Bortz, Debojyoti Ghosh, Jiun-Shyan Chen, Jonathan Belof, and Youngsoo Choi.  A comprehensive review of latent space dynamics identification algorithms for intrusive and non-intrusive reduced-order-modeling, 2024.
  * [50] Nikola B. Kovachki, Zongyi Li, Burigede Liu, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew M. Stuart, and Anima Anandkumar.  Neural operator: Learning maps between function spaces.  CoRR, abs/2108.08481, 2021.
  * [51] Zongyi Li, Nikola Kovachki, Kamyar Azizzadenesheli, Burigede Liu, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar.  Fourier neural operator for parametric partial differential equations, 2021.
  * [52] Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis.  Learning nonlinear operators via deeponet based on the universal approximation theorem of operators.  Nature Machine Intelligence, 3(3):218–229, March 2021.
  * [53] Qianying Cao, Somdatta Goswami, and George Em Karniadakis.  Lno: Laplace neural operator for solving differential equations, 2023.
  * [54] Bogdan Raonić, Roberto Molinaro, Tim De Ryck, Tobias Rohner, Francesca Bartolucci, Rima Alaifari, Siddhartha Mishra, and Emmanuel de Bézenac.  Convolutional neural operators for robust and accurate learning of pdes, 2023.
  * [55] Abouzar Choubineh, Jie Chen, David A. Wood, Frans Coenen, and Fei Ma.  Fourier neural operator for fluid flow in small-shape 2d simulated porous media dataset.  Algorithms, 16(1), 2023.
  * [56] Zhijie Li, Wenhui Peng, Yuan Zelong, and Jianchun Wang.  Fourier neural operator approach to large eddy simulation of three-dimensional turbulence.  Theoretical and Applied Mechanics Letters, 10 2022.
  * [57] Nattawat Pornthisan and Stefano Markidis.  Fast electromagnetic field pattern calculation with fourier neural operators.  In Jiří Mikyška, Clélia de Mulatier, Maciej Paszynski, Valeria V. Krzhizhanovskaya, Jack J. Dongarra, and Peter M.A. Sloot, editors, Computational Science – ICCS 2023, pages 247–255, Cham, 2023. Springer Nature Switzerland.
  * [58] Shengze Cai, Zhicheng Wang, Lu Lu, Tamer A. Zaki, and George Em Karniadakis.  Deepm&mnet: Inferring the electroconvection multiphysics fields based on operator approximation by neural networks.  Journal of Computational Physics, 436:110296, 2021.
  * [59] Somdatta Goswami, Minglang Yin, Yue Yu, and George Em Karniadakis.  A physics-informed variational deeponet for predicting crack path in quasi-brittle materials.  Computer Methods in Applied Mechanics and Engineering, 391:114587, 2022.
  * [60] Jaideep Pathak, Shashank Subramanian, Peter Harrington, Sanjeev Raja, Ashesh Chattopadhyay, Morteza Mardani, Thorsten Kurth, David Hall, Zongyi Li, Kamyar Azizzadenesheli, Pedram Hassanzadeh, Karthik Kashinath, and Animashree Anandkumar.  Fourcastnet: A global data-driven high-resolution weather model using adaptive fourier neural operators, 2022.
  * [61] Thorsten Kurth, Shashank Subramanian, Peter Harrington, Jaideep Pathak, Morteza Mardani, David Hall, Andrea Miele, Karthik Kashinath, and Anima Anandkumar.  Fourcastnet: Accelerating global high-resolution weather forecasting using adaptive fourier neural operators.  In Proceedings of the Platform for Advanced Scientific Computing Conference, PASC ’23, New York, NY, USA, 2023. Association for Computing Machinery.
  * [62] Vivek Oommen, Khemraj Shukla, Somdatta Goswami, Remi Dingreville, and George Karniadakis.  Learning two-phase microstructure evolution using neural operators and autoencoder architectures.  npj Computational Materials, 8, 09 2022.
  * [63] Gege Wen, Zongyi Li, Kamyar Azizzadenesheli, Anima Anandkumar, and Sally M. Benson.  U-fno – an enhanced fourier neural operator-based deep-learning model for multiphase flow, 2022.
  * [64] Vignesh Gopakumar, Stanislas Pamela, Lorenzo Zanisi, Zongyi Li, Anima Anandkumar, and MAST Team.  Fourier neural operator for plasma modelling, 2023.
  * [65] Olaf Ronneberger, Philipp Fischer, and Thomas Brox.  U-net: Convolutional networks for biomedical image segmentation, 2015.
  * [66] Alexey Dosovitskiy, Lucas Beyer, Alexander Kolesnikov, Dirk Weissenborn, Xiaohua Zhai, Thomas Unterthiner, Mostafa Dehghani, Matthias Minderer, Georg Heigold, Sylvain Gelly, Jakob Uszkoreit, and Neil Houlsby.  An image is worth 16x16 words: Transformers for image recognition at scale.  In International Conference on Learning Representations, 2021.
  * [67] Waleed Diab and Mohammed Al-Kobaisi.  U-deeponet: U-net enhanced deep operator network for geologic carbon sequestration, 2023.
  * [68] Zhijie Li, Wenhui Peng, Zelong Yuan, and Jianchun Wang.  Long-term predictions of turbulence by implicit U-Net enhanced Fourier neural operator.  Physics of Fluids, 35(7):075145, 07 2023.
  * [69] Oded Ovadia, Adar Kahana, Panos Stinis, Eli Turkel, and George Em Karniadakis.  Vito: Vision transformer-operator, 2023.
  * [70] Oded Ovadia, Vivek Oommen, Adar Kahana, Ahmad Peyvan, Eli Turkel, and George Em Karniadakis.  Real-time inference and extrapolation via a diffusion-inspired temporal transformer operator (ditto), 2023.
  * [71] John Guibas, Morteza Mardani, Zongyi Li, Andrew Tao, Anima Anandkumar, and Bryan Catanzaro.  Adaptive fourier neural operators: Efficient token mixers for transformers, 2022.
  * [72] Zongyi Li, Miguel Liu-Schiaffini, Nikola Kovachki, Burigede Liu, Kamyar Azizzadenesheli, Kaushik Bhattacharya, Andrew Stuart, and Anima Anandkumar.  Learning dissipative dynamics in chaotic systems, 2022.
  * [73] Vivek Oommen, Khemraj Shukla, Saaketh Desai, Remi Dingreville, and George Em Karniadakis.  Rethinking materials simulations: Blending direct numerical simulations with neural operators, 2023.
  * [74] Ingo Steinbach and F Pezzolla.  A generalized field method for multiphase transformations using interface fields.  Physica D, 134(4):385–393, 1999.
  * [75] Vittorio E Badalassi, Hector D Ceniceros, and Sanjoy Banerjee.  Computation of multiphase systems with phase field models.  J. Comput. Phys., 190(2):371–397, 2003.
  * [76] Kai Han, Yunhe Wang, Hanting Chen, Xinghao Chen, Jianyuan Guo, Zhenhua Liu, Yehui Tang, An Xiao, Chunjing Xu, Yixing Xu, Zhaohui Yang, Yiman Zhang, and Dacheng Tao.  A survey on vision transformer.  IEEE Transactions on Pattern Analysis and Machine Intelligence, 45(1):87–110, January 2023.
  * [77] Salman Khan, Muzammal Naseer, Munawar Hayat, Syed Waqas Zamir, Fahad Shahbaz Khan, and Mubarak Shah.  Transformers in vision: A survey.  ACM Computing Surveys, 54(10s):1–41, January 2022.
  * [78] Zhongkai Hao, Chang Su, Songming Liu, Julius Berner, Chengyang Ying, Hang Su, Anima Anandkumar, Jian Song, and Jun Zhu.  Dpot: Auto-regressive denoising operator transformer for large-scale pde pre-training, 2024.
  * [79] Alan John Varghese, Aniruddha Bora, Mengjia Xu, and George Em Karniadakis.  Transformerg2g: Adaptive time-stepping for learning temporal graph embeddings using transformers.  Neural Networks, 172:106086, 2024.
  * [80] Zijie Li, Kazem Meidani, and Amir Barati Farimani.  Transformer for partial differential equations’ operator learning, 2023.
  * [81] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Lukasz Kaiser, and Illia Polosukhin.  Attention is all you need, 2023.
  * [82] Lucas Beyer, Pavel Izmailov, Alexander Kolesnikov, Mathilde Caron, Simon Kornblith, Xiaohua Zhai, Matthias Minderer, Michael Tschannen, Ibrahim Alabdulmohsin, and Filip Pavetic.  Flexivit: One model for all patch sizes, 2023.
  * [83] Samad Vakili, Ingo Steinbach, and Fathollah Varnik.  On the numerical evaluation of local curvature for diffuse interface models of microstructure evolution.  Procedia Computer Science, 108:1852–1862, 2017.  International Conference on Computational Science, ICCS 2017, 12-14 June 2017, Zurich, Switzerland.
  * [84] Anh Tran, Dehao Liu, Hoang Tran, and Yan Wang.  Quantifying uncertainty in the process-structure relationship for al–cu solidification.  Modelling and Simulation in Materials Science and Engineering, 27(6):064005, jun 2019.
  * [85] Diederik P. Kingma and Jimmy Ba.  Adam: A method for stochastic optimization, 2017.
  * [86] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zach DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.  Pytorch: An imperative style, high-performance deep learning library, 2019.
  * [87] Elizabeth Herman, James A. Stewart, and Rémi Dingreville.  A data-driven surrogate model to rapidly predict microstructure morphology during physical vapor deposition.  Applied Mathematical Modelling, 88:589–603, 2020.