Open AccessArticle

Physics-Informed Neural Networks with Periodic Activation Functions for Solute Transport in Heterogeneous Porous Media

Salah A. Faroughi

^1,*

Ramin Soltanmohammadi

Pingki Datta

¹,

Seyed Kourosh Mahjour

¹ and

Shirko Faroughi

Geo-Intelligence Laboratory, Ingram School of Engineering, Texas State University, San Marcos, TX 78666, USA

Department of Mechanical Engineering, School of Engineering, Urmia University of Technology, Urmia 57561-51818, Iran

Author to whom correspondence should be addressed.

Mathematics 2024, 12(1), 63; https://doi.org/10.3390/math12010063

Submission received: 9 November 2023 / Revised: 28 November 2023 / Accepted: 30 November 2023 / Published: 24 December 2023

(This article belongs to the Special Issue Advances in Computational Fluid Dynamics)

Download

Browse Figures

Versions Notes

Abstract

Simulating solute transport in heterogeneous porous media poses computational challenges due to the high-resolution meshing required for traditional solvers. To overcome these challenges, this study explores a mesh-free method based on deep learning to accelerate solute transport simulation. We employ Physics-informed Neural Networks (PiNN) with a periodic activation function to solve solute transport problems in both homogeneous and heterogeneous porous media governed by the advection-dispersion equation. Unlike traditional neural networks that rely on large training datasets, PiNNs use strong-form mathematical models to constrain the network in the training phase and simultaneously solve for multiple dependent or independent field variables, such as pressure and solute concentration fields. To demonstrate the effectiveness of using PiNNs with a periodic activation function to resolve solute transport in porous media, we construct PiNNs using two activation functions, sin and tanh, for seven case studies, including 1D and 2D scenarios. The accuracy of the PiNNs’ predictions is then evaluated using absolute point error and mean square error metrics and compared to the ground truth solutions obtained analytically or numerically. Our results demonstrate that the PiNN with sin activation function, compared to tanh activation function, is up to two orders of magnitude more accurate and up to two times faster to train, especially in heterogeneous porous media. Moreover, PiNN’s simultaneous predictions of pressure and concentration fields can reduce computational expenses in terms of inference time by three orders of magnitude compared to FEM simulations for two-dimensional cases.

Keywords:

physics-informed neural networks; solute transport; heterogeneous porous media; advection-dispersion equation; deep learning; scientific computing

MSC:

76505

1. Introduction

Solute transport in porous media is crucial in many environmental and reservoir engineering applications, including risk and safety assessment of groundwater contamination [1], secondary recovery processes in petroleum reservoirs [2], geological storage of radioactive waste [3], hydrogen [4], and carbon dioxide [5], just to name a few [6]. It is a complicated process due to the varied temporal and spatial transport characteristics of the flow network [7]. These characteristics can be divided into flowing areas influenced primarily by advection and stagnant regions governed mainly by dispersion [8]. The process is considerably more complicated in heterogeneous porous media because all available pore space does not contribute to the flow uniformly; certain locations are dead ends, or transport is extremely slow due to inadequate connection to the main flow channels [9]. Despite a conceptual understanding of the significance of multi-component solute transport mechanisms in porous media in general, attempts to conduct reliable experimental investigations are hindered by several barriers, resulting in a dearth of published information. The major causes are the technical difficulty in creating flow conditions and pore space connections under regulated laboratory or in situ test circumstances [10]. Hence, a system of field equations, which are numerically represented by partial differential equations (PDEs), is practically the only tool to predict the behavior of such non-linear flows in interconnected systems. In solute transport studies, traditional solvers, such as finite element (FEM), finite difference (FDM), and finite volume methods (FVM), have been extensively used to solve the governing PDEs. For example, Zhang et al. [11] suggested a 1D model based on FEM for multi-component solute transport in saturated soil. Bagalkot and Suresh Kumar [12] presented a 1D numerical evaluation of the FDM for multi-species radionuclide transport in a single horizontally coupled fracture-matrix system. Mostaghimi et al. [13] and Maheshwari et al. [14] investigated the effect of pore structure heterogeneity on reactive transport using a 3D pore scale model based on FVM. Although most large-scale calculations still use FEM, FDM, or FVM, they require high-resolution meshing to describe the geometry of the domain being represented. Because the size of an elementary mesh is normally too large, discrete domain descriptions may also need the availability of efficient transport equations at the mesh scale and the consideration of unresolved subgrid-scale effects [15]. Hence, new methodologies and algorithms, especially mesh-free methods based on deep learning, are required to accelerate numerical simulations of solute transport through porous media. Deep Learning (DL) methods have the potential to offer mesh-free solvers and address some of the aforementioned challenges. Although most applications use DL (e.g., neural networks) to solve a lack of effective data modeling processes, explore vast design domains, and identify multidimensional connections [16,17], there is growing interest in using neural networks to solve PDEs [18,19,20,21]. In general, there are several main neural network frameworks to augment scientific computing: Physics-guided Neural Networks (PgNN), Physics-informed Neural Networks (PiNN), Physics-encoded Neural Networks (PeNN), and Neural Operators (NOs). The readers are referred to a recent review by Faroughi et al. [22], where different neural network frameworks are compared head-to-head, and their challenges and limitations are thoroughly discussed. In this study, we adopt PiNN because of its straightforward mechanism to integrate the underlying physics compared to other approaches. By combining a loss function composed of the residuals of the physics equations, initial conditions, and boundary restrictions, PiNN-based models adhere to the physical laws. They employ automated differentiation to differentiate the output of neural networks with respect to their input (that is, spatio-temporal coordinates and model parameters) [23]. The network can estimate the solution with high precision by minimizing the loss function [24]. Therefore, PiNN provides the foundation for a mesh-free solver that incorporates long-standing advances in mathematical physics into DL [18].

Raissi et al. [25] developed PiNN as a new computing paradigm for forward and inverse modeling in a series of studies [18,25,26]. Raissi et al. [26] developed a PiNN framework, dubbed hidden fluid mechanics (HFM), to encode the physical laws governing fluid flows, i.e., the Navier-Stokes equations. They leveraged underlying conservation laws to derive hidden quantities of objective functions such as velocity and pressure fields from spatiotemporal visualizations of a passive scalar concentration in arbitrarily complex domains. Their technique accurately predicted 2D and 3D pressure and velocity fields in benchmark problems inspired by real-world applications. Since then, PiNN and its different variants have been applied aggressively to different fields, including porous media flows. Almajid and Abu-Al-Saud [27] employed the PiNN framework to solve both the forward and inverse problems of the Buckley-Leverett PDE equation representing two-phase flow in porous media. To test their implementation, they applied the classic problem of gas drainage through a porous water-filled medium. Several cases were examined to demonstrate the significance of the connectivity between observable data and PiNNs for various parameter spaces. According to their results, PiNNs are capable of capturing the solution’s broad trend even without observed data, but their precision and accuracy improve significantly with observed data. Hanna et al. [28] employed PiNN to simulate one-dimensional (1D) and two-dimensional (2D) two-phase flow in porous media based on a new residual-based adaptive algorithm. Applying the PDE residual to build a probability density function from which additional collocation points are taken and added to the training set was fundamental to their work. The approach was applied individually to each PDE in the coupled system, considering the different collocation points for each PDE. Furthermore, the method was applied to enrich the points used to capture the initial and boundary conditions. They claimed that their approach yielded superior results with less generalization error than conventional PiNN with fixed collocation points. Haghighat et al. [29] introduced a PiNN method to solve coupled flow and deformation equations in porous media for single-phase and multiphase flow. Due to the problem’s dynamic nature, they reported a dimensionless form of the coupled governing equations for the optimizer. In addition, they presented a way for sequential training based on the stress-split algorithms of poromechanics. They demonstrated that sequential training based on stress-split performs well for a variety of problems, whereas the conventional strain-split algorithm exhibits instability comparable to that observed for FEM solvers. He et al. [30] extended a PiNN-based parameter estimation method to integrate multiphysics measurement. They studied a subsurface transport problem with sparse conductivity, a hydraulic head, and solute concentration. In their methodology, they employed the Darcy and advection-dispersion equations in conjunction with the data to train deep neural networks, reflecting space-dependent conductivity, head, and concentration fields. They proved that the proposed PiNN method considerably increased the accuracy of parameter and state estimates for sparse data compared to conventional deep neural networks trained with data alone. He and Tartakovsky [31] also suggested a discretization-free technique based on PiNN for solving coupled advection-dispersion equations and the Darcy flow equation with space-dependent hydraulic conductivity. They used PiNN for 1D and 2D forward advection-dispersion equations and compared its performance for various Peclet numbers (Pe) with analytical and numerical solutions. They found that PiNN was accurate with errors of less than 1% and outperformed other discretization-based approaches for large Pe. In addition, they proved that PiNN remained accurate for the backward advection-dispersion equations, with relative errors below 5% in the majority of instances. In another study, Vadyala et al. [32] implemented PiNN with the use of a machine learning framework such that it can be employed in reduced-order models to reduce the epistemic (model-form) ambiguity associated with the advection equation. They demonstrated that PiNN provided an accurate and consistent approximation with PDEs. Furthermore, they showed that PiNN could transform the physics simulation field by enabling real-time physics simulation and geometry optimization on large supercomputers.

Although there has been an increase in the application of PiNN to fluid flow studies in porous media [33,34,35,36], the application of PiNN to the advection–dispersion equation in heterogeneous porous media has received less attention due to the difficulty in predicting the velocity field imposed by the variation in the permeability field. The main objective of this study is to develop a PiNN model capable of accurately predicting transient solute transport in heterogeneous porous media. Due to the presence of the second-order derivative in the governing equations, selecting the activation function for PiNN becomes a critical step [37]. To this end, a PiNN is constructed using periodic activation functions, and its convergence rate and accuracy are compared against a PiNN using tanh to model transport phenomena in porous media. Although it is theoretically believed that neural networks with periodic activation functions, e.g., sinusoids, may be harder to train, it is shown in several cases that periodicity is intuitively beneficial and can learn faster and better than other activation functions [37,38]. In this study, using 1D and 2D test cases, we also show that a PiNN with a sine activation function provides more accurate predictions for solute transport in heterogeneous porous media while reducing the training time.

This paper is structured as follows: In Section 2, the underlying physics for solute transport in a porous medium is presented. Section 3 discusses PiNN’s algorithm with periodic activation functions to predict solute transport. In Section 4, the PiNN is deployed to several one-dimensional (1D) and two-dimensional (2D) case studies, and its predictions are examined by comparing them with analytical and/or FEM solutions. Finally, in Section 5, we summarize the main conclusions of this work.

2. Underlying Physics

Solute transport in porous media is governed by three main mechanisms: advection, diffusion, and dispersion [39,40,41]. The mathematical formulation of the advection-diffusion-dispersion equation is explained briefly here in terms of implementing the flow and transport factors that affect the transport process. The mass conservation equation serves as the basis for the solute transport equation [42,43] that reads as,

\nabla \cdot J = - \frac{\partial}{\partial t} (ϕ C),

(1)

where

J

[ML

^{- 2}

^{- 1}

] is the solute mass flux,

ϕ

is the porosity, and C [ML

^{- 3}

] is the solute concentration. The flux

J

contains two mechanisms,

J = - ϕ D \nabla C + u (ϕ C),

(2)

where the first term refers to the diffusion flux, and the second term refers to the advection flux [44,45]. In Equation (2), D [L

^{2}

^{- 1}

] is the molecular diffusion coefficient (

D = D_{0}

), and

u

[LT

^{- 1}

] is the velocity vector of the pore fluid flow computed using Darcy’s law,

u = q / ϕ = - \frac{k}{μ ϕ} (\nabla P + ρ g),

(3)

where

q

[LT

^{- 1}

] is the Darcy’s flux, P [ML

^{- 1}

^{- 2}

] refers to the pressure field,

g

[LT

^{- 2}

] is the gravity, k [L

^{2}

] is the permeability, and

μ

[ML

^{- 1}

^{- 1}

] and

ρ

[ML

^{- 3}]

are the fluid viscosity and density, respectively. To include the impact of dispersion, the diffusion coefficient in Equation (2) can be changed to

D = D_{0} + α U,

(4)

where

α

is known as the dynamic dispersivity, D is now called the hydrodynamic dispersion coefficient, and U is the magnitude of fluid velocity defined as

U = | u |

. In this work, we assume that longitudinal and transverse dispersion are identical; however, we allow U to vary spatially for a heterogeneous porous medium where

u

varies due to permeability changes [46,47]. The reason for this assumption is to focus on the impact of permeability and how PiNN can handle it when supplied with different activation functions. Substituting Equation (2) into Equation (1) leads to

\nabla \cdot [ϕ D \nabla C] - \nabla \cdot [u (ϕ C)] = \frac{\partial (ϕ C)}{\partial t},

(5)

which can be rewritten as,

\nabla \cdot [ϕ D \nabla C] - (\nabla \cdot u) [(ϕ C)] - u \cdot \nabla (ϕ C) = \frac{\partial (ϕ C)}{\partial t},

(6)

that incorporates diffusion, dispersion, and advection transports and serves as the foundation of solute transport in porous media [48]. If the solute transport process has no effect on fluid density (i.e., incompressible flow,

\nabla \cdot u = 0

), the velocity field can be calculated independently of the solute concentration using Darcy’s law, and Equation (6) reduces to,

\nabla \cdot [ϕ D \nabla C] - u \cdot \nabla (ϕ C) = \frac{\partial (ϕ C)}{\partial t} .

(7)

The flow incompressibility condition generates an equation for computing the pressure field as,

\nabla \cdot u = \nabla \cdot (- \frac{k}{μ ϕ} \nabla P + ρ g) = 0,

(8)

that reads as,

\frac{\partial}{\partial x} [ζ (x, y) \frac{\partial P}{\partial x}] + \frac{\partial}{\partial y} [ζ (x, y) \frac{\partial P}{\partial y}] = 0,

(9)

for two-dimensional (2D) flow problems in x and y direction assuming

g_{x} = g_{y} = 0

, where

ζ (x, y) = \frac{k (x, y)}{μ ϕ} .

(10)

Equations (7) and (9) together form the governing equations for the concentration and pressure fields in a porous medium with a heterogeneous permeability distribution [48].

3. Methodology

In this study, Physics-informed Neural Networks (PiNN) [17,18] are leveraged to resolve the solute transport in porous media, which is governed by a strong mathematical form consisting of the pressure and advection-dispersion equations as well as the relevant initial and boundary conditions. Unlike Physics-guided Neural Networks (PgNNs), PiNN is a solution learning method that does not require labeled datasets [22]. In PiNNs, the underlying physics is included outside of the neural network architecture to constrain the model during training and ensure that the outputs follow given physical laws. The most common method to emulate this process is through a weakly imposed penalty loss that penalizes the network when it does not follow the physical constraints [25,49,50].

A schematic representation of our PiNN architecture is illustrated in Figure 1. As illustrated, a PiNN consists of three elements: a neural network (NN), an automatic differentiation (AD) layer, and a feedback mechanism informed by physics [51]. A neural network is first built in order to digest the spatiotemporal features (i.e.,

x

and t) as input parameters and approximate the solution,

s (x, t)

, of a phenomenon described by PDEs with known boundary values. Assume that

N_{L} (x) : R^{d_{i n}}

→

R^{d_{o u t}}

be a L-layer neural network, or a (

L - 1

)-hidden layer neural network, with

N_{l}

neurons in the

l^{t h}

layer (

N_{0}

d_{i n}

N_{L}

d_{o u t}

) [52]. Each layer in the NN is associated with its own weight matrix

W^{l} \in R^{N_{l} \times N_{l - 1}}

and bias vector

b^{l} \in R^{N_{l}}

. The NN is recursively defined using an element-wise nonlinear activation function,

σ

, as follows, [53],

\begin{matrix} I n p u t l a y e r : N_{0} (x) & = x \in R^{d_{i n}}, \end{matrix}

(11)

\begin{matrix} H i d d e n l a y e r s : N_{l} (x) & = σ (W^{l} N_{l - 1} (x) + b^{l}) \in R^{N_{l}}, f o r 1 \leq l \leq L - 1, \end{matrix}

(12)

\begin{matrix} O u t p u t l a y e r : N_{L} (x) & = W^{l} N_{L - 1} (x) + b^{l} \in R^{d_{o u t}}, \end{matrix}

(13)

where

σ

is commonly specified as the logistic sigmoid

σ (x) = 1 / (1 - e x p (- x))

, hyperbolic tangent

σ (x) = t a n h (x)

, or the rectified linear unit,

σ (x) = m a x (x, 0)

[52]. In this work, we propose to use a sine (i.e., periodic) activation function,

σ (W^{l} N_{l - 1} (x) + b^{l}) = sin (W^{l} N_{l - 1} (x) + b^{l}),

(14)

and compare its performance against tanh, as one of the most used activation functions for PiNNs [37], in simulating solute transport in porous media. Periodic activation functions are ideally suited for representing complex physical signals and their derivatives while converging faster than baseline architectures [38,54]. This superior behavior, in comparison to tanh, is attributed to (i) the presence of simple sum and diff tasks for which a sine function can outperform hyperbolic tangent function [38], and (ii) the derivative of a network with sine activation function behaves like the network itself (i.e., the derivative of the sine is a cosine that is a phase-shifted sine) that potentially speeds up the learning process [54], especially in the solute transport problem that contains the second order derivatives in the governing equations.

The NN outputs are then fed into the AD layer, which is the central property of PiNNs. AD is employed to assess the derivatives of the network outputs with respect to the network inputs [23,52,55,56]. Consider a strong mathematical form with a PDE specified in the domain

Ω

and parameterized by

λ

f (x, t : Δ s; λ) = 0, x \in Ω,

(15)

and boundary conditions (e.g., Dirichlet, Neumann, Robin boundary condition, etc.) specified on the boundary of the domain

\partial Ω

ψ (s (x, t)) = 0, x \in \partial Ω,

(16)

where

s

is the unknown solution, and

Δ

represents the linear or nonlinear differential operator (e.g.,

\frac{\partial}{\partial t}

\frac{\partial}{\partial x}

\frac{\partial^{2}}{\partial x^{2}}

, etc.). The initial condition can be considered as a kind of Dirichlet boundary condition in the spatiotemporal domain. AD is used to apply the

Δ

and

ψ

operators to the neural network with respect to inputs, i.e.,

x

and t, and generate the required terms in the loss function in order to optimize the PDE solution,

s

. Lastly, a feedback mechanism is constructed to minimize the loss terms through optimizations. The total loss term is defined as [56,57,58],

L = w_{i} L_{I C} + w_{b} L_{B C} + w_{d} L_{D a t a} + w_{p} L_{P D E},

(17)

where

w_{i}

w_{b}

w_{d}

, and

w_{p}

, are referred to as the weights for the loss due to the initial conditions, boundary conditions, labeled data, if any, and PDEs, respectively. The individual loss terms are computed as,

\begin{matrix} L_{I C} & = \frac{1}{N_{i}} \sum_{i = 1}^{N_{i}} {(\hat{s} |_{Ω, t_{0}} - s |_{Ω, t_{0}})}^{2}, \end{matrix}

(18)

\begin{matrix} L_{B C} & = \frac{1}{N_{b}} \sum_{i = 1}^{N_{b}} {((\partial_{n} \hat{s} |_{\partial Ω} - \partial_{n} s |_{\partial Ω}) - (\hat{s} |_{\partial Ω} - s |_{\partial Ω}))}^{2}, \end{matrix}

(19)

\begin{matrix} L_{D a t a} & = \frac{1}{N_{d}} \sum_{i = 1}^{N_{d}} {(\hat{s} |_{Ω} - s |_{D a t a})}^{2}, \end{matrix}

(20)

\begin{matrix} L_{P D E} & = \frac{1}{N_{p}} \sum_{i = 1}^{N_{p}} {(f (s, \partial_{t} s, \partial_{x} s, \dots, λ))}^{2}, \end{matrix}

(21)

where

t_{0}

is the initial time, s is the neural network prediction,

\hat{s}

is the ground truth solution,

N_{i}

N_{b}

N_{d}

, and

N_{p}

are, respectively, the number of spatiotemporal points representing the initial conditions, boundary conditions, labeled data, and collocation points (i.e., spatiotemporal points within the domain where the neural network prediction,

s (x, t)

is checked against the constraints of PDEs). PiNN with weights

θ

can then be trained by optimizing,

θ^{'} = arg min_{θ} \sum_{1}^{N^{*}} L (θ),

(22)

leading to the PiNN with weights

θ^{'}

that generates reasonable predictions.

N^{*}

represents the total number of input-output pairs in the training process. PiNN training is more difficult than PgNN training, because PiNNs are composed of sophisticated non-convex and multi-objective loss functions, which may cause instability during optimization [59,60]. The selection of weights for the loss terms is ad-hoc and problem (PDE) dependent. The weights are adjusted using trials and errors in the training phase to reach the minimum minimization error or mitigate the instability of the solution [22]. In this study, we start with identical weights for all loss terms, and as needed (especially for heterogeneous cases where solutions encounter instability issues due to variation in the permeability field), we reduce the weight for PDE loss terms,

ω_{p}

, to fully respect the boundary conditions and mitigate the instability in the solution. Also, we set

w_{d} = 0

as we do not provide labeled data to PiNN.

Finally, a gradient-based optimizer such as Adam method [61] and Limited-memory Broyden–Fletcher Goldfarb–Shanno with box constraints, L-BFGS-B [62], should be used to minimize the loss function. It is found that the PiNN predictions strongly depend on which of these two algorithms is selected [53], and several studies suggested a two-step optimization algorithm that starts with the Adam method for a prescribed number of iterations, and then continues with the L-BFGS-B method until convergence [30,31]. In this work, L-BFGS-B, as a second-order, optimizer is used. It is observed that L-BFGS-B finds a satisfactory solution for smooth PDEs faster than Adam and with fewer iterations [63,64,65]. This is aligned with the goal of this work, which compares PiNN’s predictions when trained using different activation functions, without adding the complexity of alternating between optimizer schemes. Using L-BFGS-B, at each iteration, the loss is checked against

L < ϵ

, where

ϵ

is a specified tolerance. If the condition is not met, error backpropagation is implemented to update the learnable parameters (

θ

and/or

λ

). The entire cycle is repeated for a given number of iterations until the PiNN model produces learnable parameters with a loss error less than

ϵ

[66].

The proposed PiNN is deployed to model solute transport in porous media by predicting pressure and concentration fields, i.e.,

s = (P, C)

. We select the PiNN architectures using an iterative random-search hyperparameter tuning [67] (e.g., selecting the number of layers, neurons per layer, collocation points, and weights randomly within a specified range) for different case studies in this work due to differences in spatial dimension, boundary conditions, and permeability field. We ensure that the selected PiNN has enough layers and width to estimate each target field. There are two critical normalization steps in the implementation of PiNNs that can be followed to ensure faster convergence to the correct solution [68], similar to the common practice of scaling and dimensionless analysis in other mesh-based computational techniques, e.g., FVM, FEM, etc. The first step is to map the network input and output variables to the interval

[0, 1] \in R

. The second step is to scale the pressure and concentration equations such that all terms are of the same order. In this work, we use the absolute point error (APE) and the mean squared error (MSE) defined as,

\begin{matrix} A P E & = | s - \hat{s} |, \end{matrix}

(23)

\begin{matrix} M S E & = \frac{1}{N} Σ_{i = 1}^{N} {(s - \hat{s})}^{2}, \end{matrix}

(24)

as statistical measurands to assess the accuracy of PiNNs’ predictions, s, for a certain time-step against the ground truth solutions,

\hat{s}

, obtained analytically or numerically using the finite element method. Here, N is the total number of inference points (resembling the spatiotemporal mesh) to generate the solutions for comparison.

4. Computational Experiments

In this section, the proposed PiNN model is applied to solve the 1D and 2D solute transport phenomena under different conditions (i.e., a total of seven computational experiments with homogeneous and heterogeneous domains). We compare the PiNN models with sin and tanh activation functions and validate them against analytical and/or FEM solutions to examine their capability and accuracy. The comparisons are presented in terms of accuracy and training time.

4.1. Case 1: 1D Solute Transport with Constant Velocity

This test case explores solute transport in a 1D domain (

x \in [0, L]

) representing an isotropic porous medium with a constant velocity field. As the steady state velocity field is known, the governing equation, Equation (7), reduced to,

\frac{\partial C}{\partial t} + u_{x} \frac{\partial C}{\partial x} = \frac{\partial}{\partial x} (D_{x} \frac{\partial C}{\partial x}), 0 < x < L, 0 < t < t_{0},

(25)

where

L = 1

is the length of the domain,

u_{x} =

0.5 m/s is the steady state velocity field, and

D_{x}

= 0.02 m

^{2}

/s refers to the hydrodynamic dispersion coefficient in the x direction. For this case, we assume the following initial and boundary conditions,

\{\begin{matrix} C (x, t) = 0; t = 0, \\ C (x, t) = C_{0}; x = 0, \\ \frac{\partial}{\partial x} C (x, t) = 0; x = L . \end{matrix}

(26)

where

C_{0} = 1.0

kg/m

^{3}

is the injected concentration at

x = 0

(i.e., Dirichlet boundary condition). This 1D advection–dispersion can be solved analytically [69]. The analytical solution is,

C (x, t) = C_{0} (1 - 2 exp (\frac{x u_{x}}{2 D_{x}} - \frac{u_{x}^{2} t}{4 D_{x}}) \sum_{i = 1}^{\infty} \frac{β_{i} sin (\frac{β_{i} x}{L}) exp (- \frac{β_{i}^{2} D_{x} t}{L^{2}})}{β_{i}^{2} + {(\frac{u_{x} L}{2 D_{x}})}^{2} + \frac{u_{x} L}{2 D_{x}}}),

(27)

where

β_{i}

are the roots of,

β cot β + \frac{u_{x} L}{2 D_{x}} = 0 .

(28)

Equations (27) and (28) can be approximated using [42],

C (x, t) = C_{0} A (x, t), 0 < t \leq t_{0},

(29)

where

A (x, t)

is computed as,

\begin{matrix} A (x, t) & = \frac{1}{2} erfc (\frac{x - u_{x} t}{2 \sqrt{D_{x} t}}) + \frac{1}{2} exp (\frac{u_{x} x}{D_{x}}) erfc (\frac{x + u_{x} t}{2 \sqrt{D_{x} t}}) \\ + \frac{1}{2} (2 + \frac{u_{x} (2 L - x)}{D_{x}} + \frac{u_{x}^{2} t}{D_{x}}) exp (\frac{u_{x} L}{D_{x}}) erfc (\frac{2 L - x + u_{x} t}{2 \sqrt{D_{x} t}}) \\ - {(\frac{u_{x}^{2} t}{π D_{x}})}^{\frac{1}{2}} exp (\frac{u_{x} L}{D_{x}} - \frac{1}{4 D_{x} t} {(2 L - x + u_{x} t)}^{2}) . \end{matrix}

(30)

A PiNN with sin activation function and randomly distributed collocation points (with a uniform distribution scheme in the specified domain) is used to predict the solute transport in 1D domain described by Equation (25) with the initial and boundary conditions given in Equation (26). An iterative random-search hyperparameter tuning process is practiced to find the best architecture based on the MSE accuracy measures. This process yielded a network with four hidden layers with {32, 32, 16, 16} neurons per layer, 5000 spatiotemporal random points within

0 \leq x \leq 1

and

0 \leq t \leq 10

as the collocation points to enforce the PDE loss term, 5000 spatiotemporal random points to collectively represent the boundary and initial conditions informing the loss terms, and identical weights for all loss terms (

w_{i} = w_{b} = w_{p} = 1

). Another PiNN with similar architecture but using the tanh activation function is also trained and tuned for comparison purposes, see Table 1. For this PiNN model, the number of layers, neurons per layer, and collocation points are fixed, but the weights are allowed to change for better optimization. These assumptions make the comparisons fair in terms of both accuracy and training time.

Figure 2 represents the concentration field in the (

x, t

) domain predicted by PiNN using sin activation function. It also shows the PiNNs’ prediction compared with the ground truth (analytical solution given by Equation (30) at different simulation times. The comparison demonstrates a good agreement yielding

M S E = 1.15 \times 10^{- 6}

using sin activation function and

M S E = 1.21 \times 10^{- 6}

using tanh activation function for the entire spatiotemporal domain. As reported in Table 1, for this 1D case, the PiNN with tanh activation function approximates the solution with the same order of magnitude accuracy, but with nearly 25% increase in training time captured based on three runs using a system with 3.00-GHz 48-core Intel Xeon Gold 6248R CPU, Nvidia Quadro RTX 8000 GPU, and 128 GB of RAM. In the next test cases, to better differentiate the models, the PiNNs are applied to more complex problems representing the complexity of solute transport in porous media.

4.2. Case 2: 2D Solute Transport with Constant Velocity

This test case explores solute transport in a 2D rectangular domain (

x \in [0, L], y \in [0, W]

) representing an isotropic porous medium with a constant velocity field,

u = (u_{x}, u_{y})

. Since the velocity field is given, the governing equation for solute transport, Equation (7), is simplified to,

\frac{\partial C}{\partial t} + u_{x} \frac{\partial C}{\partial x} + u_{y} \frac{\partial C}{\partial y} = \frac{\partial}{\partial x} (D_{x} \frac{\partial C}{\partial x}) + \frac{\partial}{\partial y} (D_{y} \frac{\partial C}{\partial y}),

(31)

where

D_{y} = D_{x} = 0.02

^{2}

/s refers to the hydrodynamic dispersion coefficient in the y and x directions. As schematically shown in Figure 3, the following initial and boundary conditions,

\{\begin{matrix} C (x, y, t) = 0; t = 0, \\ C (x, y, t) = C_{0}; x = 0 a n d y_{1} \leq y \leq y_{2}, \\ C (x, y, t) = 0; x = 0 a n d y < y_{1} a n d y > y_{2}, \\ \frac{\partial}{\partial x} C (x, y, t) = 0; x = L, a n d 0 \leq y \leq W, \\ \frac{\partial}{\partial y} C (x, y, t) = 0; y = 0 o r y = W, a n d 0 \leq x \leq L, \end{matrix}

(32)

are considered for this test case. The blue points on the left boundary bounded by

y_{1}

and

y_{2}

represent the positions of point-source solute injection at the rate of

C_{0} = 0.2

kg/m

^{3}

. The other three sides of the domain are assigned as zero-gradient boundaries, assuming a free outflow of solute. This 2D advection-dispersion problem defined by Equations (31) and (32) can be solved analytically [69],

C (x, y, t) = C_{0} \sum_{n = 0}^{\infty} L_{n} P_{n} cos (η y) (exp (\frac{x (u_{x} - ζ)}{2 D_{x}}) erfc (\frac{x - ζ t}{2 \sqrt{D_{x} t}}) + exp (\frac{x (u_{x} + ζ)}{2 D_{x}}) erfc (\frac{x + ζ t}{2 \sqrt{D_{x} t}})),

(33)

where

L_{n}

is defined as,

L_{n} = \{\begin{matrix} \frac{1}{2}; n = 0, \\ 1; n > 0, \end{matrix}

(34)

and

P_{n}

is computed as,

P_{n} = \{\begin{matrix} (y_{2} - y_{1}) / W; n = 0, \\ (sin (η y_{2}) - sin (η y_{1})) / (n π); n > 0, \end{matrix}

(35)

where

η

and

ζ

, respectively, are,

η = \frac{n π}{W} (n = 0, 1, 2, 3, \dots), ζ = \sqrt{{u_{x}}^{2} + 4 η^{2} D_{x} D_{y}} .

(36)

First, a PiNN is employed with sin activation function and randomly distributed collocation points to predict the solute transport in the 2D domain (

L = W = 1

m) considering the dispersion only with

D_{x} = D_{y} = 0.02

^{2} /

s and

u = (u_{x}, u_{y})

= (0.0, 0.0) m/s in Equation (31) and the initial and boundary conditions given by Equation (32). The solute injection rate is set to

C_{0} =

0.2 kg/m

^{3}

between

y_{1} = 0.3

m and

y_{2} = 0.7

m. To determine the best architecture, similar to the previous case, an iterative random-search hyperparameter tuning approach is practiced. This approach produced a network with four hidden layers and {32, 16, 16, 16} neurons per layer, 8000 spatiotemporal random points within

0 \leq x, y \leq 1, 0 \leq t \leq 1

as collocation points to enforce the PDE loss term, 8000 spatiotemporal random points to collectively represent the boundary and initial conditions informing the loss terms, and identical weights for all loss terms (

w_{i} = w_{b} = w_{p} = 1

Figure 3 shows the comparison of the 2D concentration fields predicted by PiNNs with sin and tanh activation functions and the analytical solution given by Equation (33). The absolute point error is also shown, which illustrates a good agreement between the analytic solution and PiNN’s prediction. Based on the total loss plot, it can be observed that the PiNN with sin activation function is faster to converge and more accurate. The mismatch with the analytical solution in both cases is considerable at locations where extremely high concentration gradients exist (e.g., close to points

(0.0, 0.30)

and

(0.0, 0.70)

). This mismatch is due to a difficult-to-minimize approximation error, i.e., PiNN struggles to converge to ground truth solutions close to those points due to the inherent complexity that arises in high-gradient areas and the limited capacity of the network architecture and training procedure [70,71]. A PiNN with a deeper network may resolve those areas and converge to the correct solution, but that makes minimization error a harder task.

Two PiNN models with sin and tanh activation functions and the same architecture as discussed above are then employed to model the advection–dispersion solute transport in 2D domains with

D_{x} = D_{y} = 0.02

^{2} /

s and

u = (u_{x}, u_{y}) = (0.5, 0.0)

m/s. The solute is injected again at a rate of

C_{0} = 0.2

kg/m

^{3}

between

y_{1} = 0.3

m and

y_{2} = 0.7

m. Table 1 reports the results of the comparison of these two models against the ground truth. The concentration fields predicted by PiNN with sin activation function at different times (

t = 0.25

s, 0.50 s, 0.75 s, 1.00 s) are shown in Figure 4. The predictions of PiNNs for the concentration field at

t =

1.00 s are also compared with the analytical solution given by Equation (31). The absolute point error shows that there is a good match between the PiNNs’ predictions, with sin, and the analytical solution within the entire domain yielding

M S E = 1.54 \times 10^{- 6}

that is about an order of magnitude more accurate compared to the PiNN with tanh activation function while being nearly 31% faster to be trained.

4.3. Case 3: 2D Solute Transport in Homogeneous Porous Media

In this test case, the solute transport in a 2D homogeneous porous medium (

x \in [0, L], y \in [0, W]

) is investigated. For this problem, assuming

ζ (x, y) = 1

, the pressure equation defined by Equation (9) reduces to,

\frac{\partial^{2} P}{\partial x^{2}} + \frac{\partial^{2} P}{\partial y^{2}} = 0,

(37)

that must be solved to determine the velocity field using Equation (3). The governing equation for the solute concentration remains the same as Equation (31). The following boundary conditions are considered for the pressure field,

\{\begin{matrix} P (x, y) = 1.0; y = 0 o r y = W, a n d 0 \leq x \leq L, \\ P (x, y) = 0.1; x = 0 o r x = L, a n d 0 \leq y \leq W, \end{matrix}

(38)

and for the concentration field,

\{\begin{matrix} C (x, y, t) = 0; y = 0 o r y = W a n d x < x_{1} a n d x > x_{2}, \\ C (x, y, t) = C_{0}; y = 0 o r y = W a n d x_{1} \leq x \leq x_{2}, \\ \frac{\partial}{\partial x} C (x, y, t) = 0; x = 0 o r x = L, a n d 0 \leq y \leq W, \end{matrix}

(39)

with

C (x, y, t = 0) = 0

as the initial condition. The pressure value at the corner points that belong to two perpendicular sides is set to an average value of

0.55

, as additional constraints, to maintain symmetry in the solution. The ground truth solutions for pressure and concentration fields are obtained using FEM with

100 \times 100

quadrilateral elements assuming

L = W = 1 m

D_{x} = D_{y} = 0.02

^{2}

/s, and

C_{0} = 0.2

kg/m

^{3}

A PiNN with sin activation function and randomly distributed collocation points is employed to predict the solute transport in the 2D porous media described above. The PiNN’s inputs are the spatiotemporal coordinates,

(x, y, t)

, and the outputs are the pressure and concentration fields, i.e.,

s = {P, C}

in Figure 1. The pressure and concentration fields are decoupled, therefore one may use two PiNNs side-by-side to solve this problem; however, training a PiNN with multiple outputs is more desirable in this study to minimize hyperparameter tuning. The iterative random-search hyperparameter tuning process is again practiced to obtain the best PiNN architecture based on the MSE accuracy measure. This approach resulted in a network with five hidden layers and {32, 16, 16, 8, 8} neurons per layer using 12,000 spatiotemporal random points within

0 \leq x, y \leq 1, 0 \leq t \leq 1

as collocation points to enforce the pressure and concentration PDEs’ loss terms, 12,000 spatiotemporal random points to collectively represent the boundary and initial conditions loss terms (sample selected points are shown in Figure 5), and identical weights for all loss terms (

w_{i} = w_{b} = w_{p} = 1

). Again, a PiNN with similar architecture but using tanh activation function is also trained and tuned for comparison purposes, see Table 2.

Figure 5 depicts a comparison between ground truth (FEM solutions) and PiNN’s predictions, with sin activation function, for the pressure field, the velocity field in the x direction, the total velocity field, and the flow streamlines in a homogeneous porous medium. The absolute point error is also displayed for all fields, indicating good agreement between PiNN and FEM solutions. Figure 5 also demonstrates the comparison between the PiNN’s prediction and FEM solution for the concentration field at

t = 1.00

s. The absolute point error indicates the agreement between both solutions. The results reported in Table 2 reveal that the PiNN with sin activation function is almost an order of magnitude more accurate than its counterpart (the PiNN with tan activation function), while its training time is reduced by 78%. Furthermore, the faster convergence of PiNN with the sin activation function is noticeable in the total loss plot versus iteration, as illustrated in Figure 5. Based on these results, it can be inferred that the proposed PiNN is capable of simultaneously solving for the pressure and solute concentration fields in a 2D homogeneous porous medium with high accuracy. The next section of the study explores the application of PiNN and the benefit of using a periodic activation function to solve a more complex problem of solute transport imposed by permeability heterogeneity in porous media.

4.4. Case 4: 2D Solute Transport in Heterogeneous Porous Media

In this test case, the solute transport is analyzed in a 2D rectangular domain (

x \in [0, L], y \in [0, W]

) representing a heterogeneous porous medium, i.e., the permeability

k (x, y)

, and hence,

ζ (x, y)

in Equation (9) vary spatially. Three different

ζ (x, y)

fields, i.e., scaled permeability fields, are considered as shown in Figure 6. These fields are designed to exhibit structural features with distinct length scales relative to the size of the domain. As shown in Figure 6, from left to right, i.e., Case 4A, 4B, and 4C, the length scale of the structural features in the permeability field decreases, leading to more complex problems to simulate using PiNNs owing to the increase in high-frequency features. The initial and boundary conditions for these cases are identical to those reported for Case 3 (refer to Section 4.3). The pressure constraints at the corner points that correspond to the two perpendicular sides are not maintained due to variations in permeability. For each case, the ground truth solutions for the pressure and concentration fields are determined using FEM with

100 \times 100

quadrilateral elements considering

L = W = 1 m

α = 0

D_{x} = D_{y} = 0.02

^{2}

/s,

μ ϕ = 0.001

, and

C_{0} =

0.2 kg/m

^{3}

Two PiNN models with five hidden layers and {32, 32, 16, 16, 16} neurons per layer are trained using sin and tanh activation functions, respectively, to resolve the solute transport in these three heterogeneous cases (i.e., predict the pressure and concentration fields). For each case, due to changes in the permeability field, the hyperparameter tuning process leads to different randomly distributed collocation points. For Case 4A, a total of 30,000 spatiotemporal random points (15,000 points within

0 \leq x, y \leq 1, 0 \leq t \leq 1

and 15,000 points on the boundaries) are used as collocation points to enforce the PDE, initial, and boundary loss terms. For Case 4B and Case 4C, the hyperparameter tuning process leads to a total of 36,000 and 40,000 spatiotemporal random points as collocation points, respectively. The weights for the loss terms are also tuned for both PiNN models to achieve solution convergence while meeting the boundary condition requirements satisfactorily. Table 3 reports the collocation points within the domain and the weights of the loss terms used to predict the pressure and concentration fields. It also reports the comparison of the PiNN models, in terms of accuracy and training time, to predict the ground truth solution. Figure 7 and Figure 8 illustrate the comparison between the PiNNs’ predictions with the ground truth (FEM solutions) obtained for solute transport in heterogeneous porous media. The comparisons are shown for the steady state pressure field (Figure 7) and the concentration field at

t = 1.00 s

(Figure 8).

Figure 7 shows a comparison of PiNN models using sin and tanh activation functions in predicting the pressure field in heterogeneous porous media represented by Case 4A, 4B, and 4C scaled permeability fields (Figure 6). The absolute point error for each case, compared to the ground truth obtained using FEM with

100 \times 100

quadrilateral elements, shows that the PiNN with the sin activation function provides a more accurate prediction for the pressure field. Table 3 reports the comparison in terms of MSE (Equation (24)). Based on the MSE reported, the PiNN model with sin activation function is an order of magnitude more accurate in predicting the pressure field. As can be seen in Figure 7, the PiNN model with the tanh activation function encounters more difficulty in converging to the ground truth solution as the variation of permeability in the domain increases (i.e., as the length-scale of structural features decreases). In contrast, this issue is less pronounced for the PiNN model with sin activation function.

Figure 8 shows a comparison of PiNN models using sin and tanh activation functions in predicting the concentration field in heterogeneous porous media represented by Case 4A, 4B, and 4C scaled permeability fields (Figure 6). The PiNNs’ predictions are validated against the FEM solution with a

100 \times 100

quadrilateral element grid. The comparison is illustrated for the concentration field at

t = 1.00

s considering

D_{x} = D_{y} = 0.02

^{2}

/s,

α = 0

μ ϕ = 0.001

, and

C_{0} = 0.2

kg/m

^{3}

. The absolute point error is also shown for each PiNN model to illustrate the mismatch between the prediction and ground truth solution. Table 3 reports the MSE calculated using Equation (24), and the training time. It is clearly observed that the mismatch for both PiNN models increases as the length-scale of the structural features in the permeability field decreases (from left to right in Figure 8), which is attributed to the rapid variation of the velocity field caused by heterogeneities. As shown, the PiNN model with sin activation function encounters less difficulty in converging to the ground truth solution, and this is more noticeable as the variation of permeability in the domain increases. In contrast, the mismatch is smaller for the PiNN model with sin activation function, making it two orders of magnitude more accurate and nearly 2x faster to train compared to the PiNN model using tanh activation function. This order of magnitude for accuracy is consistent for all time steps. Figure 9 shows an example comparison between the PiNN’s predictions (with sin activation function) and ground truth solutions for transient solute transport in a heterogeneous porous medium with a permeability field described by Case 4B. The MSE is calculated using Equation (24) for the concentration fields obtained at each given time. These comparative findings clearly show that a PiNN with a periodic activation function can solve for the pressure and solute concentration fields in 2D heterogeneous media with a higher degree of accuracy compared to tanh activation function. The results also show that a PiNN with the periodic activation function reduces the training time.

Table 4 compares the inference time for PiNN with the sin activation function and FEM resolving solute transport in 2D homogeneous (case 3) and heterogeneous (Case 4) porous media. The speed-up factor achieved is around three orders of magnitude and can considerably escalate when a higher-resolution mesh (i.e., a higher number of elements) is used by the FEM solver; a situation that is often faced in cases with higher levels of heterogeneities. In fact, in such cases, the number of collocation points used by PiNN also increases, which could negatively impact the cost of training but not the inference time. Therefore, the computational efficiency of PiNN is certainly superior to that of a traditional solver (e.g., FEM) in domains with complex geometry and heterogeneity. It is also notable that PiNN comes with several drawbacks. The two most common drawbacks that were also faced in this study are: (i) PiNN’s training can face gradient vanishing problems, which can prohibitively slow down the learning process; and (ii) PiNN did not always converge due to competing non-linear loss terms (e.g., PDE loss terms related to the pressure and concentration fields in addition to the loss terms related to initial/boundary conditions for both of the fields). It requires trials and errors or other adaptive techniques to adjust the loss terms’ weight functions to mitigate the instability. The lack of a theoretical condition or constraint (e.g., Courant number in traditional computational fluid dynamics [72]) to ensure convergence is an open research area for investigation. Overall, considering the effectiveness of PiNN in combining scientific computing and DL and accelerating computation, the addition of PiNNs to the simulations of porous media flow should be further investigated. For example, a future study could be the assessment of PiNNs to model non-isothermal reactive flows with anisotropic dispersion coefficient in highly heterogeneous and deformable porous media.

5. Conclusions

In this study, we demonstrated the use of physics-informed neural networks (PiNNs) and the advantages of employing a periodic activation function in addressing the solute transport problem in both homogeneous and heterogeneous porous media, as governed by the advection-dispersion equation. We constructed PiNNs using the sin and tanh activation functions to predict pressure and concentration fields within a given spatiotemporal domain. To evaluate the capabilities of PiNN models, we conducted seven case studies (1D and 2D) with varying degrees of permeability heterogeneity. We utilized an iterative random-search hyperparameter tuning method to determine the optimal architecture for each PiNN model in the test cases. The accuracy of the PiNNs’ predictions was evaluated using absolute point error and mean square error metrics, and compared to the ground truth solutions obtained analytically or through FEM numerical methods. Our results showed that the PiNN with sin activation function was capable of accurately predicting the behavior of the solute transport under a variety of conditions. The PiNN model employing the sin activation function was found to be up to two orders of magnitude more accurate and up to two times faster to train than the commonly used tanh activation function when applied to 2D homogeneous and heterogeneous porous media exhibiting structural features with distinct length scales relative to the domain size. Furthermore, the inference speed-up results showed that the PiNN’s simultaneous predictions of pressure and concentration fields can reduce computational time by three orders of magnitude compared to FEM simulations.

Author Contributions

Conceptualization, S.A.F. and S.F.; investigation, S.A.F., P.D., R.S. and S.K.M.; writing—original draft preparation, P.D., S.A.F. and S.K.M.; writing—review and editing, S.A.F. and S.F.; funding acquisition, S.A.F.; All authors have read and agreed to the published version of the manuscript.

Funding

S.A.F. would like to acknowledge the support by the Texas Sate University’s International Research Accelerator Grant (award no. 9000003039).

Data Availability Statement

Data will be provided upon request.

Conflicts of Interest

The authors declare no conflict of interests.

References

Cook, P.G.; Böhlke, J.K. Determining timescales for groundwater flow and solute transport. In Environmental Tracers in Subsurface Hydrology; Springer: Berlin/Heidelberg, Germany, 2000; pp. 1–30. [Google Scholar] [CrossRef]
Gaus, I.; Audigane, P.; André, L.; Lions, J.; Jacquemet, N.; Durst, P.; Czernichowski-Lauriol, I.; Azaroual, M. Geochemical and solute transport modelling for CO₂ storage, what to expect from it? Int. J. Greenh. Gas Control. 2008, 2, 605–625. [Google Scholar] [CrossRef]
Pruess, K.; Yabusaki, S.; Steefel, C.; Lichtner, P. Fluid flow, heat transfer, and solute transport at nuclear waste storage tanks in the Hanford vadose zone. Vadose Zone J. 2002, 1, 68–88. [Google Scholar] [CrossRef]
Bienert, G.P.; Schjoerring, J.K.; Jahn, T.P. Membrane transport of hydrogen peroxide. Biochim. Biophys. Acta (BBA)-Biomembr. 2006, 1758, 994–1003. [Google Scholar] [CrossRef] [PubMed]
Kristensen, E.; Hansen, K. Transport of carbon dioxide and ammonium in bioturbated (Nereis diversicolor) coastal, marine sediments. Biogeochemistry 1999, 45, 147–168. [Google Scholar] [CrossRef]
Li, X.; Li, D.; Xu, Y.; Feng, X. A DFN based 3D numerical approach for modeling coupled groundwater flow and solute transport in fractured rock mass. Int. J. Heat Mass Transf. 2020, 149, 119179. [Google Scholar] [CrossRef]
Hasan, S.; Niasar, V.; Karadimitriou, N.K.; Godinho, J.R.; Vo, N.T.; An, S.; Rabbani, A.; Steeb, H. Direct characterization of solute transport in unsaturated porous media using fast X-ray synchrotron microtomography. Proc. Natl. Acad. Sci. USA 2020, 117, 23443–23449. [Google Scholar] [CrossRef] [PubMed]
Faraji, M.; Mazaheri, M. Mathematical model of solute transport in rivers with storage zones using nonlinear dispersion flux approach. Hydrol. Sci. J. 2022, 67, 1656–1668. [Google Scholar] [CrossRef]
Yang, X.R.; Wang, Y. Ubiquity of anomalous transport in porous media: Numerical evidence, continuous time random walk modelling, and hydrodynamic interpretation. Sci. Rep. 2019, 9, 4601. [Google Scholar] [CrossRef]
Zhao, Z.; Jing, L.; Neretnieks, I.; Moreno, L. Numerical modeling of stress effects on solute transport in fractured rocks. Comput. Geotech. 2011, 38, 113–126. [Google Scholar] [CrossRef]
Zhang, Z.h.; Zhang, J.p.; Ju, Z.y.; Zhu, M. A one-dimensional transport model for multi-component solute in saturated soil. Water Sci. Eng. 2018, 11, 236–242. [Google Scholar] [CrossRef]
Bagalkot, N.; Suresh Kumar, G. Effect of nonlinear sorption on multispecies radionuclide transport in a coupled fracture-matrix system with variable fracture aperture: A numerical study. ISH J. Hydraul. Eng. 2015, 21, 242–254. [Google Scholar] [CrossRef]
Mostaghimi, P.; Liu, M.; Arns, C.H. Numerical simulation of reactive transport on micro-CT images. Math. Geosci. 2016, 48, 963–983. [Google Scholar] [CrossRef]
Maheshwari, P.; Ratnakar, R.; Kalia, N.; Balakotaiah, V. 3-D simulation and analysis of reactive dissolution and wormhole formation in carbonate rocks. Chem. Eng. Sci. 2013, 90, 258–274. [Google Scholar] [CrossRef]
Noetinger, B.; Roubinet, D.; Russian, A.; Le Borgne, T.; Delay, F.; Dentz, M.; De Dreuzy, J.R.; Gouze, P. Random walk methods for modeling hydrodynamic transport in porous and fractured media from pore to reservoir scale. Transp. Porous Media 2016, 115, 345–385. [Google Scholar] [CrossRef]
Im, S.; Lee, J.; Cho, M. Surrogate modeling of elasto-plastic problems via long short-term memory neural networks and proper orthogonal decomposition. Comput. Methods Appl. Mech. Eng. 2021, 385, 114030. [Google Scholar] [CrossRef]
Karniadakis, G.E.; Kevrekidis, I.G.; Lu, L.; Perdikaris, P.; Wang, S.; Yang, L. Physics-informed machine learning. Nat. Rev. Phys. 2021, 3, 422–440. [Google Scholar] [CrossRef]
Raissi, M.; Perdikaris, P.; Karniadakis, G.E. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. J. Comput. Phys. 2019, 378, 686–707. [Google Scholar] [CrossRef]
Han, J.; Jentzen, A.; E, W. Solving high-dimensional partial differential equations using deep learning. Proc. Natl. Acad. Sci. USA 2018, 115, 8505–8510. [Google Scholar] [CrossRef]
Berg, J.; Nyström, K. A unified deep artificial neural network approach to partial differential equations in complex geometries. Neurocomputing 2018, 317, 28–41. [Google Scholar] [CrossRef]
Sirignano, J.; Spiliopoulos, K. DGM: A deep learning algorithm for solving partial differential equations. J. Comput. Phys. 2018, 375, 1339–1364. [Google Scholar] [CrossRef]
Faroughi, S.A.; Pawar, N.; Fernandes, C.; Das, S.; Kalantari, N.K.; Mahjour, S.K. Physics-Guided, Physics-Informed, and Physics-Encoded Neural Networks in Scientific Computing. arXiv 2022, arXiv:2211.07377. [Google Scholar] [CrossRef]
Van Merriënboer, B.; Breuleux, O.; Bergeron, A.; Lamblin, P. Automatic differentiation in ML: Where we are and where we should be going. Adv. Neural Inf. Process. Syst. 2018, 31. [Google Scholar] [CrossRef]
Baydin, A.G.; Pearlmutter, B.A.; Radul, A.A.; Siskind, J.M. Automatic differentiation in machine learning: A survey. J. Marchine Learn. Res. 2018, 18, 1–43. [Google Scholar] [CrossRef]
Raissi, M.; Perdikaris, P.; Karniadakis, G.E. Physics informed deep learning (part i): Data-driven solutions of nonlinear partial differential equations. arXiv 2017, arXiv:1711.10561. [Google Scholar] [CrossRef]
Raissi, M.; Yazdani, A.; Karniadakis, G.E. Hidden fluid mechanics: A Navier-Stokes informed deep learning framework for assimilating flow visualization data. arXiv 2018, arXiv:1808.04327. [Google Scholar] [CrossRef]
Almajid, M.M.; Abu-Al-Saud, M.O. Prediction of porous media fluid flow using physics informed neural networks. J. Pet. Sci. Eng. 2022, 208, 109205. [Google Scholar] [CrossRef]
Hanna, J.M.; Aguado, J.V.; Comas-Cardona, S.; Askri, R.; Borzacchiello, D. Residual-based adaptivity for two-phase flow simulation in porous media using Physics-informed Neural Networks. Comput. Methods Appl. Mech. Eng. 2022, 396, 115100. [Google Scholar] [CrossRef]
Haghighat, E.; Amini, D.; Juanes, R. Physics-informed neural network simulation of multiphase poroelasticity using stress-split sequential training. Comput. Methods Appl. Mech. Eng. 2022, 397, 115141. [Google Scholar] [CrossRef]
He, Q.; Barajas-Solano, D.; Tartakovsky, G.; Tartakovsky, A.M. Physics-informed neural networks for multiphysics data assimilation with application to subsurface transport. Adv. Water Resour. 2020, 141, 103610. [Google Scholar] [CrossRef]
He, Q.; Tartakovsky, A.M. Physics-Informed neural network method for forward and backward advection-dispersion equations. Water Resour. Res. 2021, 57, e2020WR029479. [Google Scholar] [CrossRef]
Vadyala, S.R.; Betgeri, S.N.; Betgeri, N.P. Physics-informed neural network method for solving one-dimensional advection equation using PyTorch. Array 2022, 13, 100110. [Google Scholar] [CrossRef]
Rodriguez-Torrado, R.; Ruiz, P.; Cueto-Felgueroso, L.; Green, M.C.; Friesen, T.; Matringe, S.; Togelius, J. Physics-informed attention-based neural network for hyperbolic partial differential equations: Application to the Buckley–Leverett problem. Sci. Rep. 2022, 12, 7557. [Google Scholar] [CrossRef] [PubMed]
Zhang, W.; Al Kobaisi, M. On the Monotonicity and Positivity of Physics-Informed Neural Networks for Highly Anisotropic Diffusion Equations. Energies 2022, 15, 6823. [Google Scholar] [CrossRef]
Fuks, O.; Tchelepi, H.A. Limitations of physics informed machine learning for nonlinear two-phase transport in porous media. J. Mach. Learn. Model. Comput. 2020, 1, 19–37. [Google Scholar] [CrossRef]
Zhang, Z.; Yan, X.; Liu, P.; Zhang, K.; Han, R.; Wang, S. A physics-informed convolutional neural network for the simulation and prediction of two-phase Darcy flows in heterogeneous porous media. J. Comput. Phys. 2023, 477, 111919. [Google Scholar] [CrossRef]
Jagtap, A.; Karniadakis, G.E. How important are activation functions in regression and classification? A survey, performance comparison, and future directions. J. Mach. Learn. Model. Comput. 2022, 4, 21–75. [Google Scholar] [CrossRef]
Parascandolo, G.; Huttunen, H.; Virtanen, T. Taming the waves: Sine as activation function in deep neural networks. In Proceedings of the ICLR 2017 Conference Track, Toulon, France, 24–26 April 2017. [Google Scholar]
Zhang, C.; Kaito, K.; Hu, Y.; Patmonoaji, A.; Matsushita, S.; Suekane, T. Influence of stagnant zones on solute transport in heterogeneous porous media at the pore scale. Phys. Fluids 2021, 33, 036605. [Google Scholar] [CrossRef]
Khan, S.; Alhazmi, S.E.; Alotaibi, F.M.; Ferrara, M.; Ahmadian, A. On the Numerical Approximation of Mobile-Immobile Advection-Dispersion Model of Fractional Order Arising from Solute Transport in Porous Media. Fractal Fract. 2022, 6, 445. [Google Scholar] [CrossRef]
Zhao, X.; Toksoz, M.N. Solute transport in heterogeneous porous media. Mass. Inst. Technol. Earth Resour. Lab. 1994, 145, 151–177. [Google Scholar] [CrossRef]
Van Genuchten, M.T.; Alves, W. Analytical Solutions of the One-Dimensional Convective-Dispersive Solute Transport Equation; Technical Bulletin (USA); United States Department of Agriculture: Washington, DC, USA, 1982. [Google Scholar] [CrossRef]
Sun, L.; Qiu, H.; Wu, C.; Niu, J.; Hu, B.X. A review of applications of fractional advection–dispersion equations for anomalous solute transport in surface and subsurface water. Wiley Interdiscip. Rev. Water 2020, 7, e1448. [Google Scholar] [CrossRef]
Haigh, M.; Sun, L.; McWilliams, J.C.; Berloff, P. On eddy transport in the ocean. Part II: The advection tensor. Ocean. Model. 2021, 165, 101845. [Google Scholar] [CrossRef]
Lou, S.; Chen, S.s.; Lin, B.x.; Yu, J.; Yan, C. Effective high-order energy stable flux reconstruction methods for first-order hyperbolic linear and nonlinear systems. J. Comput. Phys. 2020, 414, 109475. [Google Scholar] [CrossRef]
Talon, L. On the statistical properties of fluid flows with transitional power-law rheology in heterogeneous porous media. J. Non-Newton. Fluid Mech. 2022, 304, 104789. [Google Scholar] [CrossRef]
Baioni, E.; Mousavi Nezhad, M.; Porta, G.M.; Guadagnini, A. Modeling solute transport and mixing in heterogeneous porous media under turbulent flow conditions. Phys. Fluids 2021, 33, 106604. [Google Scholar] [CrossRef]
Berkowitz, B.; Klafter, J.; Metzler, R.; Scher, H. Physical pictures of transport in heterogeneous media: Advection-dispersion, random-walk, and fractional derivative formulations. Water Resour. Res. 2002, 38, 9-1–9-12. [Google Scholar] [CrossRef]
Jagtap, A.D.; Kharazmi, E.; Karniadakis, G.E. Conservative physics-informed neural networks on discrete domains for conservation laws: Applications to forward and inverse problems. Comput. Methods Appl. Mech. Eng. 2020, 365, 113028. [Google Scholar] [CrossRef]
Jagtap, A.D.; Karniadakis, G.E. Extended Physics-informed Neural Networks (XPINNs): A Generalized Space-Time Domain Decomposition based Deep Learning Framework for Nonlinear Partial Differential Equations. Commun. Comput. Phys. 2020, 8, 2002–2041. [Google Scholar] [CrossRef]
Cuomo, S.; Di Cola, V.S.; Giampaolo, F.; Rozza, G.; Raissi, M.; Piccialli, F. Scientific Machine Learning through Physics-Informed Neural Networks: Where we are and What’s next. arXiv 2022, arXiv:2201.05624. [Google Scholar] [CrossRef]
Depina, I.; Jain, S.; Mar Valsson, S.; Gotovac, H. Application of physics-informed neural networks to inverse problems in unsaturated groundwater flow. Georisk Assess. Manag. Risk Eng. Syst. Geohazards 2022, 16, 21–36. [Google Scholar] [CrossRef]
Lu, L.; Meng, X.; Mao, Z.; Karniadakis, G.E. DeepXDE: A deep learning library for solving differential equations. SIAM Rev. 2021, 63, 208–228. [Google Scholar] [CrossRef]
Sitzmann, V.; Martel, J.N.P.; Bergman, A.W.; Lindell, D.B.; Wetzstein, G. Implicit Neural Representations with Periodic Activation Functions. arXiv 2020, arXiv:2006.09661. [Google Scholar] [CrossRef]
Rumelhart, D.E.; Hinton, G.E.; Williams, R.J. Learning representations by back-propagating errors. Nature 1986, 323, 533–536. [Google Scholar] [CrossRef]
Guo, Y.; Cao, X.; Liu, B.; Gao, M. Solving partial differential equations using deep learning and physical constraints. Appl. Sci. 2020, 10, 5917. [Google Scholar] [CrossRef]
Cai, S.; Mao, Z.; Wang, Z.; Yin, M.; Karniadakis, G.E. Physics-informed neural networks (PINNs) for fluid mechanics: A review. Acta Mech. Sin. 2022, 37, 1727–1738. [Google Scholar] [CrossRef]
Strelow, E.L.; Gerisch, A.; Lang, J.; Pfetsch, M.E. Physics informed neural networks: A case study for gas transport problems. J. Comput. Phys. 2023, 481, 112041. [Google Scholar] [CrossRef]
Cuomo, S.; Giampaolo, F.; Izzo, S.; Nitsch, C.; Piccialli, F.; Trombetti, C. A physics-informed learning approach to Bernoulli-type free boundary problems. Comput. Math. Appl. 2022, 128, 34–43. [Google Scholar] [CrossRef]
Shah, K.; Stiller, P.; Hoffmann, N.; Cangi, A. Physics-Informed Neural Networks as Solvers for the Time-Dependent Schrödinger Equation. arXiv 2022, arXiv:2210.12522. [Google Scholar] [CrossRef]
Kingma, D.P.; Ba, J. Adam: A method for stochastic optimization. arXiv 2014, arXiv:1412.6980. [Google Scholar]
Byrd, R.H.; Lu, P.; Nocedal, J.; Zhu, C. A limited memory algorithm for bound constrained optimization. SIAM J. Sci. Comput. 1995, 16, 1190–1208. [Google Scholar] [CrossRef]
Bengio, Y. Gradient-based optimization of hyperparameters. Neural Comput. 2000, 12, 1889–1900. [Google Scholar] [CrossRef]
Kylasa, S.; Roosta, F.; Mahoney, M.W.; Grama, A. GPU accelerated sub-sampled Newton’s method for convex classification problems. In Proceedings of the 2019 SIAM International Conference on Data Mining, SIAM, Calgary, AB, Canada, 2–4 May 2019; pp. 702–710. [Google Scholar] [CrossRef]
Richardson, A. Seismic full-waveform inversion using deep learning tools and techniques. arXiv 2018, arXiv:1801.07232. [Google Scholar] [CrossRef]
Olmo, A.; Zamzam, A.; Glaws, A.; King, R. Physics-Driven Convolutional Autoencoder Approach for CFD Data Compressions. arXiv 2022, arXiv:2210.09262. [Google Scholar] [CrossRef]
Escapil-Inchauspé, P.; Ruz, G.A. Hyper-parameter tuning of physics-informed neural networks: Application to Helmholtz problems. Neurocomputing 2023, 561, 126826. [Google Scholar] [CrossRef]
Rasht-Behesht, M.; Huber, C.; Shukla, K.; Karniadakis, G.E. Physics-Informed Neural Networks (PINNs) for Wave Propagation and Full Waveform Inversions. J. Geophys. Res. Solid Earth 2022, 127, e2021JB023120. [Google Scholar] [CrossRef]
Zhou, J.G. A lattice Boltzmann method for solute transport. Int. J. Numer. Methods Fluids 2009, 61, 848–863. [Google Scholar] [CrossRef]
Yu, J.; Lu, L.; Meng, X.; Karniadakis, G.E. Gradient-enhanced physics-informed neural networks for forward and inverse PDE problems. Comput. Methods Appl. Mech. Eng. 2022, 393, 114823. [Google Scholar] [CrossRef]
Chiu, P.H.; Wong, J.C.; Ooi, C.; Dao, M.H.; Ong, Y.S. CAN-PINN: A fast physics-informed neural network based on coupled-automatic–numerical differentiation method. Comput. Methods Appl. Mech. Eng. 2022, 395, 114909. [Google Scholar] [CrossRef]
Atmakidis, T.; Kenig, E.Y. A study on the Kelvin-Helmholtz instability using two different computational fluid dynamics methods. J. Comput. Multiph. Flows 2010, 2, 33–45. [Google Scholar] [CrossRef]

Figure 1. A schematic architecture of Physics-informed Neural Networks (PiNNs). The network digests spatiotemporal coordinates, (x,t), as inputs to predict a solution set,

s

, as an approximate to the ground truth solution,

\hat{s}

. The automatic differentiation (AD) is then used to generate the derivatives of the predicted solution

s

with respect to inputs. These derivatives are used to formulate the residuals of the governing equations in the loss function weighted by different coefficients.

θ

and

λ

are the learnable parameters for weights/biases and unknown PDE parameters, respectively, that can be learned simultaneously while minimizing the loss function.

Figure 1. A schematic architecture of Physics-informed Neural Networks (PiNNs). The network digests spatiotemporal coordinates, (x,t), as inputs to predict a solution set,

s

, as an approximate to the ground truth solution,

\hat{s}

. The automatic differentiation (AD) is then used to generate the derivatives of the predicted solution

s

with respect to inputs. These derivatives are used to formulate the residuals of the governing equations in the loss function weighted by different coefficients.

θ

and

λ

are the learnable parameters for weights/biases and unknown PDE parameters, respectively, that can be learned simultaneously while minimizing the loss function.

Figure 2. One-dimensional solute transport with a constant velocity field,

u_{x} = 0.5

m/s. The upper panel shows the concentration field predicted by PiNN with sin activation function within the spatiotemporal domain

(0 \leq x \leq 1, 0 \leq t \leq 10)

at an injection rate of

C_{0} =

1 kg/m

^{3}

. The lower panels show the comparison of PiNNs’ prediction with the ground truth obtained analytically using Equation (29) at three different times. The comparison demonstrates a good agreement yielding

M S E = 1.15 \times 10^{- 6}

using sin activation function and

M S E = 1.21 \times 10^{- 6}

using tanh activation function for the entire spatiotemporal domain.

Figure 2. One-dimensional solute transport with a constant velocity field,

u_{x} = 0.5

m/s. The upper panel shows the concentration field predicted by PiNN with sin activation function within the spatiotemporal domain

(0 \leq x \leq 1, 0 \leq t \leq 10)

at an injection rate of

C_{0} =

1 kg/m

^{3}

M S E = 1.15 \times 10^{- 6}

using sin activation function and

M S E = 1.21 \times 10^{- 6}

using tanh activation function for the entire spatiotemporal domain.

Figure 3. A comparison between the predictions of PiNNs with sin and tanh activation functions and ground truth for solute transport in a 2D domain representing an isotropic porous medium considering

D_{x} = D_{y} = 0.02

^{2}

/s,

u = (u_{x}, u_{y}) = (0.0, 0.0)

m/s, and a solute injection rate of

C_{0}

= 0.2 kg/m

^{3}

between

y_{1} = 0.3

m and

y_{2}

= 0.7 m. The upper panels show the domain with boundary conditions and distributions of randomly selected points on which different terms in the loss function are evaluated. The red dots represent sample collocation points inside the domain corresponding to the loss term associated with the 2D solute transport PDE, and the blue, green, and black dots represent the points on the boundary of the domain corresponding to the loss terms associated with the boundary conditions. The lower panel shows the comparison between the PiNNs’ predictions for the concentration field at

t = 1.00

s and the ground truth obtained analytically using Equation (33). The absolute point error shows the mismatch between the solutions, and the total loss plot depicts the difference in the convergence rate of PiNNs.

D_{x} = D_{y} = 0.02

^{2}

/s,

u = (u_{x}, u_{y}) = (0.0, 0.0)

m/s, and a solute injection rate of

C_{0}

= 0.2 kg/m

^{3}

between

y_{1} = 0.3

m and

y_{2}

t = 1.00

Figure 4. A comparison between the prediction of PiNN with sin and tanh activation functions and ground truth obtained analytically for the concentration field considering

D_{x} = D_{y} = 0.02

^{2}

/s,

u = (u_{x}, u_{y}) = (0.5, 0.0)

m/s, and an injection rate of

C_{0} = 0.2

kg/m

^{3}

between

y_{1} = 0.3

m and

y_{2} = 0.7

m. The upper panels show the predictions of the PiNN with sin activation function at

t =

0.25 s, 0.50 s, 0.75 s, and

1.00

s. At

t = 1.00

s, the PiNNs’ predictions are compared with the analytical solution obtained using Equation (33). The absolute point error shows the mismatch between the solutions yielding the MSEs reported in Table 1.

Figure 4. A comparison between the prediction of PiNN with sin and tanh activation functions and ground truth obtained analytically for the concentration field considering

D_{x} = D_{y} = 0.02

^{2}

/s,

u = (u_{x}, u_{y}) = (0.5, 0.0)

m/s, and an injection rate of

C_{0} = 0.2

kg/m

^{3}

between

y_{1} = 0.3

m and

y_{2} = 0.7

m. The upper panels show the predictions of the PiNN with sin activation function at

t =

0.25 s, 0.50 s, 0.75 s, and

1.00

s. At

t = 1.00

Figure 5. A comparison between the prediction of the PiNN with sin activation function and the ground truth for the pressure field in a 2D homogeneous porous medium considering

D_{x} = D_{y} = 0.02

^{2} /

s and

C_{0} = 0.2

kg/m

^{3}

. The upper panels show the domain with pressure and concentration boundary conditions, as well as the comparison of total loss vs. iteration for PiNNs using sin and tanh activation functions. The lower panels show the comparison between the PiNN’s predictions and ground truth solutions obtained using FEM with

100 \times 100

quadrilateral elements for the pressure field, the velocity field in the x direction, and the total velocity field under steady state conditions as well as the concentration field at

t = 1.00

s. The flow direction of the pore fluid is also represented by the streamlines. The absolute point error is also shown for each field to illustrate the mismatch between the solutions.

Figure 5. A comparison between the prediction of the PiNN with sin activation function and the ground truth for the pressure field in a 2D homogeneous porous medium considering

D_{x} = D_{y} = 0.02

^{2} /

s and

C_{0} = 0.2

kg/m

^{3}

100 \times 100

quadrilateral elements for the pressure field, the velocity field in the x direction, and the total velocity field under steady state conditions as well as the concentration field at

t = 1.00

s. The flow direction of the pore fluid is also represented by the streamlines. The absolute point error is also shown for each field to illustrate the mismatch between the solutions.

Figure 6. Scaled permeability fields selected to possess structural features with different length-scales compared to the domain size. From left to right, i.e., Case 4A, 4B, and 4C, the length-scale of the structural features decreases, leading to more intricate problems to simulate using PiNNs owing to the existence of high-frequency features.

Figure 7. A comparison of PiNN models using sin and tanh activation functions in predicting the pressure field in heterogeneous porous media, validated against the FEM solution with a

100 \times 100

quadrilateral element grid. The comparison is illustrated for the pressure field considering

D_{x} = D_{y} = 0.02

^{2}

/s,

α = 0

μ ϕ = 0.001

, and

C_{0} = 0.2

kg/m

^{3}

. The absolute point error is also shown for each PiNN model to illustrate the mismatch between the prediction and ground truth solution. The PiNN model with the tanh activation function encounters more difficulty in converging to the ground truth solution as the variation of permeability in the domain increases (i.e., as the length-scale of structural features decreases). In contrast, this issue is less pronounced for the PiNN model with sin activation function.

Figure 7. A comparison of PiNN models using sin and tanh activation functions in predicting the pressure field in heterogeneous porous media, validated against the FEM solution with a

100 \times 100

quadrilateral element grid. The comparison is illustrated for the pressure field considering

D_{x} = D_{y} = 0.02

^{2}

/s,

α = 0

μ ϕ = 0.001

, and

C_{0} = 0.2

kg/m

^{3}

Figure 8. A comparison of PiNN models using sin and tanh activation functions to predict the concentration field in heterogeneous porous media, validated against the FEM solution with a

100 \times 100

quadrilateral element grid. The comparison is illustrated for the concentration field at

t = 1.00

s considering

D_{x} = D_{y} = 0.02

^{2} /

α = 0

μ ϕ = 0.001

, and

C_{0} = 0.2

kg/m

^{3}

. The absolute point error is also shown for each PiNN model to illustrate the mismatch between the predictions and ground truth solutions. The PiNN model with tanh activation function encounters more difficulty in converging to the ground truth solution as the variation of permeability in the domain increases (from left to right). This disparity is attributed to the rapid variation of the velocity field caused by heterogeneities (i.e., as the length-scale of structural features). In contrast, the mismatch is less pronounced for the PiNN model with sin activation function, making it two orders of magnitude more accurate than the PiNN model with tanh activation function.

Figure 8. A comparison of PiNN models using sin and tanh activation functions to predict the concentration field in heterogeneous porous media, validated against the FEM solution with a

100 \times 100

quadrilateral element grid. The comparison is illustrated for the concentration field at

t = 1.00

s considering

D_{x} = D_{y} = 0.02

^{2} /

α = 0

μ ϕ = 0.001

, and

C_{0} = 0.2

kg/m

^{3}

Figure 9. A comparison between the PiNN’s predictions (with sin activation function) and ground truth solutions for transient solute transport in a heterogeneous porous medium with the given permeability field, Case 4B. The MSE is calculated using Equation (24) for the concentration fields obtained at each given time.

Table 1. A comparison of PiNN models using sin and tanh activation functions in terms of accuracy and training time to predict the concentration field in homogeneous 1D and 2D porous media, validated against the analytical solutions. The MSE is calculated using Equation (24) in (t,x) domain for Case 1, and in (x, y) domain at

t = 1.0

(s) for Case 2 assuming

u = (u_{x}, u_{y}) = (0.0, 0.50)

m/s. The training time of the PiNNs is measured through three runs on a system equipped with a 3.00-GHz 48-core Intel Xeon Gold 6248R CPU, Nvidia Quadro RTX 8000 GPU, and 128 GB of RAM. Additionally, the number of collection points in the domain and the weights assigned to the loss terms are reported for each PiNN model.

t = 1.0

(s) for Case 2 assuming

u = (u_{x}, u_{y}) = (0.0, 0.50)

	Case 1		Case 2
	sin	tanh	sin	tanh
Collocation Points, $N_{p}$	5000	5000	8000	8000
Weights: $w_{i}$ , $w_{b}$ , $w_{p}$	1.0, 1.0, 1.0	1.0, 1.0, 0.9	1.0, 1.0, 1.0	1.0, 1.0, 0.8
Accuracy (MSE)	1.15e-06	1.21e-06	1.54e-06	2.37e-05
Training Time (s)	$86 \pm 7$	$108 \pm 3$	$388 \pm 24$	$509 \pm 15$

Table 2. A comparison of PiNN models using sin and tanh activation functions in terms of accuracy and training time to predict the pressure and concentration fields in homogeneous 2D porous media, validated against the ground truth solutions obtained using FEM with

100 \times 100

quadrilateral elements. The MSE is calculated in (x,y) domain using Equation (24) for the pressure field at steady state, and for the concentration field at

t = 1.00

(s). The training time of the PiNNs is measured through three runs on a system equipped with a 3.00-GHz 48-core Intel Xeon Gold 6248R CPU, Nvidia Quadro RTX 8000 GPU, and 128 GB of RAM. Additionally, the number of collection points in the domain and the weights assigned to the loss terms are reported for each PiNN model.

100 \times 100

quadrilateral elements. The MSE is calculated in (x,y) domain using Equation (24) for the pressure field at steady state, and for the concentration field at

t = 1.00

	Case 3
	sin	tanh
Collocation Points, $N_{p}$	12,000	12,000
Initial and Boundary Loss Weights: $w_{i}$ , $w_{b}$	1.0, 1.0	1.0, 1.0
Pressure and Concentrations PDE Loss Weights: $w_{p}^{P}$ , $w_{p}^{C}$	1.0, 1.0	0.7, 0.9
Pressure MSE	2.84e-06	4.36e-05
Concentration MSE	1.22e-06	1.51e-05
Training Time (s)	$946 \pm 81$	$1687 \pm 53$

Table 3. A comparison of PiNN models using sin and tanh activation functions to predict the pressure and concentration fields in heterogeneous porous media, validated against the FEM solution with a

100 \times 100

quadrilateral element grid. The MSE is calculated in (x, y) domain using Equation (24) for the pressure field at steady state and for the concentration field at

t = 1.00

s. The training time of the PiNNs is measured through three runs on a system equipped with a 3.00-GHz 48-core Intel Xeon Gold 6248R CPU, Nvidia Quadro RTX 8000 GPU, and 128 GB of RAM. Additionally, the number of collection points in the domain and the weights assigned to the loss terms are reported for each PiNN model.

100 \times 100

quadrilateral element grid. The MSE is calculated in (x, y) domain using Equation (24) for the pressure field at steady state and for the concentration field at

t = 1.00

	Case 4A		Case 4B		Case 4C
	sin	tanh	sin	tanh	sin	tanh
$N_{p}$	15,000	15,000	18,000	18,000	20,000	20,000
$w_{i}$ , $w_{b}$	1.0, 1.0	1.0, 1.0	1.0, 1.0	1.0, 1.0	1.0, 1.0	1.0, 1.0
$w_{p}^{P}$ , $w_{p}^{C}$	0.20, 0.70	0.25, 0.65	0.15, 0.80	0.25, 0.70	0.15, 0.80	0.3, 0.90
Pressure MSE	1.13e-05	1.38e-04	2.48e-05	1.96e-04	5.38e-05	2.08e-04
Concentration MSE	1.12e-06	1.10e-04	3.62e-06	1.25e-04	9.86e-06	2.72e-04
Training Time (s)	$1345 \pm 112$	$2710 \pm 73$	$1566 \pm 131$	$3038 \pm 92$	$2042 \pm 179$	$4113 \pm 143$

Table 4. A comparison of the inference time for the PiNN with sin activation function and FEM resolving the solute transport in the 2D homogeneous (Case 3) and heterogeneous (Case 4) porous media. The computational time is measured through three runs on a system equipped with a 3.00-GHz 48-core Intel Xeon Gold 6248R CPU, Nvidia Quadro RTX 8000 GPU, and 128 GB of RAM. Note: This comparison only shows the inference speed-up factor, not considering the PiNN training time in the calculation.

	FEM	PiNN (sin)	Speed-Up Factor
	(s)	(s)	(—)
2D Homogeneous	$183.8 \pm 0.1$	$0.126 \pm 0.01$	1458.7×
2D Heterogeneous-Case 4A	$363.2 \pm 0.3$	$0.259 \pm 0.02$	1402.3×
2D Heterogeneous-Case 4B	$363.9 \pm 0.3$	$0.260 \pm 0.01$	1399.6×
2D Heterogeneous-Case 4C	$364.5 \pm 0.2$	$0.260 \pm 0.03$	1401.9×

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

© 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Faroughi, S.A.; Soltanmohammadi, R.; Datta, P.; Mahjour, S.K.; Faroughi, S. Physics-Informed Neural Networks with Periodic Activation Functions for Solute Transport in Heterogeneous Porous Media. Mathematics 2024, 12, 63. https://doi.org/10.3390/math12010063

AMA Style

Faroughi SA, Soltanmohammadi R, Datta P, Mahjour SK, Faroughi S. Physics-Informed Neural Networks with Periodic Activation Functions for Solute Transport in Heterogeneous Porous Media. Mathematics. 2024; 12(1):63. https://doi.org/10.3390/math12010063

Chicago/Turabian Style

Faroughi, Salah A., Ramin Soltanmohammadi, Pingki Datta, Seyed Kourosh Mahjour, and Shirko Faroughi. 2024. "Physics-Informed Neural Networks with Periodic Activation Functions for Solute Transport in Heterogeneous Porous Media" Mathematics 12, no. 1: 63. https://doi.org/10.3390/math12010063

APA Style

Faroughi, S. A., Soltanmohammadi, R., Datta, P., Mahjour, S. K., & Faroughi, S. (2024). Physics-Informed Neural Networks with Periodic Activation Functions for Solute Transport in Heterogeneous Porous Media. Mathematics, 12(1), 63. https://doi.org/10.3390/math12010063

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

Physics-Informed Neural Networks with Periodic Activation Functions for Solute Transport in Heterogeneous Porous Media

Abstract

1. Introduction

2. Underlying Physics

3. Methodology

4. Computational Experiments

4.1. Case 1: 1D Solute Transport with Constant Velocity

4.2. Case 2: 2D Solute Transport with Constant Velocity

4.3. Case 3: 2D Solute Transport in Homogeneous Porous Media

4.4. Case 4: 2D Solute Transport in Heterogeneous Porous Media

5. Conclusions

Author Contributions

Funding

Data Availability Statement

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI