Fundamentals of Physics-Informed Neural Networks Applied to Solve the Reynolds Boundary Value Problem

Abstract

This paper presents a complete derivation and design of a physics-informed neural network (PINN) applicable to solve initial and boundary value problems described by linear ordinary differential equations. The objective with this technical note is not to develop a numerical solution procedure which is more accurate and efficient than standard finite element- or finite difference-based methods, but to give a fully explicit mathematical description of a PINN and to present an application example in the context of hydrodynamic lubrication. It is, however, worth noticing that the PINN developed herein, contrary to FEM and FDM, is a meshless method and that training does not require big data which is typical in machine learning.

1. Introduction

There are various categories of artificial neural networks (ANN) and a physics-informed neural network (PINN), see [1] for a recent review on the matter, is a neural network trained to solve both supervised and unsupervised learning tasks while satisfying some given laws of physics, which may be described in terms of nonlinear partial differential equations (PDE). For example, the balance of momentum and conservation laws in solid- and fluid mechanics and various types of initial value problems (IVP) and boundary value problems (BVP), see e.g., [2,3]. The application of a PINN (of this type) to solve differential equations, renders meshless numerical solution procedures [4], and an important feature from a machine learning perspective, is that it is not a data-driven approach requiring a large set of training data to learn the solution.

In fluid mechanics, under certain assumptions, i.e., that the fluid is incompressible, iso-viscous, the balance of linear momentum and the continuity equation, for flows in narrow interfaces reduces to the classical Reynolds equation [5]. For more recent work establishing lower-dimensional models in a similar manner, see e.g., [6,7,8]. The present work describes how a PINN can be adapted and trained to solve both initial and boundary value problems, described by ordinary differential equations, numerically. The theoretical description starts by presenting the neural network’s architecture and it is first applied to solve an initial value problem, which is described by a first order ODE, which can be solved analytically so that the validity of the solution can be thoroughly assessed. Thereafter, it is used to obtain a PINN for the classical one-dimensional Reynolds equation, which is a boundary value problem governing, e.g., the flow of lubricant between the runner and the stator in a 1D slider bearing. The novelty and originality of the present work lays in the explicit mathematical description of the cost function, which constitutes the physics-informed feature of the ANN, and the associated gradient with respect to the networks weights and bias. Important features of this particular numerical solution procedure, that is publicly available here: [9], are that it is not data driven, i.e., no training data need to be provided and that it is a meshless method [4].

2. PINN Architecture

Knowing the characteristics of the solution to the differential equation under consideration is very helpful when designing the PINN architecture, including structure, number of hidden layers, activation function, etc. For this reason, the PINN developed here has one input node x (the independent variable representing the spatial coordinate), one hidden layer consisting of N nodes and one output node y (the dependent variable representing pressure). Figure 1 depicts a graphical illustration of the present architecture, which when trained solves both the IVP example and the Reynolds BVP considered here.

The Sigmoid function, i.e.,

$$\varphi \left(\xi \right)=\frac{1}{1+e^{-\xi }},$$

(1)

which is mapping $\mathbb{R}$ to $[0,1]$ and exhibits the property

$$\begin{array}{c}\hfill {\varphi }^{\prime }\left(\xi \right)=\varphi \left(\xi \right)\left(1-\varphi \left(\xi \right)\right).\end{array}$$

(2)

is employed as activation function for the hidden layer. This means that the neural network has $3N+1$ trainable parameters. That is, the weights $w_i^{\left(0\right)}$ and bias $b_i^{\left(0\right)}$ for the nodes in the hidden layer and the weights $w_i^{\left(1\right)}$, $i=1\dots N$, for each synapses connecting them with the output node, plus the bias $b^{\left(1\right)}$ applied there.

Based on this particular architecture, the output $z_i$ of each node in the first hidden layer is,

$$z_i\left(x\right)=\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right).$$

(3)

The output value is then given by applying the Sigmoid activation function scaled by the weight from the node in the second layer and yields

$$\begin{array}{cc}\hfill y\left(x\right)& =b^{\left(1\right)}+\sum _{i=1}^N{w}_i^{\left(1\right)}{z}_i\left(x\right)=b^{\left(1\right)}+\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right).\hfill \end{array}$$

(4)

Let us now construct the cost function which the network will be trained to minimise. While the cost function appearing in a typical machine learning procedure is just the quadratic difference between the predicted and the target values, it will here be defined by means of the operators $\mathcal{L}$ and $\mathcal{B}$. The cost function applied here reads

$$\begin{array}{c}\hfill l=〈{\left(\mathcal{L}y-f\right)}^2〉+{\left(\left(\mathcal{B}y-\mathbf{b}\right)·{\mathbf{e}}_{\mathbf{1}}\right)}^2+{\left(\left(\mathcal{B}y-\mathbf{b}\right)·{\mathbf{e}}_{\mathbf{2}}\right)}^2,\end{array}$$

(5)

where $〈f〉$ defines the average value of f, and this is exactly the feature that makes an ANN “physics informed”, i.e., a PINN.

Since $\mathcal{L}y$ is a differential operator the cost function contains derivatives of the network output (4). In order to obtain an expression of the cost function, in terms of the input x, the weights w and bias b, the network output (4), must be differentiated twice with respect to (w.r.t.) x. This can be accomplished by some kind of automatic differentiation (AD) (also referred to as algorithmic differentiation, computer differentiation, auto-differentiation or simply autodiff), which is a computerised methodology based on the chain rule, which can be applied to efficiently and accurately evaluate derivatives of numeric functions, see e.g., [10,11]. The present work instead applies symbolic differentiation to clearly explain all the essential details of the PINN. Indeed, differentiating one yield

$$\begin{array}{cc}\hfill & y^{\prime }\left(x\right)=\frac{\partial }{\partial x}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}{z}_i\left(x\right)\right)+b^{\left(1\right)}\right)=\frac{\partial }{\partial x}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=\hfill \\ \hfill & =\sum _{i=1}^N{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=\sum _{i=1}^N{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\left(1-\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right),\hfill \end{array}$$

(6)

and, because of (2), a consecutive differentiation then yields

$$\begin{array}{cc}\hfill & y^{\prime \prime }\left(x\right)=\frac{\partial }{\partial x}\left(\sum _{i=1}^N{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)=\sum _{i=1}^N{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=\hfill \\ \hfill & =\sum _{i=1}^N{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\left(1-2\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)=\hfill \\ \hfill & =\sum _{i=1}^N{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\left(1-\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)\left(1-2\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right).\hfill \end{array}$$

(7)

Moreover, finding the set of weights and bias minimising the cost function requires its partial derivatives w.r.t. to each weight and bias defining the PINN. In the subsections below, we will present how to achieve this, by first considering a first order differential equation with an analytical solution, and, thereafter, we will consider the classical Reynolds equation which is a second order (linear) ODE that describes laminar flow of incompressible and iso-viscous fluids in narrow interfaces.

3. A First Order ODE Example

Let us consider the first order ODE, describing the initial value problem (IVP) given by

$$\begin{array}{cc}& \mathcal{L}y-f=y^{\prime }+2xy=0,x>0\hfill \end{array}$$

(8a)

$$\begin{array}{cc}& \mathcal{B}y-\mathbf{b}=y\left(0\right)-1=0,\hfill \end{array}$$

(8b)

with the exact solution $y=e^{-x^2}$. By means of (6), a cost function suitable for solving (8) may be generated by

$$\begin{array}{cc}\hfill l=& 〈\left[\sum _{i=1}^N{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\left(1-\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)\right.+\hfill \\ \hfill & {\left.2x\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)\right]}^2〉+{\left[y\left(0\right)-1\right]}^2\hfill \end{array}$$

(9)

The solution of (8) can be obtained by implementing a training routine which iteratively finds the set of weights w and bias b that minimises (9) (and similarly for (19) minimising (17)). The most well-known of these is the Gradient Decent method attributed to Cauchy, who first suggested it in 1847 [12]. For an overview, see, e.g., [13].

As mentioned in the previous section, the derivatives of (4) w.r.t. to the weights w and bias b are required to find them, and automatic differentiation is, normally, employed to perform the differentiation. However, here we carry out symbolic differentiation to demonstrate exactly the explicit expressions that constitute the gradient of the cost function. Indeed, by taking the partial derivatives we obtain

$$\begin{array}{cc}\hfill & \frac{\partial y}{\partial w_i^{\left(0\right)}}=\frac{\partial }{\partial w_i^{\left(0\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=w_i^{\left(1\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)x,\hfill \end{array}$$

(10a)

$$\begin{array}{cc}\hfill & \frac{\partial y}{\partial w_i^{\left(1\right)}}=\frac{\partial }{\partial w_i^{\left(1\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(10b)

$$\begin{array}{cc}\hfill & \frac{\partial y}{\partial b_i^{\left(0\right)}}=\frac{\partial }{\partial b_i^{\left(0\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=w_i^{\left(1\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(10c)

$$\begin{array}{cc}\hfill & \frac{\partial y}{\partial b^{\left(1\right)}}=1.\hfill \end{array}$$

(10d)

Moreover, the derivatives of the cost function (5) w.r.t. to the weights and bias are also required. For the derivative w.r.t. $w_i^{\left(0\right)}$ for the first order ODE (8), this means that

$$\begin{array}{c}\hfill 〈2\left(y^{\prime }+2xy\right)\left(\frac{\partial y^{\prime }}{\partial w_i^{\left(0\right)}}+2x\frac{\partial y}{\partial w_i^{\left(0\right)}}\right)〉+2\left(y\left(0\right)-1\right)\frac{\partial y\left(0\right)}{\partial w_i^{\left(0\right)}}.\end{array}$$

(11)

To complete the analysis, we also need expressions for the derivatives of $y^{\prime }$ w.r.t. $w_i^{\left(0\right)}$, $w_i^{\left(1\right)}$, $b_i^{\left(0\right)}$ and $b^{\left(1\right)}$. By the chain rule, the following expressions can be obtained, viz.

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime }}{\partial w_i^{\left(0\right)}}=\frac{\partial }{\partial w_i^{\left(0\right)}}\sum _{i=1}^N{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=\hfill \end{array}$$

$$\begin{array}{cc}\hfill & =w_i^{\left(1\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)+x{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(12a)

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime }}{\partial w_i^{\left(1\right)}}=\frac{\partial }{\partial w_i^{\left(1\right)}}{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=w_i^{\left(0\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(12b)

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime }}{\partial b_i^{\left(0\right)}}=\frac{\partial }{\partial b_i^{\left(0\right)}}\sum _{i=1}^N{w}_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=w_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(12c)

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime }}{\partial b^{\left(1\right)}}=0.\hfill \end{array}$$

(12d)

What remains now is to obtain expressions for $y\left(0\right)$ and the partial derivatives of $y\left(0\right)$, w.r.t. to the weights and bias. Let us start with $y\left(0\right)$. With $y\left(x\right)$ given by (4) we directly have

$$\begin{array}{cc}\hfill y\left(0\right)& =\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)},\hfill \end{array}$$

(13)

which, in turn, means that

$$\begin{array}{cc}\hfill & \frac{\partial y\left(0\right)}{\partial w_i^{\left(0\right)}}=0,\hfill \end{array}$$

(14a)

$$\begin{array}{cc}\hfill & \frac{\partial y\left(0\right)}{\partial w_i^{\left(1\right)}}=\frac{\partial }{\partial w_i^{\left(1\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=\varphi \left(b_i^{\left(0\right)}\right),\hfill \end{array}$$

(14b)

$$\begin{array}{cc}\hfill & \frac{\partial y\left(0\right)}{\partial b_i^{\left(0\right)}}=\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=w_i^{\left(1\right)}{\varphi }^{\prime }\left(b_i^{\left(0\right)}\right)\hfill \end{array}$$

(14c)

$$\begin{array}{cc}\hfill & \frac{\partial y\left(0\right)}{\partial b^{\left(1\right)}}=1.\hfill \end{array}$$

(14d)

The PINN (following the architecture presented above) was implemented as computer program in MATLAB. The program was employed to obtain a numerical solution to the IVP in (8), using the parameters in

Table 1. Parameters used to defined the PINN to for the IVP in (8).

Parameter	Description	Value
$N_i$	# of grid points for the solution domain $[0,2]$	41
$N_e$	# of training batches (# or corrections during 1 Epoch)	1000
$T_b$	# of Epochs (1 Epoch contains $T_b$ training batches)	100
$L_r$	Learning rate coefficient (relaxation for the update)	0.01
N	# of nodes/neurons in the hidden layer	10

.

The weights $w_i^{\left(0\right)}$ and bias $b_i^{\left(0\right)}$ were initialised using randomly generated and uniformly distributed numbers in the interval $[-2,2]$, while the weights $w_i^{\left(1\right)}$ were initially set to zero and the bias $b^{\left(1\right)}$ to one, to ensure fulfilment of the initial condition ($y\left(0\right)=1$).

Table 2. Parameters used to defined the PINN for the IVP (8).

Node	${\mathit{w}}^{\left(0\right)}$	${\mathit{b}}^{\left(0\right)}$	${\mathit{w}}^{\left(1\right)}$	${\mathit{b}}^{\left(1\right)}$
1	1.8500	−0.5946	−3.5805	0.3055
2	1.8588	1.5974	0.9712
3	0.3025	1.9241	0.8921
4	1.4546	0.3742	−0.9955
5	0.5065	1.2535	−0.1430
6	−1.0898	−1.0199	−1.1067
7	−0.8302	0.3519	−1.1668
8	0.3789	1.6502	0.1754
9	2.5012	0.7657	1.2955
10	2.2743	1.4172	1.2787

lists the weights an bias corresponding to the solution presented in Figure 2. We note that, with the weights and bias given by Table 2, the trained network’s prediction exhibits the overall error

$$\begin{array}{cc}\hfill & \frac{1}{N_i}\sqrt{\sum _{k=1}^{N_i}{\left(e^{-x_k^2}-y\left(x_k\right)\right)}^2}=5.8\times {10}^{-4},\hfill \end{array}$$

(15)

and $1-y\left(0\right)=2.2\times {10}^{-4}$, when comparing against the initial condition.

4. A PINN for the Classical Reynolds Equation

The Reynolds equation for a one-dimensional flow situation, where the lubricant is assumed to be incompressible and iso-viscous, is a second order Boundary Value Problem (BVP), which in dimensionless form can be formulated as

$$\begin{array}{cc}& \frac{d}{dx}\left(c\left(x\right)\frac{dy}{dx}\right)=f\left(x\right),0<x<1,\hfill \end{array}$$

(16a)

$$\begin{array}{cc}& y\left(0\right)=0,y\left(1\right)=0,\hfill \end{array}$$

(16b)

where $c\left(x\right)=H^3$, $f\left(x\right)=dH/dX$ and H is the dimensionless film thickness, if it is assumed that the pressure y at the boundaries is zero. For the subsequent analysis it is, however, more suitable work with a condensed form which can be obtained by defining the operators $\mathcal{L}$ and $\mathcal{B}$ as

$$\begin{array}{cc}& \mathcal{L}y=c\left(x\right)y^{\prime \prime }+c^{\prime }\left(x\right)y^{\prime },\hfill \end{array}$$

(17a)

$$\begin{array}{cc}& \mathcal{B}y=\left[\begin{array}{c}y\left(0\right)\\ y\left(1\right)\end{array}\right].\hfill \end{array}$$

(17b)

The Reynolds BVP given by (16) can then be presented as

$$\begin{array}{cc}& \mathcal{L}y-f=0,0<x<1,\hfill \end{array}$$

(18a)

$$\begin{array}{cc}& \mathcal{B}y-\mathbf{b}=\mathbf{0},\hfill \end{array}$$

(18b)

where $\mathbf{b}=\mathbf{0}$.

For the Reynolds BVP, the cost function (5) becomes

$$\begin{array}{c}\hfill l=〈{\left(c\left(x\right)y^{\prime \prime }+c^{\prime }\left(x\right)y^{\prime }-f\right)}^2〉+y^2\left(0\right)+y^2\left(1\right),\end{array}$$

(19)

and from the analysis presented for the IVP in Section 3 above, we have all the “ingredients” except for the partial derivatives of $y^{\prime \prime }$ and $y\left(1\right)$ w.r.t. to the weights and bias. For $y^{\prime \prime }$, based on (7) and (12), we obtain

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime \prime }}{\partial w_i^{\left(0\right)}}=\frac{\partial }{\partial w_i^{\left(0\right)}}\sum _{i=1}^N{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=\hfill \end{array}$$

$$\begin{array}{cc}\hfill & =2w_i^{\left(1\right)}{w}_i^{\left(0\right)}{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)+x{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(20a)

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime \prime }}{\partial w_i^{\left(1\right)}}=\frac{\partial }{\partial w_i^{\left(1\right)}}\sum _{i=1}^N{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)={\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(20b)

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime \prime }}{\partial b_i^{\left(0\right)}}=\frac{\partial }{\partial b_i^{\left(0\right)}}\sum _{i=1}^N{w}_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right)=w_i^{\left(1\right)}{\left(w_i^{\left(0\right)}\right)}^2{\varphi }^{\prime \prime \prime }\left(w_i^{\left(0\right)}x+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(20c)

$$\begin{array}{cc}\hfill & \frac{\partial y^{\prime \prime }}{\partial b^{\left(1\right)}}=0,\hfill \end{array}$$

(20d)

where the third derivative of the Sigmoid function (1) is required. It yields

$$\begin{array}{cc}\hfill & \frac{d}{d\xi }\left({\varphi }^{\prime \prime }\left(\xi \right)\right)=\frac{d}{d\xi }\left({\varphi }^{\prime }\left(\xi \right)\left(1-2\varphi \left(\xi \right)\right)\right)={\varphi }^{\prime \prime }\left(\xi \right)\left(1-2\varphi \left(\xi \right)\right)-2{\left({\varphi }^{\prime }\left(\xi \right)\right)}^2=\hfill \\ \hfill & =\varphi \left(\xi \right)\left(1-\varphi \left(\xi \right)\right){\left(1-2\varphi \left(\xi \right)\right)}^2-2{\left(\varphi \left(\xi \right)\left(1-\varphi \left(\xi \right)\right)\right)}^2=\hfill \\ \hfill & =\varphi \left(\xi \right){\left(1-\varphi \left(\xi \right)\right)}^2\left(1-3\varphi \left(\xi \right)\right).\hfill \end{array}$$

(21)

For $y\left(1\right)$ we obtain

$$\begin{array}{cc}\hfill & \frac{\partial y\left(1\right)}{\partial w_i^{\left(0\right)}}=\frac{\partial }{\partial w_i^{\left(0\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=w_i^{\left(1\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(22a)

$$\begin{array}{cc}\hfill & \frac{\partial y\left(1\right)}{\partial w_i^{\left(1\right)}}=\frac{\partial }{\partial w_i^{\left(1\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=\varphi \left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(22b)

$$\begin{array}{cc}\hfill & \frac{\partial y\left(1\right)}{\partial b_i^{\left(0\right)}}=\frac{\partial }{\partial b_i^{\left(0\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=w_i^{\left(1\right)}{\varphi }^{\prime }\left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right),\hfill \end{array}$$

(22c)

$$\begin{array}{cc}\hfill & \frac{\partial y\left(1\right)}{\partial b^{\left(1\right)}}=\frac{\partial }{\partial b^{\left(1\right)}}\left(\left(\sum _{i=1}^N{w}_i^{\left(1\right)}\varphi \left(w_i^{\left(0\right)}+b_i^{\left(0\right)}\right)\right)+b^{\left(1\right)}\right)=1,\hfill \end{array}$$

(22d)

and we now have all the “ingredients” required to fully specify (19). To test the performance of the PINN, a Reynolds BVP was specified for a linear slider with dimensionless film thickness defined by

$$H\left(x\right)=1+K-Kx.$$

(23)

This means that $c\left(x\right)={(1+K-Kx)}^3$ and $f\left(x\right)=dH/dx=-K$ and that the exact solution is

$$y_{exact}\left(x\right)=\left[\frac{1}{K}\left(\frac{1}{1+K-Kx}-\frac{1+K}{2+K}\frac{1}{{\left(1+K-Kx\right)}^2}-\frac{1}{2+K}\right)\right],$$

(24)

see, e.g., [14].

The PINN (following the architecture suggested herein) was implemented in MATLAB and a numerical solution to (16) was obtained using the parameters in

Table 3. Parameters used to defined the ANN to for the Reynolds equation.

Parameter	Description	Value
$N_i$	# of grid points for the solution domain $[0,1]$	21
K	Slope parameter for the Reynolds equation	1
$N_e$	# of training batches (# or corrections during 1 epoch)	2000
$T_b$	# of Epochs (1 epoch contains $T_b$ training batches)	600
$L_r$	Learning rate coefficient (relaxation for the update)	0.005
N	# of nodes/neurons in the hidden layer	10

. As for the IVP, addressed in the previous section, the weights $w_i^{\left(0\right)}$ and bias $b_i^{\left(0\right)}$ were, again, initialised using randomly generated numbers, uniformly distributed in $[-2,2]$, while the weights $w_i^{\left(1\right)}$ and the bias $b^{\left(1\right)}$ was initially set to zero, to ensure fulfilment of the boundary conditions.

Figure 3 depicts solution predicted by the PINN (red line with circle markers) and the exact solution obtained by integration (blue continuous line).

Table 4. Parameters used to defined the ANN.

Node	${\mathit{w}}^{\left(0\right)}$	${\mathit{b}}^{\left(0\right)}$	${\mathit{w}}^{\left(1\right)}$	${\mathit{b}}^{\left(1\right)}$
1	0.0557	1.9808	−0.2186	−0.0641
2	−6.3047	6.1664	0.1220
3	−9.3674	11.4571	0.3843
4	−4.5473	3.3266	0.0305
5	−2.4464	−1.9884	0.1188
6	−0.1365	−0.1674	0.4155
7	0.8581	0.5253	0.5089
8	1.0901	2.0858	0.3348
9	0.2085	0.2523	−0.2024
10	−3.2168	5.9722	−0.9899

lists the weights an bias corresponding to the solution presented in Figure 3.

We note that, with these weights and bias, the trained network’s prediction of the solution to the Reynolds BVP exhibits the overall error

$$\begin{array}{cc}\hfill & \frac{1}{N_i}\sqrt{\sum _{k=1}^{N_i}{\left(y_{exact}\left(x_k\right)-y\left(x_k\right)\right)}^2}=6.2\times {10}^{-5},\hfill \end{array}$$

(25)

while $y\left(0\right)=4.1\times {10}^{-4}$ and $y\left(1\right)=-4.0\times {10}^{-4}$.

Remark 1.

The formulation of the PINN presented here is applicable as a numerical solution procedure for the Reynolds BVP (16) and it does not consider the effect of cavitation. Exactly how the effect of cavitation can be included is, however, out of the scope of this paper.

5. Concluding Remarks

A physics-informed neural network (PINN) applicable to solve initial and boundary value problems has been established. The PINN was applied to solve an initial value problem described by a first order ordinary differential equation and to solve the Reynolds boundary value problem, described by a second order ordinary differential equation. Both these problems were selected since they can be solved analytically, and the error analysis showed that the predictions returned by the PINN was in good agreement with the analytical solutions for the specifications given. The advantage of the present approach is, however, neither accuracy nor efficiency when solving these linear equations, but that it presents a meshless method and that it is not data driven. This concept may, of course, be generalised, and it is hypothesised that future research in this direction may lead to more accurate and efficient in solving related but nonlinear problems than currently available routines.