Physics-informed neural networks for unsteady incompressible flows with
  time-dependent moving boundaries
Physics-informed neural networks (PINNs) employed in fluid mechanics deal
primarily with stationary boundaries. This hinders the capability to address a
wide range of flow problems involving moving bodies. To this end, we propose a
novel extension, which enables PINNs to solve incompressible flows with
time-dependent moving boundaries. More specifically, we impose Dirichlet
constraints of velocity at the moving interfaces and define new loss functions
for the corresponding training points. Moreover, we refine training points for
flows around the moving boundaries for accuracy. This effectively enforces the
no-slip condition of the moving boundaries. With an initial condition, the
extended PINNs solve unsteady flow problems with time-dependent moving
boundaries and still have the flexibility to leverage partial data to
reconstruct the entire flow field. Therefore, the extended version inherits the
amalgamation of both physics and data from the original PINNs. With a series of
typical flow problems, we demonstrate the effectiveness and accuracy of the
extended PINNs. The proposed concept allows for solving inverse problems as
well, which calls for further investigations.
Yongzheng Zhu, Weizhen Kong, Jian Deng, Xin Bian
Physics-informed neural networks for unsteady incompressible flows with time-dependent moving boundaries

###### Abstract

Physics-informed neural networks (PINNs) employed in fluid mechanics deal primarily with stationary boundaries. This hinders the capability to address a wide range of flow problems involving moving bodies. To this end, we propose a novel extension, which enables PINNs to solve incompressible flows with time-dependent moving boundaries. More specifically, we impose Dirichlet constraints of velocity at the moving interfaces and define new loss functions for the corresponding training points. Moreover, we refine training points for flows around the moving boundaries for accuracy. This effectively enforces the no-slip condition of the moving boundaries. With an initial condition, the extended PINNs solve unsteady flow problems with time-dependent moving boundaries and still have the flexibility to leverage partial data to reconstruct the entire flow field. Therefore, the extended version inherits the amalgamation of both physics and data from the original PINNs. With a series of typical flow problems, we demonstrate the effectiveness and accuracy of the extended PINNs. The proposed concept allows for solving inverse problems as well, which calls for further investigations.

keywords: physics-informed neural networks; unsteady incompressible flows; time-dependent moving boundaries.

## 1 Introduction

In recent years, machine learning has achieved enormous accomplishments in part due to the technological innovation in graph processing units and healthy development of software ecosystems. It finds applications in a variety of disciplines, including physics, linguistics, biology, and many others, as well as in scenarios such as image recognition, natural language processing, and autonomous driving [1]. Inspired by these developments, applied mathematicians have swiftly proposed various ingenious frameworks based on machine learning, aiming to shift the traditional paradigm of computational science and engineering [2; 3]. The novel frameworks may be data-driven, physics-informed, partial differential equations (PDEs)-regularized, in a variational form, optimization of discrete loss and so on [4; 5; 6; 7; 8]. In particular, the physics-informed neural networks (PINNs) have stimulated a wide range of research efforts [3]. The key idea behind PINNs is to incorporate the residual of PDEs into the loss function of the neural network and leverage its powerful nonlinear fitting abilities to approximate the solutions [7]. Given the impressive
performance, researchers have extended PINNs to solve various types of differential equations [9, 10, 11]. Many algorithmic aspects have also been considered to facilitate the training of PINNs to improve the accuracy and/or convergence. These include but are not limited to adaptive activation functions [12], adaptive weights [13], gradient-enhanced approaches [14], artificial-viscosity regularized loss [15], residual-based and/or solution-gradient based adaptive sampling [16, 17], among others. Furthermore, PINNs are also extended to deal with multi-fidelity data and multiscale physics [18, 19]. To tackle various differential equations steadily, an open-source code known as DeepXDE has been provided as a platform for stream-lined training of PINNs [20].

Most known physics laws are described by PDEs and therefore, PINNs provides a new paradigm for solving problems in various fields. Moreover, the solution strategy is universal for both forward and inverse problems, whether with partial data or without data. In fluid mechanics, Navier-Stokes (NS) equations are regarded as the fundamental PDEs. Raissi et al. [7] applied PINNs to infer an unknown coefficient of the NS equations from flow data around a stationary cylinder. Furthermore, they adjusted PINNs for flow visualization, where concentration fields were observed from artwork or medical images to reconstruct the complete flow fields [21]. Rao et al. [22] utilized PINNs as a direct numerical simulation (DNS) tool to study flow around a stationary cylinder at low Reynolds number. Jin et al. [23] evaluated the accuracy, convergence, and computational cost of PINNs on laminar flows such as flow around a cylinder and turbulent flows in pipeline using DNS dataset as benchmark. Wang et al. [24] proposed an algorithm to enable PINNs to deal with 1D two-phase Stefan problems with free boundaries. Cai et al. [25] explored PINNs for heat transfer problems, in which they utilized sparse measuring points of temperature to infer both temperature and velocity fields of the entire domain. Despite so many successful applications of PINNs in fluid mechanics, flows with moving boundaries are rarely considered, which are frequently encountered in nature and industry, such as flow over a moving blunt body, bird/insect hovering and fish swimming [26, 27, 28, 29, 30, 31]. One exception is the work of Raissi et al. [32], where they adopted a transformed coordinate attached to the cylinder to resolve the flow for vortex-induced vibration. The coordinate transformation is elegant, but it is not suitable for general moving boundaries such as multiple blunt bodies in flows. In contrast, computational fluid dynamics (CFD) methods have been employed routinely to solve flow problems with moving boundaries [33, 34, 35], although they may require complex meshing techniques such as morphing mesh [36] and overset mesh [37], or interpolation techniques such as the immersed boundary method (IBM) [38]. To bridge this gap and allow PINNs to handle flow problems with moving boundaries, we introduce new loss functions proper for training points at moving interfaces and around moving boundaries in the entire temporal domain. Therefore, we are still able to leverage the power and flexibility of the classical PINNs to infer all flow fields, whether with partial data or without any exogenous data. Interactions between fluid and a moving body may be categorized into two types: one-way coupling, where the motion of the body is prescribed and its influence on the fluid dynamics needs to be computed; two-way coupling, where interface between the fluid and the body is part of the solution itself and dynamics of both parts must be computed as well. In this work, we focus on the first type, as in reality it is rather straightforward to extract the trajectories of blunt bodies and reconstruction of the flow fields around moving bodies often remain obscure.

This rest of the paper is organized as follows. In Section 2, we review PINNs briefly and present the no-slip velocity condition for time-dependent moving bodies in detail in Section 3. Direct numerical simulation results of PINNs are validated by three examples in Section 4. In Section 5, we further reconstruct the whole flow fields with partial data in three representative cases. Discussions and conclusions are made in Section 6. Auxiliary data from CFD methods are presented in the Appendix.
## 2 Physics-informed neural networks

Consider the parametric nonlinear PDEs in a general form expressed as follows:

\[f\left(\mathbf{x};\frac{\partial s}{\partial x_{1}},\ldots,\frac{\partial s}{ \partial x_{d}};\frac{\partial^{2}s}{\partial x_{1}\partial x_{1}},\ldots, \frac{\partial^{2}s}{\partial x_{1}\partial x_{d}};\ldots;\lambda\right)=0, \quad\mathbf{x}\in\Omega, \tag{1}\]

with the boundary condition and initial condition:

\[\mathcal{B}(s,\mathbf{x})=0,\quad\mathbf{x}\in\partial\Omega, \tag{2}\] \[\mathcal{I}(s,\mathbf{x})=0,\quad\mathbf{x}\in\Gamma_{i}. \tag{3}\]

Here \(\mathbf{x}=[x_{1},x_{2},\ldots,x_{d}]\) are the independent variables, \(f\) denotes the linear or nonlinear differential operators, \(s\) is the solution, and \(\lambda=[\lambda_{1},\lambda_{2},\ldots]\) are the parameters for combining each component. \(\Omega\) and \(\partial\Omega\) represent the computational domain and the boundaries, respectively. \(\Gamma_{i}\) stands for the space of \(\mathbf{x}\) at the initial snapshot.

Within the framework of PINNs [7], a fully connected neural network (FCNN) is exploited as a highly nonlinear function \(\hat{s}(\mathbf{x})\) to approximate the solution \(s\) of the PDEs. It is composed of an input layer, multiple hidden layers, and an output layer, while each layer contains several neurons with weights, biases, and nonlinear activation function \(\sigma(\cdot)\). Let \(\mathbf{x}\) be the input and the implicit variable of the \(i\)th hidden layer be \(\mathbf{y}^{i}\), then a FCNN with \(L\) layers can be expressed as:

\[\begin{cases}\mathbf{y}^{0}=\mathbf{x},\\ \mathbf{y}^{i}=\sigma\left(\mathbf{W}^{i}\mathbf{y}^{i-1}+\mathbf{b}^{i}\right),\quad 1\leq i \leq L-1,\\ \mathbf{y}^{i}=\mathbf{W}^{i}\mathbf{y}^{i-1}+\mathbf{b}^{i},\quad i=L,\end{cases} \tag{4}\]

where \(\mathbf{W}^{i}\) denotes the weight matrix and \(\mathbf{b}^{i}\) is the bias vector to be trained at \(i\)th layer, respectively. \(\sigma(\cdot)\) functions element-wise.

Since \(\hat{s}\) is the approximate solution, each component in the PDEs can be obtained by taking derivatives of \(\hat{s}\) with respect to \(x_{i}\) one or more times by the automatic differentiation (AD) techniques combined with the chain rule [39]. Subsequently, we define a composite loss function as the sum of the residuals of the equation, boundary condition, initial condition, and labeled data:

\[\mathcal{L}(\mathcal{T})=w_{f}\mathcal{L}_{f}\left(\mathcal{T}_{f}\right)+w_{ b}\mathcal{L}_{b}\left(\mathcal{T}_{b}\right)+w_{i}\mathcal{L}_{i}\left( \mathcal{T}_{i}\right)+w_{d}\mathcal{L}_{data}\left(\mathcal{T}_{data}\right), \tag{5}\]

where \(w_{f}\), \(w_{b}\), \(w_{i}\) and \(w_{d}\) are the corresponding weights and \(\mathcal{T}_{f}\), \(\mathcal{T}_{b}\), \(\mathcal{T}_{i}\) and \(\mathcal{T}_{data}\) are the corresponding sets of training points, respectively.

With the target of minimizing the loss function, the neural network is trained by optimizing the weights and biases through the back-propagation algorithm with an optimizer such as Adam [40]. The flexibility of PINNs lies in the capability of switching between supervised, weakly supervised, and unsupervised learning approaches. For a supervised learning or data-driven approach, the loss function is guided only by a labeled dataset, that is, \(\mathcal{T}_{data}\), and PINNs degenerates to be FCNNs. For an unsupervised learning, the loss function is defined by the equation, the boundary condition, and the initial condition. The training dataset \(\mathcal{T}_{f}\) corresponds to collocation points, where solutions need to be inferred, and they are constrained by the points within \(\mathcal{T}_{b}\) and \(\mathcal{T}_{i}\). More frequently, PINNs are employed with all the four types of loss functions of training points.
## 3 PINNs for incompressible flows with time-dependent moving boundaries

### Governing equations

The NS equations are expressed as continuity and momentum equations in two dimensional coordinates \(\mathbf{x}=(x,y)\) as follows:

\[\nabla\cdot\mathbf{u} =0,\] \[\frac{\partial\mathbf{u}}{\partial t}+(\mathbf{u}\cdot\nabla)\mathbf{u} =-\nabla p+\frac{1}{Re}\nabla^{2}\mathbf{u}, \tag{6}\]

where \(t\) is the time, \(\mathbf{u}=(u,v)\) is the velocity vector, \(p\) is the pressure, and \(Re=UD/v\) is the Reynolds number defined by the reference velocity \(U\), the characteristic length \(D\), and the kinematic viscosity \(v\).

If there are solid or truncated boundaries at stationary, one or more constraints for velocity and/or pressure may be specified as follows:

\[\mathbf{u}=\mathbf{u}_{\Gamma}(\mathbf{x}),\quad\mathbf{x}\in\Gamma_{D}, \tag{7}\] \[p=p_{\Gamma}(\mathbf{x}),\quad\mathbf{x}\in\Gamma_{D},\] (8) \[\frac{\partial\mathbf{u}}{\partial\mathbf{n}}=0,\quad\mathbf{x}\in\Gamma_{N},\] (9) \[\frac{\partial p}{\partial\mathbf{n}}=0,\quad\mathbf{x}\in\Gamma_{N}, \tag{10}\]

where \(\mathbf{n}\) is the normal vector at the boundary; \(\Gamma_{D}\) and \(\Gamma_{N}\) denote the Dirichlet and Neumann boundaries, respectively.

### Boundary conditions of moving bodies

Firstly, we take a cylinder that travels through the fluid constantly with time, as an example. Other moving solid bodies can be accomplished similarly. As the interface between the fluid and solid is time-dependent, we may "excavate a tunnel" occupied by the cylinder in the spatial-temporal domain, as illustrated by Fig. 1(a). We denote the fluid domain outside of the tunnel as \(\Omega\). Accordingly, an interface \(\partial\Omega\) is

Figure 1: Schematic of a moving cylinder in spatial-temporal domain and its corresponding training points. (a) Training points are randomly distributed in \(\Omega\), which excludes the ”tunnel” occupied by the travelling cylinder with \(\partial\Omega\) as the interface. (b) Distribution of training points in \(x-y\) plane at time \(t=t_{1}\), with finer points sampled at the interface and around the cylinder. (c) Aerial view in the y-direction of all training points sampled at the moving interface.
formed as the cylinder travels through the fluid and it is given as a function of time:

\[\partial\Omega:=f(\mathbf{\delta}(t),\alpha(t),t), \tag{11}\]

where \(\mathbf{\delta}=(\zeta,\eta)\) is the spatial coordinate for the center of the cylinder. Moreover, \(\alpha\) is its rotation angle vector, which shall be useful for a general-shaped moving body. The Dirichlet velocity constraint according to the no-slip condition is as follows:

\[\mathbf{u}=\frac{\partial\mathbf{\delta}}{\partial t}+\mathbf{r}\times\frac{\partial \alpha}{\partial t},\quad(\mathbf{x},t)\in\partial\Omega, \tag{12}\]

where \((\mathbf{x},t)=(x,y,t)\) is a spatial-temporal coordinate at the interface and \(\mathbf{r}\) is the corresponding radial vector to the center of the cylinder.

Traditional CFD methods may address this relatively complex issue, for example, by regenerating computational mesh frequently. In contrast, PINNs does not have a time-stepping, but instead optimizes the flow prediction in a continuous spatial-temporal domain via minimizing the composite loss of training points. It is simple to see that _the time-dependent moving boundary problem in the two-dimensional coordinate of space corresponds to a stationary problem in the three dimensional coordinate of space and time._ In this way, one or more moving boundaries can be dealt with steadily. Therefore, we firstly generate sufficient training points randomly in the entire \(\Omega\), as illustrated in Fig. 1(a). Furthermore, to take into account the stiff flow variation in the vicinity of the moving boundary, we distribute finer training points around the moving boundary, as shown in Fig. 1(b). In Fig. 1(c), an aerial view of the training points that are sampled at the moving interface, are looked upon in the \(y\) direction. Meanwhile, additional training points at the interface \(\partial\Omega\) are imposed with Dirichlet velocity constraint.

Theoretically, this strategy can be extended to three spatial dimensions as well, which just needs to be elevated into four dimensional coordinate of space and time. Next we shall introduce the technical details of PINNs' architecture and elaborate each component of the composite loss function.

### Each loss function and technical details for PINNs

Based on the principles of PINNs in Section 2, we take spatial and temporal coordinates \((x,y,t)\) as inputs and \((u,v,p)\) as outputs of the neural networks to predict the unsteady velocity and pressure fields. Fig. 2 illustrates a schematic of PINNs for solving unsteady flow problems with moving boundaries. The derivatives of \(u\), \(v\), and \(p\) with respect to the inputs are calculated using the AD. The composite loss of PINNs is repeated as follows:

\[\mathcal{L}(\theta)=w_{f}\mathcal{L}_{f}(\theta)+w_{b}^{s}\mathcal{L}_{b}^{s} (\theta)+w_{b}^{m}\mathcal{L}_{b}^{m}(\theta)+w_{i}\mathcal{L}_{i}(\theta)+w_ {data}\mathcal{L}_{data}(\theta). \tag{13}\]

The difference from Eq. (5) is that we replace \(w_{b}\mathcal{L}_{b}\) with \(w_{b}^{s}\mathcal{L}_{b}^{s}+w_{b}^{m}\mathcal{L}_{b}^{m}\) to differentiate the sub-loss for stationary and moving boundaries, respectively. Morever, we adopt the Glorot scheme to randomly initialize all weights and biases denoted by \(\theta\)[41], which are to be optimized. The definitions for all components in
Eq. (13) are expressed in mean squared errors (MSE) as:

\[\mathcal{L}_{f} =\frac{1}{N_{f}}\sum_{(\mathbf{x},t)\in\Omega}\|\nabla\cdot\mathbf{u}\|_{2}^ {2}+\frac{1}{N_{f}}\sum_{(\mathbf{x},t)\in\Omega}\bigg{\|}\frac{\partial\mathbf{u}}{ \partial t}+(\mathbf{u}\cdot\nabla)\mathbf{u}+\nabla p-\frac{1}{Re}\nabla^{2}\mathbf{u} \bigg{\|}_{2}^{2}, \tag{14}\] \[\mathcal{L}_{b}^{s} =\frac{1}{N_{b}^{s}}\sum_{\mathbf{x}\in\Gamma_{D}}\Big{(}\|\mathbf{u}-\bm {u}_{\Gamma}(\mathbf{x})\|_{2}^{2}+\|p-p_{\Gamma}(\mathbf{x})\|_{2}^{2}\Big{)}+\frac{1 }{N_{b}^{s}}\sum_{x\in\Gamma_{N}}\bigg{(}\bigg{\|}\frac{\partial\mathbf{u}}{ \partial\mathbf{u}}\bigg{\|}_{2}^{2}+\bigg{\|}\frac{\partial p}{\partial\mathbf{u}} \bigg{\|}_{2}^{2}\bigg{)},\] (15) \[\mathcal{L}_{b}^{m} =\frac{1}{N_{b}^{m}}\sum_{(\mathbf{x},t)\in\partial\Omega}\bigg{(} \bigg{\|}\mathbf{u}-\frac{\partial\mathbf{\delta}}{\partial t}-\mathbf{r}\times\frac{ \partial\alpha}{\partial t}\bigg{\|}_{2}^{2}\bigg{)},\] (16) \[\mathcal{L}_{i} =\frac{1}{N_{i}}\sum_{(\mathbf{x},t)\in\Gamma_{i}}\Big{(}\|\mathbf{u}-\bm {u}_{i}(t)\|_{2}^{2}+\|p-p_{i}(t)\|_{2}^{2}\Big{)},\] (17) \[\mathcal{L}_{data} =\frac{1}{N_{data}}\sum_{n=1}^{N_{data}}\bigg{(}\Big{|}\mathbf{u}( \mathbf{x}_{u}^{n},t_{u}^{n})-\mathbf{u}_{data}^{n}\big{|}^{2}+\big{|}p(\mathbf{x}_{p}^{n},t_{p}^{n})-p_{data}^{n}\big{|}^{2}\bigg{)}, \tag{18}\]

where \(N_{f}\), \(N_{b}^{s}\), \(N_{b}^{m}\), \(N_{i}\) and \(N_{data}\) denote the number of training points corresponding to each loss term. \(\Omega\), \(\Gamma_{D}/\Gamma_{N}\), \(\partial\Omega\) and \(\Gamma_{i}\) indicate the training domain, the stationary boundaries, the moving boundaries and the spatial-temporal coordinates for initial condition, respectively. Each component of the loss function is also illustrated in Fig. 2 for a clear overview.

The numbers of hidden layers and neurons in each layer of the neural network are usually selected to meet specific needs of the problem. After trial-and-error for the sensitivity of the results, we employ 8 hidden layers and each layer contains 40 neurons universally. We choose the continuous and differentiable hyperbolic tangent function _tanh_ as the activation function \(\sigma(\cdot)\)[1]. Typically we minimize the total loss

Figure 2: Schematic of PINNs for solving unsteady flow problems with moving boundaries: four types of loss functions and the one for boundary condition consists of stationary and moving parts.
function using the Adam optimizer [40] with an adaptive stochastic objective function combined with the L-BFGS optimizer [42] based on the Quasi-Newton method. The Adam optimizer works with a decreasing learning rate schedule. More specifically, we train the model with a learning rate of \(10^{-3}\) for the first 100,000 epochs, \(5\times 10^{-4}\) for the subsequent 30,000 epochs, and \(10^{-4}\) for the last 30,000 epochs. The L-BFGS optimizer follows for about 50,000 epochs to further diminish the residuals. We obtain approximate solutions when the loss function reaches a sufficiently low level. To evaluate the accuracy of the prediction, we select the relative \(L_{2}\) error as a metric:

\[\varepsilon_{V}=\frac{\|V-V^{*}\|_{2}}{\|V^{*}\|_{2}}, \tag{19}\]

where \(V\) represents the predicted solutions \((u,v,p)\), and \(V^{*}\) denotes the corresponding reference solutions.

All the implementations are based on the framework of DeepXDE [20], which further employs TensorFlow [43] in the backend. High-fidelity DNS data computed by OpenFOAM (FVM-based) [44] are taken as reference in Section 4 and also leveraged as partial data for flow reconstruction in Section 5. The mesh designed for each case is presented in A.

## 4 PINNs as a direct numerical simulation solver

In this section, we employ PINNs only with initial data and validate the proposed extension as direct numerical simulation for three challenging flow problems, which involve moving boundaries.

### In-line oscillating cylinder

The interaction of an oscillating cylinder with a fluid at rest is of particular interest in fields such as marine engineering. It is a classical case of fluid-structure interaction problem and is well-documented in literature [45; 46]. The flow manifests itself by a complex vortex-structure interaction phenomena. Two key dimensionless parameters are Reynolds number \(Re\) and Keulegan-Carpenter number \(KC\), defined as:

\[Re=\frac{U_{max}D}{\nu},\quad KC=\frac{U_{max}}{fD}, \tag{20}\]

where \(U_{max}\) is the maximum velocity of the moving cylinder, \(D\) is the cylinder diameter, \(\nu\) is the kinematic viscosity of the fluid, and \(f\) is the characteristic frequency of the oscillation. The geometry of this flow is shown in Fig. 3(a). The translational motion of the cylinder is prescribed by a simple harmonic oscillation:

\[x(t)=-A\cdot\sin(2\pi f\cdot t), \tag{21}\]

where \(x\) denotes the location of the cylinder's center, and \(A\) is the amplitude of the oscillation. Thus, the Keulegan-Carpenter number can also be written as:

\[KC=\frac{2\pi A}{D}. \tag{22}\]

Referring to the experiments of Dutsch et al. [45] and the numerical simulations of Liu et al. [46], the Reynolds and Keulegan-Carpenter numbers are set to \(Re=125\) and \(KC=5\), respectively. Accordingly, \(\rho=1\), \(\nu=0.008\), \(D=1\), \(U_{max}=1\) and oscillation period \(T=1/f=5\). The far-field velocities are subjected to Neumann velocity conditions, as illustrated in Fig. 3(a). Furthermore, Dirichlet pressure condition is applied on the left and right sides, while Neumann pressure condition is imposed on the top and bottom sides. To capture the cylinder oscillation accurately, we enforce the Dirichlet velocity condition
on the cylinder boundary, wherein \(u(t)=-2\pi f\cdot A\cdot\cos(2\pi f\cdot t)\) and \(v(t)=0\). As for the initial conditions, we adopt the high-fidelity velocity and pressure data obtained from the FVM at \(t^{\prime}=15\) (phase \(0^{\circ}\)) as the starting point \(t=0\) of PINNs.

We first define the spatial computational domain as a rectangular region \([x_{min},x_{max}]\times[y_{min},y_{max}]=[-5,5]\times[-5,5]\). For the time domain, a single oscillation period from \(0.0\) to \(5.0\) is partitioned into two separate half-periods, \([0.0,2.5]\) and \([2.5,5.0]\), to be solved individually. Next, we excavate a tunnel formed by the oscillation of the cylinder in the spatial-temporal domain. A visual representation of the tunnel corresponding to one oscillation period is shown in Fig. 3(b). Note that the computational domain excludes the interior of the moving cylinder. Therefore, only training points obtained from locations in the domain other than the interior of the tunnel are considered. A number of \(20,120\) data points from the FVM containing the velocity and pressure data of the flow field are used as the initial conditions. We randomly sample \(100,000\) training points for the fluid, \(20,000\) training points at the cylinder, and \(10,000\) training points at the truncated rectangular boundary. We further enhance the sampling density by adding an additional \(50,000\) points of finer resolution in close proximity to the surface of the cylinder, as well as \(20,000\) finer points within a rectangular region \([-2,-1]\times[2,1]\). Here, a snapshot of the distribution of the points sampled within the temporal neighborhood \(\epsilon=0.1\) at \(t\)=1.25 is shown in Fig. 3(c). In Fig. 3(d), an aerial view of the training points that are sampled at the moving interface in two half-periods, are looked upon in the \(y\) direction. The convergence of the total loss, PDEs loss, moving boundary loss, stationary boundary loss, and initial con

Figure 3: Problem setup and training strategy for an in-line oscillating cylinder in fluid at rest. (a) The geometry of the computational domain and boundary conditions. (b) Training points at the moving interface, around the oscillating cylinder boundary and within the computational domain in space and time. (c) Snapshots of spatial training points sampled within the temporal neighborhood \(\epsilon=0.1\) at \(t\)=1.25. The points are sampled with finer resolution in the region close to the cylinder and within the in-line oscillation range. (d) Aerial view of the training points sampled at the moving interface in the two half-period time domains from the \(y\)-direction, respectively. (e) The training losses versus the number of optimization epochs in the first (left) and second (right) half periods.
dition loss, are depicted in Fig. 3(e), where we observe an effective reduction of the residuals by the Adam optimizer followed by the L-BFGS optimizer.

In Fig. 4, the predicted vorticity and pressure contours and isolines of the flow at four typical phases (\(0^{\circ}\), \(96^{\circ}\), \(192^{\circ}\), and \(288^{\circ}\)) by PINNs are presented, exhibiting excellent consistency with the results of FVM. This indicates that the method is capable of capturing the symmetrical vortex shedding characteristics and pressure distribution during the oscillation of the cylinder. The first row of Fig. 5 shows the velocity \(u\), \(v\), and pressure \(p\) at phase \(180^{\circ}\) predicted by PINNs, where the maximum velocity of the circular cylinder is located. These results are compared to those of FVM in the second row, and meanwhile point-wise absolute errors are presented in the third row of Fig. 5, where the errors are well below 5% overall. For selected time instants, we present relative \(L_{2}\) errors for the whole field in Fig. 6. We note that the grid nodes in the FVM

Figure 4: Results of PINNs and FVM for vorticity and pressure contours, and isolines in several typical phases. Four distinct phases and the corresponding time: (a) \(0^{\circ}\) (\(t\)=0); (b) \(94^{\circ}\) (\(t\)=1.3); (c) \(194^{\circ}\) (\(t\)=2.7); (d) \(288^{\circ}\) (\(t\)=4). Solid and dotted lines denote positive and negative contours, respectively.
are employed as the coordinate points for PINNs' predictions and validations.

The velocity distributions in the \(y\)-directional cross-section are further shown in Fig. 7 for \(u\) and \(v\) in three phases (180\({}^{\circ}\), 210\({}^{\circ}\), and 330\({}^{\circ}\)) of the oscillating cylinder, and compared with FVM results. The comparison demonstrates that the proposed model captures the finer flow distribution with great accuracy.

### Multiple cylinders translating along a circle

We further consider the translational motion of multiple bodies in fluid to evaluate the capability and universality of the proposed extension on PINNs. The key objective is to examine the relation between

Figure 5: The comparison between PINNs and FVM for the velocity \(u\), \(v\) and pressure \(p\) at phase 180\({}^{\circ}\), where the velocity of the cylinder is maximum. The point-wise absolute errors are depicted in the bottom panels.

Figure 6: Relative \(L_{2}\) errors of PINNs for the in-line oscillating cylinder in fluid at rest: the predicted velocity \(u\), \(v\) and pressure \(p\).
the number of moving boundaries and the accuracy of model training in the simulation of multiple moving bodies. Inspired by the numerical simulations of Xu et al. [47], we consider 2, 3 and 4 cylinders moving clockwise in fluid along a circle of radius \(R\) as an orbit. These cylinders start from rest in a sinusoidal pattern and gradually accelerate to a steady maximum velocity. The geometries of the setups, featuring 2, 3, and 4 cylinders at the maximum velocity, are presented in Fig. 8(a).

The key dimensionless parameter in this flow is the Reynolds number \(Re\), defined as:

\[Re=\frac{U_{max}D}{\nu}, \tag{23}\]

where \(U_{max}\) is the maximum velocity of the moving cylinders, \(D\) is the diameter of the cylinders, \(\nu\) is the kinematic viscosity of the fluid. The motion of the cylinders is given by:

\[\begin{cases}x=R\cdot\cos(\theta(t)),\\ y=R\cdot\sin(\theta(t)),\end{cases} \tag{24}\]

where \(\theta(t)\) represents the angle of the center of the cylinder in a polar coordinate, which can be expressed

Figure 7: PINNs results for velocity contours and selected profiles for \(u\) and \(v\) along the \(y\)-direction are compared with those of FVM at three typical phases of the in-line moving cylinder. Solid and dotted lines denote PINNs and FVM results along four profiles indicated with colors: red: \(x=-0.6D\); green: \(x=0D\); blue: \(x=0.6D\); orange: \(x=1.2D\).
as:

\[\theta(t)=\begin{cases}\varphi+\dfrac{A}{R}\cdot\cos(2\pi f\cdot t), \quad 0\leq t\leq 2.5,\\ \varphi-\dfrac{U_{\max}}{R}\cdot(t-2.5),\quad t\geq 2.5,\end{cases} \tag{25}\]

where \(A\) and \(f\) represent the amplitude and frequency, respectively, of the cosine acceleration, \(R\) is the radius of the circle along which the cylinders are translating. The variable \(\varphi\) represents the phase angle, given in polar coordinates, at which the cylinders reach their maximum linear velocity. For the two-cylinder system, the phase angles of cylinders A and B are \(\varphi_{A}=\pi/2\) and \(\varphi_{B}=-\pi/2\), respectively; For the three-cylinder system, the phase angles of cylinders A, B and C are \(\varphi_{A}=\pi/2\), \(\varphi_{B}=-\pi/6\) and \(\varphi_{C}=-5\pi/6\), respectively; For the four-cylinder system, the phase angles of cylinders A, B, C and D are \(\varphi_{A}=\pi/2\), \(\varphi_{B}=0\), \(\varphi_{C}=-\pi/2\) and \(\varphi_{D}=-\pi\), respectively.

In all three systems, the Reynolds number is set to \(Re=100\). The parameters are set as follows, \(\rho=1\), \(\nu=0.01\), \(D=1\), \(R=2\), \(U_{max}=1\) and acceleration parameters \(A=1.592\), \(f=0.1\). The far-field boundaries are subjected to Neumann velocity conditions and constant pressure condition \(p=1\) on four sides, as illustrated in Fig. 8(a). To accurately describe the motion of the cylinder, Dirichlet velocity conditions are imposed on each boundary of the cylinder, where \(u(t)=-R\cdot\sin(\theta(t))\cdot\frac{d\theta(t)}{dt}\) and \(v(t)=R\cdot\cos(\theta(t))\cdot\frac{d\theta(t)}{dt}\). The initial condition \(t=0\) for these systems is that both the fluid and the cylinders are stationary, i.e. \(u=v=p=0\).

Figure 8: Problem setup and training points for multiple cylinders translating along a circle. (a) The geometry of the computational domain and boundary conditions. The left, middle, and right panels depict the geometric configurations of 2, 3, and 4 cylinders, respectively. (b) Scattered training points at four moving interfaces, around four moving boundaries, and in the fluid domain in space and time. (c) Snapshots of spatial training points sampled within a temporal neighborhood of \(\epsilon=0.2\) at \(t=0\). (d) Aerial view of the training points sampled at the interfaces of four moving cylinders from the \(y\) direction.
The spatial computational domain is defined as a rectangular region \([-4,4]\times[-4,4]\). The time interval [0, 5] is partitioned into two stages: [0, 2.5] for acceleration of the cylinders, and [2.5, 5] for maintaining maximum velocity through cylinder translation. And the two stages are treated as separate intervals and solved independently. Similarly, we excavate the tunnels that are formed as a result of the movements of multiple cylinders throughout the space-time domain. A visual representation of these moving boundary tunnels for the four-cylinder system is displayed in Fig. 8(b). We randomly sample 70,000 points in the fluid domain, 10000 initial points, 10,000 rectangular boundary points, 15,000 points at each cylinder interface. We have enhanced the sampling density by adding an additional 40,000 points of higher resolution, in close proximity to the boundary surface of each cylinder. In addition, 50,000 points are sampled within the orbit region \((R-D/2<r<R-D/2)\) of the cylinder motion. We illustrate the spatial distribution of the training points within the temporal neighborhood where \(\epsilon=0.2\) at \(t\)=0 for a system composed of four cylinders, as shown in Fig. 8(c). Additionally, an aerial view of the training points sampled at the interfaces of four moving cylinders from the \(y\) direction is depicted in Fig. 8(d).

Fig. 9 illustrates the convergence of the total loss during training for the three multi-cylinder systems, as well as the individual components of the loss, namely PDEs loss, moving boundary loss, stationary boundary loss, and initial boundary loss. The predicted vorticity and pressure contours of the flow for three multi-cylinder systems are presented in Fig. 10. The results demonstrate remarkable consistency with the FVM simulations, thereby confirming PINNs' effectiveness in solving multi-body flow problems. Fig. 11 provides the point-wise relative \(L_{2}\) errors of predicted velocity and pressure solutions in three multi-cylinder systems. It should be noted that the grid nodes in the FVM have been selected as the coordinate points for PINNs' validations. The results from the relative errors show that the model is able to solve the flow

Figure 9: The training losses versus the number of optimization epochs for (a) two-cylinder system, (b) three-cylinder system, and (c) four-cylinder system. The training losses for the acceleration and steady velocity stages are depicted in the top and bottom panels, separately.
problem containing multiple rigid bodies quite accurately. As the number of moving bodies increases, there is a slight increment of the prediction error on average. Nevertheless, the magnitude of the error remains within an acceptable range. We postulate that this error increases due to the increased complexity of the flow as a result of hydrodynamic interactions from multiple moving objects. To address this, it may become imperative to supplement PINNs with additional training points and/or explore a new network architecture to capture finer details. The primary objective of this study is to interrogate the feasibility of the new extension of PINNs to solve multi-body flow problems. Further research endeavors may potentially verify these hypotheses.
### Flow around a flapping wing

We study a flapping wing that undergoes both translation and rotation. This motion pattern of the wing is similar to that of insect or bird wings during hovering. Here, we provide only the initial condition data of the flow field and to predict the entire flow field by PINNs. We consider a rigid wing with an elliptical shape, characterized by a chord length \(C\) and an aspect ratio \(E\). Fig. 12(a) depicts a sinusoidal motion of the wing cross-section in the chord direction, which is prescribed as follows:

\[\left\{\begin{aligned} x(t)&=\frac{A}{2}\left[\cos \left(\frac{2\pi t}{T}\right)+1\right]\cos\beta,\\ y(t)&=\frac{A}{2}\left[\cos\left(\frac{2\pi t}{T} \right)+1\right]\sin\beta,\\ \alpha(t)&=\alpha_{0}\left[1-\sin\left(\frac{2\pi t }{T}+\varphi\right)\right].\end{aligned}\right. \tag{26}\]

Figure 11: Relative \(L_{2}\) errors of velocity \(u\), \(v\) and pressure \(p\) by PINNs for multiple cylinders translating along a circle.

Figure 10: PINNs-predicted and FVM results for vorticity and pressure contours in three multi-cylinder flow problems.
Here \(x\) and \(y\) denote the coordinates of the wing's center, \(\alpha\) represents the attack angle, \(T\) is the flapping period, \(A\) is the flapping displacement, \(\beta\) is the angle of inclination of the stroke plane, and \(\varphi\) is the phase angle. The velocity amplitude of the flapping wing motion is, therefore, \(U_{max}=\pi\cdot A/T\). The key dimensionless parameter in this flow is the Reynolds number \(Re\) defined as

\[Re=\frac{U_{max}C}{\nu}=\frac{\pi AC}{\nu T}, \tag{27}\]

where \(\nu\) is the kinematic viscosity of the fluid.

The Reynolds number for is selected as \(Re\)=157, following the simulation undertaken by Wang et al. [26]. The various simulation parameters are set as: \(C=1\), \(E=4\), \(A=2.5C\), \(\alpha_{0}=\pi/4\), \(\beta=\pi/3\), \(\varphi=0\), \(T=5,\rho=1\), and \(\nu=0.01\). The spatial computational domain is defined as a rectangular region \([-3,4]\times[-4,4]\). We employ Neumann velocity conditions and a constant pressure condition \(p=1\) at all far-field boundaries, as illustrated in Fig. 12(a). To precisely describe the wing's motion, the boundary of the wing is subject to Dirichlet velocity conditions, which can be computed using Eq. (12). As for the initial conditions, we adopt the high-fidelity velocity and pressure data obtained from the FVM at \(t^{\prime}=20\) as the starting point \(t=20\) of PINNs, and we solve the flow within the time domain T = [20.0, 20.6]. We use the time sequence scheme with a subdomain size of \(\Delta T=0.2\). Fig. 12(b) displays a visual representation of the tunnel of a wing in space and time. 100,000 training points are randomly sampled in the domain composed of a rectangular region of width \(7C\) and height \(8C\) and the time subdomain. 30,000 and 20,000 points are

Figure 12: Problem setup and training strategy of flow around a flapping wing. (a) The geometry of the computational domain and boundary conditions. The solid ellipses indicate the down-stroke phase while the dashed ellipses represent the upstroke phase. (b) Visual representation of scattered training points at interface, around boundary and in the fluid domain. (c) Snapshots of spatial training points sampled within the temporal neighborhood \(\epsilon=0.04\) at \(t\)=20. The training points are refined with higher resolution in the region close to the wing and within the flapping range. (d) Aerial view of all the points sampled at the wing surface from the \(x\) (left) and \(y\) (right) directions. (e) The training losses versus the number of optimization epochs.
sampled on the wing surface and rectangular boundaries respectively. 30,000 finer sampled points are added to a rectangular region \([-1.3,1.3]\times[2.5,3.5]\). In addition, 80,000 finer sampled points are added to a circular region with a radius of \(1.5C\) centered on the wing. There are about 8,967 data points from FVM used as initial conditions. Fig. 12(c) presents a snapshot of the distribution of sampled points except for initial data points within the temporal neighborhood \(\epsilon=0.04\) at \(t\)=20.0. We additionally present, in Fig. 12(d), the aerial views of the training points sampled at the wing surface, are looked upon in the \(x\) and \(y\) direction. The convergence of various training losses is depicted in Fig. 12(e).

Fig. 13 illustrates PINNs' results of velocity \(u\), \(v\) and pressure \(p\) contours in three time frames \(t=20.2\), 20.4, 20.6 within a flapping period. Fig. 14 presents the snapshots of the vorticity contours from \(t=20.0\) to 20.6. The model captures the typical and subtle variations in the vortex structure around the wing as it rotates and translates. Fig. 15 provides the point-wise relative \(L_{2}\) errors at each time from \(t=20.0\) to 20.6. The comparison between the predictions by PINNs and results obtained by the FVM demonstrates the capability of the former to accurately predict changes in the flow field based on the boundary conditions of the flapping wing motion.

Fig. 13: PINNs’ results for flow around a flapping wing: velocity \(u\), \(v\) and pressure \(p\) contours at several time frames: (a) \(t\)=20.2; (b) \(t\)=20.4; (c) \(t\)=20.6.
## 5 Reconstructing the flow fields via partial data

In this section, we evaluate the capability of the proposed extension in reconstructing the entire flow field when only partial data is available. This is also one of the best properties of PINNs, which is unavailable in other computational frameworks.

### Transversely oscillating cylinder in steady flow

We consider a transversely oscillating cylinder immersed in a steady flow. In this particular instance, the existence of asymmetrical vortex shedding introduces an elevated level of complexity in comparison to the preceding case lacking a steady flow, thus presenting an additional challenge for PINNs. The two key dimensionless parameters are Reynolds number \(Re\) and Strouhal number \(St\), defined as:

\[Re=\frac{U_{\infty}D}{\nu},\quad St=\frac{fD}{U_{\infty}}, \tag{28}\]

where \(U_{\infty}\) is the velocity of uniform stream, \(D\) is the diameter of the cylinder, \(\nu\) is the kinematic viscosity of the fluid, and \(f\) is the vortex shedding frequency. The geometry is illustrated in Fig. 16(a), where setups

Figure 14: PINNs’ results for flow around a flapping wing: vorticity contours at several time frames: (a) \(t\)=20.1; (b) \(t\)=20.2; (c) \(t\)=20.3; (d) \(t\)=20.4; (e) \(t\)=20.5; (f) \(t\)=20.6.

Figure 15: Relative \(L_{2}\) error of velocity \(u\), \(v\) and pressure \(p\) by PINNs for flow around a flapping wing.
for inlet, outlet, and far-field boundaries are also given. A prescribed motion of the cylinder is expressed as:

\[y(t)=-A\cdot\cos(2\pi f_{c}\cdot t), \tag{29}\]

where \(y\) denotes the cross stream location of the cylinder's center, \(A\) and \(f_{c}\) are the amplitude and characteristic frequency of the oscillation, respectively.

The Reynolds number is set to \(Re=185\), consistent with the numerical simulations reported by Guilmineau et al. [48]. Accordingly, the simulation parameters are set as follows: \(D=1\), the fluid density \(\rho=1\), \(\nu=0.01\), \(U_{\infty}=1.85\), the oscillation amplitude \(A/D=0.2\), and the frequency ratio \(f_{c}/f_{0}=0.8\), where \(f_{0}=St\cdot U_{\infty}/D\) is the vortex shedding frequency of the stationary cylinder at \(St=0.195\). Neumann velocity and pressure boundary conditions are used at the far-field boundaries on the top and bottom sides. At the inlet boundary, a steady horizontal free stream is imposed with the velocity Dirichlet boundary condition \(u=U_{\infty}\) and \(v=0\). And the zero pressure condition and Neumann velocity conditions are applied to the outlet. Dirichlet velocity conditions are imposed on the cylinder boundary to capture the no-slip condition due to oscillation. Specifically, we set \(u(t)=0\) and \(v(t)=2\pi f_{c}\cdot A\cdot\sin(2\pi f_{c}\cdot t)\). We take the velocity and pressure data on 21,535 points obtained from the FVM at \(t^{\prime}=2.0\) as the initial conditions to train the model from \(t=2.0\).

We first define the computational spatial domain, i.e., a rectangular region \([-5,10]\times[-5,5]\). As vortex shedding takes time to develop, we adapt a time sequence training scheme in which the time domain [2, 10]

Figure 16: Problem setup and training strategy of a transversely oscillating cylinder in steady flow. (a) The geometry of the computational domain and boundary conditions. (b) Visual representation of scattered training points around the oscillating cylinder boundary and in the fluid domain in space and time. (c) Snapshots of spatial training points sampled within the temporal neighborhood \(\epsilon=0.07\) at \(t\)=0.86. The training points are sampled with higher resolution in the vicinity of the cylinder and in the wake region where vortex shedding is prominent. (d) Aerial view of the training points sampled at the moving cylinder surface in the time domain [2, 10] from the \(x\)-direction. (e) The training losses versus the number of optimization epochs in two randomly selected time domains, [2, 3] and [5, 6], respectively.
is discretized into \(n\) domains as: [\(T_{0}=2\), \(T_{1}\)], [\(T_{1}\), \(T_{2}\)], \(\cdots\), [\(T_{n-1}\), \(T_{n}=10\)], to be solved individually. We set the range of the time sub-domain to \(\Delta T=2\). Next, we excavate the tunnel formed by the oscillation of the cylinder in the whole spatial-temporal domain. A visual representation of the tunnel corresponding to about two oscillation periods (one period is 3.465) is shown in Fig. 16(b). For each sub-domain [\(T_{n-1}\), \(T_{n}\)], a total of 452,235 FVM data points, corresponding to 21 time snapshots, scatter in space and time, with a time interval of \(\Delta t\) = 0.1 between two consecutive snapshots. Among these FVM data points, a total of 150,000 data points are randomly selected to reconstruct the entire field. In each sub-domain [\(T_{n-1}\), \(T_{n}=T\)], we randomly sample 50,000 domain points, 10,000 cylinder boundary points, and 10,000 rectangular boundary points. Moreover, we additionally sample 30,000 finer points in the vicinity of the cylinder and 30,000 points in the wake region where vortex shedding is prominent. Fig. 16(c) displays a snapshot of the distribution of sampled points from the above-mentioned processes within the temporal neighborhood \(\epsilon\) = 0.07 at \(t\)=4.32. Furthermore, we provide an aerial perspective of the training points that are sampled on the moving cylinder surface in \(t\)=[2, 10], from the \(x\) direction, as shown in Fig. 16(d).

Fig. 16(e) illustrates the convergence of the total loss and its individual components during training, including the PDEs loss, moving boundary loss, stationary boundary loss, and initial boundary loss. Fig. 17 shows the predicted velocity and pressure at time instances \(t\)=1.0, 4.0, and 7.0. Additionally, Fig. 18 illustrates the vorticity predictions of the trained model for the time period of t=[0, 7] and compared with

Figure 17: PINNs-predicted and FVM results for velocity \(u\), \(v\) and pressure \(p\) contours for several time frames. (a) \(t\)=1.0; (b) \(t\)=4.0; (c) \(t\)=7.0.
the FVM results. The comparison of vorticity indicates that the model successfully reproduces the vortex shedding phenomenon in alignment with the FVM calculation. Despite the slightly larger dissipation in the wake region, the predicted solution accurately and comprehensively captures the process of vortex shedding from the cylinder and its development in the wake region. Fig. 19 provides the point-wise relative \(L_{2}\) errors of inferred velocity and pressure solutions at each time. The relative \(L_{2}\) errors illustrate the rapid accumulation of errors in time series model predictions over time. Offering initial conditions with high precision may mitigate errors to a certain extent in the early stages of the simulation, but these errors will surge rapidly in the subsequent stages. The velocity distributions in the \(x\)-directional cross-section are further shown in Fig. 20 for \(u\) and \(v\) in three time frames (\(t\) = 3.0, 5.0, and 7.0) of the oscillating cylinder, and compared with FVM results. The comparison shows that the proposed model accurately captures finer flow distributions. Overall, our strategy for PINNs is capable of reproducing the phenomenon of vortex shedding with an acceptable margin of error.

Figure 19: Relative \(L_{2}\) errors of PINNs for the transversely oscillating cylinder in steady flow: the inferred velocity \(u\), \(v\) and pressure \(p\).

Figure 18: Results of PINNs and FVM for vorticity contours in several time frames. (a) \(t\)=0.1; (b) \(t\)=1.0; (c) \(t\)=2.0; (d) \(t\)=3.0; (e) \(t\)=4.0; (f) \(t\)=5.0; (g) \(t\)=6.0; (h) \(t\)=7.0. A significant vortex shedding phenomenon is observed.
### Flow around a flapping wing

In this case, we still consider a flapping wing that undergoes both translational and rotational motions. We solve this problem to demonstrate the flexibility of our proposed framework in the moving boundary scenario, that is, reconstruct the entire flow field (velocity \(v\) and pressure \(p\)) based on the partial data (velocity \(u\)) of the spatial-temporal domain of the flow field.

Here, we employ the velocity \(u\) data at the grid points of the FVM as Dirichlet constraint to train the neural network. Meanwhile, we take the velocity and pressure data on \(8,967\) points obtained from the FVM at \(t^{\prime}=20.0\) as the initial conditions to train the model from \(t=20.0\). A single neural network is used to reconstruct the flow field in the time domain [20.0, 20.4]. A total of 367,738 FVM data points, corresponding to 41 time snapshots, scatter in space and time, with a time interval of \(\Delta t\) = 0.01 between two consecutive snapshots. Among these FVM data points, a total of 120,000 data points are randomly selected to reconstruct the entire field. In addition, 50,000 collocation points are randomly sampled in the domain to complement the resolution of the computational domain. 25,000 and 10,000 points are sampled on the wing surface and rectangular boundaries, respectively. 40,000 finer sampled points are sampled within a circular region with a radius of \(1.5C\) centered on the wing. 10,000 finer sampled points are collected within a rectangular region \([-1.3,1.3]\times[2.5,3.5]\). Fig. 21(a) depicts a snapshot at \(t\)=20.0 demonstrating

Figure 20: The PINNs-predicted velocity contours and the velocity \(u\), \(v\) distribution in the \(y\)-direction for several frames of the cylinder, and are compared with the FVM results. Solid and dotted lines denote PINNs and FVM results in four profiles(indicated with colors: red: \(y=-0.6D\); green: \(y=0D\); blue: \(y=0.6D\); orange: \(y=1.2D\)).
the distribution of data points and sampled collocation points.

And the convergence of various training losses is given (see Fig. 21(b)). Fig. 22 shows the snapshots of inferred velocity \(v\) and pressure \(p\) contours in two time frames \(t=20.2\), \(20.4\). Fig. 23 provides the point-wise relative \(L_{2}\) errors at each time from \(t=20.0\) to \(20.4\). Obviously, providing partial data for the entire flow field can significantly reduce the relative error in velocity and especially pressure compared to providing data only of the initial conditions. This demonstrates the flexibility of our proposed strategy in solving moving boundary scenarios, which can rely on different types or amount of data for the prediction and reconstruction of the entire flow field.

Figure 21: Training strategy for the reconstruction of the flow field around a flapping wing. (a) Snapshots of data points randomly selected and training points sampled within the temporal neighborhood \(\epsilon=0.02\) at \(t\)=20. (b) The training losses versus the number of optimization epochs.

Figure 22: Flow around a flapping wing: inferred velocity \(v\) and pressure \(p\) contours at two time frames: (a) \(t\)=20.2; (b) \(t\)=20.4.
### A cylinder settling under gravity

In the last case, we conduct simulations pertaining to the free settling of a cylinder under gravity. This instance can be classified as a two-way fluid-structure interaction problem. In this case, the settling trajectory of the cylinder is extracted from the FVM results, making this problem degenerate into a flow problem with prescribed motion to be solved using our proposed strategy.

The objective is to employ partial flow field information for the reconstruction of the complete flow

Figure 23: Relative \(L_{2}\) errors of PINNs for flow around a flapping wing: the inferred velocity \(u\), \(v\) and pressure \(p\).

Figure 24: Problem setup and training strategy of a cylinder settling under gravity. (a) The geometry of the computational domain and boundary conditions. (b) Visual representation of scattered training points around the moving boundary and grid points in space and time. (c) Snapshots of the distribution of data points at \(t=4.4\). The spatial resolution of the points is gradually enhanced in the proximity of the cylinder region. (d) Aerial view of all the points sampled at the moving cylinder surface from the \(x\) (left) and \(y\) (right) directions. (e) The training losses versus the number of optimization epochs.
field, with the pressure \(p\) inferred from the known velocities \(u\) and \(v\) obtained via the FVM. In the present study, a computational domain is defined by intercepting a portion of the whole region in which the cylinder is settling. The geometric description and boundary conditions of the domain are provided in Fig. 24(a). The simulation parameters are set as follows: the gravity \(g=1\), the cylinder diameter \(D=0.2\), the fluid density \(\rho=1\), and the kinematic viscosity \(\nu=0.001\). In this computational region, the maximum velocity of the settling cylinder is approximately \(U_{max}=1.79\), resulting in a Reynolds number of \(Re=U_{max}D/\nu=358\) for this problem. The Neumann velocity conditions and the zero pressure condition are imposed on the left and right sides of the far-field boundaries, as depicted in Fig. 24(a). We impose Dirichlet velocity conditions on the boundary of the cylinder, where the velocities \(u\) and \(v\) are determined based on the positions \(x\) and \(y\) of the cylinders. Notably, we employ the velocities \(u\) and \(v\) at the grid points of the FVM as labeled data to train the neural network.

The cylinder settles from rest at \(y=0\). The selected spatial domain for simulation corresponds to a rectangular region bounded by \([-1,1]\times[-6,-3]\). And the time domain is [4.0, 4.4]. The location and velocity of the settling cylinder are obtained by fitting the Fourier series to the cylinder trajectory extracted from the FVM results. We excavate the tunnel formed by the settling cylinder in the whole spatial-temporal domain. A visual representation of this tunnel corresponding to the settling trajectory is shown in Fig. 24(b). A total of 229,433 data points, corresponding to 21 time snapshots, scattering in space and time, have been utilized to perform inference of the pressure field. The time interval between two consecutive snapshots is \(\Delta t=0.02\). Fig. 24(c) depicts a snapshot at \(t\)=4.4 demonstrating the distribution of these data points. In addition to the 1575 sparse data points in the 21 time snapshots that are at the cylinder boundary, we randomly sample an additional 10,000 data points on cylinder surface as a supplement. Fig. 24(d) provides an aerial perspective of all the points located at the moving cylinder surface, viewed from both the \(x\) and

Figure 25: A cylinder settling under gravity: PINNs-inferred velocity \(u\), \(v\) fields and pressure \(p\) field contours at two time frames are in the first columns for (a) \(t\)=4.1; (b) \(t\)=4.4. The corresponding results of FVM are in the second columns. The point-wise absolute errors for each frame are in the third columns.
\(y\) directions. The convergence of the total loss and its individual components during training, including PDEs loss, moving boundary loss, stationary boundary loss, and velocity \((u,v)\) data loss, are presented in Fig. 24(e).

Fig. 25 shows snapshots of inferred velocity \(u\), \(v\) and pressure \(p\) contours, as well as a visual comparison with the FVM results at two time frames. The proposed framework is capable of accurately reconstructing the pressure. Fig. 26 provides the point-wise relative \(L_{2}\) errors from \(t\)= 4.0 to 4.4. The network demonstrates the capability to accurately reconstruct the complete flow field.

## 6 Discussions and conclusions

We presented a novel extension to incorporate moving boundary conditions in fluid mechanics into the loss functions of physics-informed neural networks (PINNs). As we interpreted the time-dependent moving boundaries in the spatial domain as stationary boundaries in a higher dimension of the spatial-temporal domain, we can effortlessly distribute fine training points at and around the interfaces for accuracy. Therefore, we are able to solve one class of unsteady flow problems with moving bodies, where the time-dependent interfaces are available via other means. We have validated the extended PINNs on a variety of classical flow problems involving moving bodies. For instance, we solved the velocity, vorticity, and pressure fields for a complete cycle of an in-line oscillating cylinder in fluid at rest using only data for the initial condition. Additionally, we predicted the flow field surrounding two, three, or four cylinders translating along a circle and demonstrated the effectiveness of this approach for multiple moving bodies. Furthermore, for flows around a flapping wing, we illustrated the flexibility of the proposed framework first by relying data from the initial conditions to predict the later flow field in Section 4 and subsequently using data on the velocity \(u\) to infer the velocity \(v\) and pressure for a reconstruction of the entire flow field in Section 5. In the case of a transversely oscillating cylinder in flow, we also reproduced the phenomenon of rapid vortex shedding due to oscillation with pressure data. Lastly, for a cylinder settling under gravity, we relied on both the trajectory of the cylinder and the velocity data to infer the pressure field. Comparisons with reference solutions for these flow problems demonstrated the effectiveness and accuracy of the extended PINNs.

The extension proposed here may be insightful for flow problems solved by other machine learning frameworks and should be also applicable to physics problems with time-dependent moving boundaries beyond fluid mechanics. As PINNs exploits a fully connected neural network to approximate the solution of PDEs, currently we cannot achieve the same order of accuracy as the CFD methods within the same wall-time window. For example, the in-line oscillating cylinder case took about 3 hours to solve for a full oscillating period on a Nvidia RTX 4090 graphic card. Improving the accuracy further would require an even larger number of training points and a more significant amount of epochs for optimization, which both

Figure 26: Relative \(L_{2}\) errors of PINNs for a cylinder settling under gravity: the inferred velocity \(u\), \(v\) and pressure \(p\).
result in a higher computational cost. Due to this exact reason, it does not seem necessary to further extend PINNs to handle the two-way coupling problems, which may be more feasible by CFD methods currently. For flow problems with access to the trajectories of moving interfaces, yet with limited measurements on the flow fields, the proposed method can be applied steadily. A potential revolution to the existing neural network architecture and training strategy may also fundamentally reform the status and allow PINNs for addressing more computational challenges in fluid mechanics.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgments

Support from the grant of Innovative Research Foundation of Ship General Performance under contract number 31422121 is gratefully acknowledged. X. Bian also received the starting grant from 100 talents program of Zhejiang University.

## Appendix A Technical details of FVM simulations

The simulations by FVM in this study have been executed within OpenFOAM. Morphing mesh and overset mesh techniques are constructed for handling the moving boundaries. We employed the pimple-Foam solver coupled with the morphing mesh technique to solve flow problems in Section 4.1, 5.1, and used the overPimpleDyMFoam solver coupled with the overset mesh technique to solve flow problems in Section 4.2, 4.3, 5.2, 5.3. The time step is automatically adjusted within the limit of the Courant number being 0.4. Other technical details are given in the captions of the figures as follows.
Figure 28: Multiple cylinders translating along a circle: a uniform background mesh (\(80\times 80\)) with finer resolution overset mesh (each contains a total of 2,250 cells and 4,650 nodes) around each cylinder is adopted at time t=2.5. And the left, middle, and right panels depict the configurations of 2, 3, and 4 cylinders, respectively.

Figure 29: Flow around a flapping wing: a uniform background mesh (\(70\times 80\)) with a finer resolution overset mesh (contains a total of 13,564 cells and 17,720 nodes) around wing is adopted at time t=20. Adaptive mesh refinement technology is adopted to improve the resolution of the mesh near the wing.
Figure 31: A cylinder settling under gravity: a uniform background mesh (\(80\times 120\)) with a finer resolution overset mesh (contains a total of 2,250 cells and 4,650 nodes) around cylinder is adopted.

Figure 30: A transversely oscillating cylinder in steady flow: a morphing mesh with a finer resolution zone of (Diameter=3D) around cylinder is adopted at time t=7. It contains a total of 21,200 cells and 43,070 nodes.