Physics-informed Neural Networks approach to solve the Blasius function
Deep learning techniques with neural networks have been used effectively in
computational fluid dynamics (CFD) to obtain solutions to nonlinear
differential equations. This paper presents a physics-informed neural network
(PINN) approach to solve the Blasius function. This method eliminates the
process of changing the non-linear differential equation to an initial value
problem. Also, it tackles the convergence issue arising in the conventional
series solution. It is seen that this method produces results that are at par
with the numerical and conventional methods. The solution is extended to the
negative axis to show that PINNs capture the singularity of the function at
$\eta=-5.69$
Greeshma Krishna, Malavika S Nair, Pramod P Nair, Anil Lal S
# Physics-informed Neural Networks approach to solve the Blasius function

###### Abstract

Deep learning techniques with neural networks have been used effectively in computational fluid dynamics (CFD) to obtain solutions to nonlinear differential equations. This paper presents a physics-informed neural network (PINN) approach to solve the Blasius function. This method eliminates the process of changing the non-linear differential equation to an initial value problem. Also, it tackles the convergence issue arising in the conventional series solution. It is seen that this method produces results that are at par with the numerical and conventional methods. The solution is extended to the negative axis to show that PINNs capture the singularity of the function at \(\eta=-5.69\).

 Blasius equation, physics-informed neural networks, automatic differentiation, computational fluid dynamics, boundary layer flow, singularity.

## I Introduction

Machine learning helps learn patterns from the data to make predictions. The three basic requirements for machine learning are data, theory, and hardware to overcome computational difficulties using better GPUs and new approaches. Computational science is an essential tool that we can use to incorporate physical invariances into learning. For example, the laws that govern the conservation of momentum, mass, and energy. To quote Dr.Tinsley Oden - "Computational Science can analyze past events and look into the future. It can explore the effects of thousands of scenarios instead of actual experiments and be used to study events beyond the reach of expanding the boundaries of experimental science". Deep learning can be quite useful in the real world, recognising various types of cancer like skin cancer [1], lung cancer [2] and many more by just passing the necessary images. Also has various applications in agriculture [7].

Neural networks(NNs) are the most used Machine Learning technique. It offers a set of robust tools for solving variously supervised and unsupervised problems in pattern recognition, data processing, and non-linear control, which can be regarded as complementary to the conventional approaches [3]. It can be considered a functional approximator as it maps a set of input variables to a set of output variables. A set of parameters called weights manages this mapping. These weights are updated to train the network that fetches the required output. A polynomial can be viewed as a function that transforms a single input variable into a single output variable. The coefficients in the polynomial are comparable to the weights in a neural network. Determining these coefficients helps in evaluating the solution.

### _Physics Informed Neural Networks_

Physics-informed neural networks (PINNs) play a huge role where fewer data and the system's conventional physics are known to get the solutions to a differential equation. NN can approximate the error to ground truth function, generalize to unseen data and train the model. The contribution of PINNs is to obtain a neural network that knows about the physics hidden behind the equation and thus efficiently solves differential equations. The earlier numerical methods for solving them are finite difference, finite elements, spectral elements, and finite volume. When working with PINNs, there is data efficiency from a machine learning perspective since it is regularized heavily with physics. The required derivatives can be calculated using automatic differentiation at the end of the network. It is possible to train a nonlinear neural network using this approach and its layers of hidden nodes [4]. It then finds a loss function corresponding to the differential equation and boundary conditions.

### _Boundary Layer Theory_

The boundary layer of a flowing fluid is a thin layer near a solid surface, and the flow near the solid surface is known as the boundary layer flows. Ludwig Prandtl is credited for developing the boundary layer theory. In 1904, he published a paper titled "On the Motion of a Fluid with Very Small Viscosity" [5]. He laid out the mathematical foundation for flows and condensed the two-dimensional Navier-Stokes equations (NSE) into the boundary layer equations. This
publication made understanding fluid motion physics possible, which is regarded as the beginning of contemporary fluid mechanics [6]. The primary disadvantage of using a CFD solver for turbomachinery optimization is the amount of time needed to complete the numerous computationally demanding simulations [7].

The introduction of the similarity variables and the transformation of the PDEs into nonlinear ODEs in one coordinate allowed the successful resolution of the problem. Still, specific methods were needed to handle the unbounded boundary conditions. The Blasius benchmark problem has been resolved [8] using the trial function method put out by Lagaris [9] or a hybrid approach [10]. The solution function is seen to have a singularity on the negative real axis at approximately -5.69. The Blasius equation is a single equation that solely models the viscous boundary layer and is one of CFD's fundamental models.

In this paper, we propose solving the Blasius equation using PINNs and comparing the solution we obtain with the best-known solutions available in the literature. We shall also extend our solution to the negative real axis to locate the singularity existing for the function. In the next section of the paper, we shall present a brief literature review on the Blasius equation and PINNs. In Section III, we present the methodology that has been adopted to solve the Blasius Equation. Section IV discusses the results and analysis of the proposed algorithm. The conclusions are discussed in the last section.

## II Literature review

### _Physics Informed Neural Networks_

Neural networks can be used to attain solutions to ODEs and PDEs by reducing them to an optimization problem instead of numerically solving the equations. The application of NNs in fluid mechanics began in the 1990s. Raissi and his team [11, 12] developed a technique called Physics informed neural network (PINN) where the loss function defined in the corresponding NN is extracted from the physics behind the PDE and related equations. In earlier years, Computational fluid dynamics (CFD) had been a great relief in numerically solving the compressible and incompressible Naive Stokes equations. Minimizing the loss function is challenging in PINNs, as it can be highly complex. Nevertheless, PINNs have been proven to be more accurate than conventional CFD techniques with lesser computations [13].

Although numerical discretization of the Navier-Stokes equations (NSE) has made significant progress in simulating flow problems over the past 50 years, the present algorithms cannot solve governed by high-parametrized NSE. Additionally, it is expensive to solve inverse flow problems due to their complexity, expensive formulations, and need for new algorithms. PINNs are expanded to fractional PINNs (fPINNs) that explore their convergence methodically to solve space-time fractional advection-diffusion equations (fractional ADEs) [14]. A hybrid method for building the residual in the loss function that utilizes both automated differentiation for the integer-order operators and numerical discretization for the fractional operators is a novel component of the fPINNs.

Wide varieties of PINNs have been found in the literature in recent times. can-PINNs [15] link derivative terms with nearby support points, which generally apply to Taylor series expansion-based numerical systems. Apart from demonstrating good dispersion and dissipation characteristics, they are highly trainable and require four to sixteen times fewer collocation points than original PINNs. Auxiliary PINNs (A-PINN) is a technique for solving forward and inverse problems of non-linear integrodifferential equations [16]. ViscoelasticNet is another PINN framework for stress discovery and model selection [17]. PINNs can also be used to solve Reynolds-averaged Navier-Stokes equations [18], full waveform seismic inversions in 2D acoustic media, and wave propagation as it seamlessly handles boundary conditions and physical constraints [19]. In addition to addressing ill-posed problems beyond the scope of conventional computing techniques, PINNs can also close the discrepancy between computational and experimental heat transfer [20].

### _Blasius Equation_

Blasius equation is a third order non-linear ordinary differential equation of the form \(f^{\prime\prime\prime}+\frac{1}{2}ff^{\prime\prime}=0\) with the boundary conditions \(f(0)=0\), \(f^{\prime}(0)=0\), \(f^{\prime}(\infty)=1\). It governs the boundary layer flow over a semi-infinite flat plate. Suppose that the \(u-\)velocity, the velocity parallel to the surface, is much greater than the \(v-\)velocity, perpendicular to the surface, and the changes in the perpendicular direction to the surface are much greater than changes parallel to the surface. The boundary layer equations consists of conservation of mass (II.1), conservation of x-momentum (II.2), and conservation of y-momentum. In a flat plate boundary layer, the pressure gradient term appearing in the x-momentum equation becomes zero (II.2). This leads to the hydrodynamic solution for the flat plate boundary layer in a laminar flow called the Blasius solution.

\[\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}=0\] (II.1)

\[u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}=\nu\frac{ \partial^{2}u}{\partial y^{2}}\] (II.2)

Blasius' analysis focus on the laminar boundary layer forming on a flat plate. The main aspect of the work is the transformation of the PDE for a flat plate boundary layer with zero pressure gradient into a single ordinary differential equation(ODE) by considering the velocity components that satisfy equation II.1.

\[u\equiv\frac{\partial\psi}{\partial y}\ \ \ \ v\equiv-\frac{\partial\psi}{ \partial x}\] (II.3)

The stream function \(\psi=U_{\infty}\sqrt{\frac{\nu\kappa}{U_{\infty}}}f(\eta)\), is directly proportional to the function \(f(\eta)\) called the Blasius function. Here, \(U_{\infty}\) is the free stream velocity, and \(\eta\) is a transformed coordinate called the similarity parameter. Here, the velocity components are proportional to the first derivative of \(f(\eta)\)
and the second and third derivatives of \(f\) are proportional to the first and second derivatives of velocity. Substituting these relations into the momentum equation (II.2), the final form of the Blasius boundary layer equation \(f^{\prime\prime\prime}+\frac{1}{2}ff^{\prime\prime}=0\), for a flat plate can be obtained. The first and second derivatives of \(f(\eta)\) are given by \(f^{\prime}=\frac{u}{U_{\infty}}\) (Non-dimensional velocity profile) and \(f^{\prime\prime}=\frac{1}{U_{\infty}}\sqrt{\frac{\nu\nu\nu}{U_{\infty}}}\frac{ \partial u}{\partial y}\) (quantity related to shear stress). The boundary conditions are set considering the laminar flow on a flat plate, the no-slip condition, and free-stream velocity outside the boundary layer. Hence we have

\[f(0)=0\hskip 14.226378ptf^{\prime}(0)=0\hskip 14.226378ptf^{\prime}(\infty)=1\] (II.4)

Since the value of \(f^{\prime\prime}\) at \(\eta=0\) is unknown, it cannot be considered an initial value problem.

A power series solution to the boundary layer equation of flow across a flat plate was proposed by Blasius [21]. Schmidt and Beckmann [22] and Ostrach [23] conducted the most important work on the topic, conducting theoretical and experimental research on the free convection flow of air around a vertical flat plate under the influence of gravity. Boyd solved the equation using an analytical series solution technique [24]. Nowadays, with the availability of computers, we can obtain a numerical solution to this equation and calculate it with a very high degree of accuracy. To solve the Blasius equation numerically, a possible method is to use the shooting algorithm to find what value would satisfy the boundary condition at \(\eta=0\). The first step is to guess a value at the wall. Then solving the ODE along the non-dimensional coordinate until the first derivative of \(f\) stops changing. This checks if the first derivative satisfies the given boundary condition; instead, adjust the guessed value to decrease or increase if the first derivative is higher or lower than one. The algorithm is repeated until the boundary condition is finally satisfied to reach the final solution. Howarth [25] found the solution using a numerical method that accurately predicted the value of \(f^{\prime\prime}(0)\) to be 0.332.

This differential equation is a direct representation of the velocity profile inside the boundary layer. Once the solution is obtained, the velocity components can be calculated using these relations. Hence the boundary layer thickness can be calculated located at the point where the velocity is 99\(\%\) of the free stream velocity. It was found that this occurs approximately at \(\eta\) equal to 5. The displacement and momentum thickness can also be estimated using the solution. The wall shear stress based on the velocity gradient at the wall can be estimated, and the friction and drag coefficients can be calculated. Note that the expression is only accounting for one side of the plate. Comparing the exact solution with the approximation solution from the integral analysis is interesting. The relations are close to each other. Indeed the integral analysis is within 10\(\%\) of the exact solution; unlike Blasius' solution, these values were obtained without complicated math. The Blasius solution provides a self-similar solution meaning that the solution is the same if the independent and dependent variables of the governing equations are appropriately scaled. This can be seen in comparing experimental results with the Blasius solution. Different Reynolds numbers (\(Re\)) provide the same profile when the variables on the two axes are appropriately scaled.

The solution to the Blasius equation can be found in [26, 27], where the accurate benchmark results of the Blasius boundary layer problem using a leaping Taylors series that converges for all real values. There have been different methods to solve the Blasius equation, such as the Topber transformation, which is executed using inverse transformation, [28, 29, 30]. Runge-Kutta, incorporated with the shooting method, finds the solution numerically, while Adomian Decomposition Method [31] finds the solution analytically.

## III Proposed methodology

The ODE representing the Blasius equation derived from the PDE is considered in the present work. We follow a similar methodology proposed for the general form of generalized non-linear ODE by Raissi [11]. The workflow diagram of the present method is shown in Fig.1.

The input values \(\eta\) are discrete and are equally distributed between the boundary points \(\eta_{0}\) and \(\eta_{m}\). The beauty of PINNs lies in loss function as it incorporates the boundary conditions and the differential equation and thus includes the physics information. The loss function is the total loss from given ODE \(L_{o}\), initial conditions \(L_{i}\), and boundary \(L_{b}\) conditions.

\[L=L_{o}+L_{b}+L_{i}\] (III.1)

\[L_{o}=\sum_{\eta}[\hat{f}^{\prime\prime\prime}+\frac{1}{2}\hat{f}\hat{f}^{ \prime\prime}]^{2}\] (III.2)

\[L_{i}=[\hat{f}(\eta_{0})-f(\eta_{0})]^{2}+[\hat{f}^{\prime}(\eta_{0})-f^{ \prime}(\eta_{0})]^{2}\] (III.3)

\[=[\hat{f}(\eta_{0})-0]^{2}+[\hat{f}^{\prime}(\eta_{0})-0]^{2}\]

\[L_{b}=[\hat{f}(\eta_{m})-f(\eta_{m})]^{2}=[\hat{f}(\eta_{m})-1]^{2}\] (III.4)

From the definition of the loss function, it is clear that the total loss will be zero if the function \(f(\eta)\) is exact. Here we update the weights of the neural network in each iteration such that the loss function is minimized. Supplying the independent input values representing discrete spatial coordinates \(\eta\) ranging from zero to \(\eta_{m}\) into the neural network is sufficient to solve the ODE via PINNs. The neural network maps the input \(\eta\) to \(f(\eta)\)

Fig. 1: Workflow diagram to solve for \(f(\eta)\)
which is the estimated solution to the stream function. Unlike standard neural network techniques that approximate the value of \(f(\eta)\) in a heuristic manner from sample output values, PINNs obtain a solution function that minimizes the loss function, which is a combination of the differential equation and boundary values.

Several combinations of the number of hidden layers and the number of nodes in each layer were tested on a trial-and-error basis. We see that the results obtained are at par with that of [27] when initializing the neural network using two fully connected hidden layer structures and setting each hidden layer's width to 100 neurons. The learning rate set in this case is 0.96. Some of the relatively good results obtained for various combinations of hidden layers and nodes, along with the changes in the learning rate of the algorithms, are presented in table I. The solution \(f_{i}(\eta)\) for each of the cases in table I is graphically presented in Fig.2.

In earlier methods proposed by Lagaris [9], all the boundary conditions had to be found before moving to find the solution. In this method, it is only required to set an appropriate finite value for the \(\eta\) at \(\infty\), which we assumed to be \(\eta_{m}=8\) in our case. Although higher values could be set, we see from the literature and our results that \(f^{\prime}(\eta)\) takes the value very close to one when \(\eta=5\). Hence the justification for taking the value of \(\eta\) at \(\infty\) as \(\eta_{m}=8\). Another benefit of using PINNs is that we can considerably reduce the number of collocation points and still attain the same level of accuracy. In our method, we have considered a total of 100 equidistant collocation points of \(\eta\) between 0 and 8. The network trains with this input data and finds the derivatives using the automatic differentiation tool to calculate the loss \(L_{o}\).

The Adam optimizing algorithm is incorporated into the PINNs methodology to set the model's adaptive learning rates. It adds momentum as the estimate of the first-order moment of the gradient and includes bias corrections to the estimates of first and second-order moments. The second-order method to train the network was the Limited memory Broyden Fletcher Goldfarb Shanno (L-BFGS) algorithm. The cost of memory has been decreased by avoiding the hessian approximation of BFGS and replacing it with an identity matrix. The loss function is optimized using ADAM and L-BFGS until it converges.

The method was extended to the negative axis setting \(\eta_{0}\) to -5.69 and \(\eta_{m}\) to 7. The value of the \(f^{\prime\prime}(0)\) found from the previous program is incorporated into the existing loss function. Our purpose in doing this was to check if our method could capture the singularity of the solution function mentioned in the literature.

## IV Results and discussion

The proposed PINNs method can reduce computational time in solving the Blasius equation over conventional CFD techniques as it reduces the burden of finding the value of \(f^{\prime\prime}(0)\) and then solving the initial value problem. The neural network constructed with two hidden layers with 100 neurons each is considered for our results and discussions. The input of \(\eta\) ranges from 0 to 8 (\(\eta_{m}\)). The learning rate for adam optimizer is assigned the value of 0.96. From Fig.3, it is clear that the results obtained from PINNS perfectly match with numerical techniques' results [26].

In the present method, the value of \(f^{\prime\prime}(\eta)\) at 0 is not considered for finding the solution. Further, once \(f(\eta)\) is obtained, we can also find the functions corresponding to \(f^{\prime}(\eta)\) and \(f^{\prime\prime}(\eta)\). The value of the second derivative at the wall

Fig. 3: Estimated results \(f(\eta),f^{\prime}(\eta),f^{\prime\prime}(\eta)\) obtained from PINNs and those obtained by numerical methods

Fig. 2: Solutions obtained for various widths and depths of NN
is observed from the PINNs method to be 0.33165, which aligns approximately with those found by Howarth [25]. The loss generated from this method is 1.67x10\({}^{-6}\) that is of order \(10^{-6}\) indicating the results to be accurate with sufficiently less error value. As proposed by Gaurav Pandey and Ambedkar Dukkipati [32], increasing the width (number of neurons) and depth (number of hidden layers) of the network can be beneficial. But in our case, we see that the function \(\hat{f}(\eta)\) tends to become more vulnerable to a high amount of bias by choosing a deeper network. Thus, the initial conditions must be chosen much more accurately for the loss function to converge to the global minimum. Further, fixing the width and depth of the network to 100 and 2, respectively, depicts accurate results; hence, the need for a more complex network was overridden.

The velocity components are directly proportional to \(f^{\prime}(\eta)\). Thus, the proposed method's solution helps in mapping non-dimensional coordinate \(\eta\) with \(f^{\prime}(\eta)\). This is the velocity profile inside the boundary layer. The values of velocity components can also be calculated from \(f^{\prime}(\eta)\) and \(f^{\prime\prime}(\eta)\) From the graph, it is visually appealing that when \(\eta\) approaches five, the \(f^{\prime}(\eta)\) goes to 1. Hence the boundary layer thickness is calculated where the velocity is 99\(\%\) of free stream velocity.

In 1999, Boyd [24] suggested that convergence of his power series method to solve the Blasius equation is limited by a singularity on the negative \(\eta-\)axis at \(\eta=-5.6900380545\). Anil and Milan [27] considered a leaping Taylor's series solution for the Blasius equation to overcome this singularity. In the second part of our investigation, we extended the collocation points to the negative portion to check for the solution of the Blasius equation on the negative \(\eta-\)axis. The solution obtained from PINNs also shows that the function \(f(\eta)\) increases rapidly near \(\eta=-5.7\) signifying the presence of a singularity near the point. A graph of \(f(\eta)\) including the negative \(\eta-\)axis is shown in Fig.4.

## V Conclusion

PINNs use knowledge of the governing equation in deep learning and find a solution to the differential equation by minimizing a loss function, including the physics information. The proposed methodology using PINNs is free from mesh generation, which is an integral part of conventional CFD techniques. Here, we have obtained the solutions to the Blasius equations, which agree with all the numerical methods in the literature. Further, the solution captures the singularity mentioned while solving the differential equation analytically.
