Turbulence model augmented physics informed neural networks for mean flow reconstruction

Experimental measurements and numerical simulations of turbulent flows are characterised by a trade-off between accuracy and resolution. In this study, we combine accurate sparse pointwise mean velocity measurements with the Reynolds-Averaged Navier-Stokes (RANS) equations using data assimilation methods. Importantly, we bridge the gap between data assimilation (DA) using Physics-Informed Neural Networks (PINNs) and variational methods based on classical spatial discretisation of the flow equations, by comparing both approaches on the same turbulent flow case. Firstly, by constraining the PINN with sparse data and the under-determined RANS equations without closure, we show that the mean flow is reconstructed to a higher accuracy than a RANS solver using the Spalart-Allmaras (SA) turbulence model. Secondly, we propose the SA turbulence model augmented PINN (PINN-DA-SA), which outperforms the former approach by up to 73% reduction in mean velocity reconstruction error with coarse measurements. The additional SA physics constraints improve flow reconstructions in regions with high velocity and pressure gradients and separation. Thirdly, we compare the PINN-DA-SA approach to a variational data assimilation using the same sparse velocity measurements and physics constraints. The PINN-DA-SA achieves lower reconstruction error across a range of data resolutions. This is attributed to discretisation errors in the variational methodology that are avoided by PINNs. We demonstrate the method using high fidelity measurements from direct numerical simulation of the turbulent periodic hill at Re = 5600.


I. INTRODUCTION
Whilst the Navier-Stokes equations describe accurately the evolution of fluid flow in space and time, their direct numerical simulation (DNS) is intractable in turbulent regimes.As a compromise, industrial simulation is dominated by Reynolds-Averaged Navier-Stokes simulations (RANS), which govern the time-averaged flow quantities instead of unsteady time-varying values.However, the RANS equations are not in closed form and solving for these first-order statistics (mean flow quantities) requires knowledge of the second-order velocity statistics, known as the Reynolds stresses.This closure problem is tackled by use of turbulence models, most commonly Boussinesq linear eddy viscosity models (LEVM), including the Spalart-Allmaras (SA) model [1], k − ω [2] and k − ε [3] models.Whilst solutions of the RANS equations are computationally tractable, they are also less accurate than scale-resolving simulations.Although the turbulence models fully close the system of equations for the mean flow quantities, they typically rely on empirical estimates and parameter tuning, lacking accuracy and generality, as quantified in a review by Xiao et al. [4].On the other hand, experimental methods may provide valuable information about real-world flows, which are nevertheless generally limited in terms of spatio-temporal resolution and do not give access to a full flow description.Error from sources, such as imperfect test parts and sensor noise may also impact measurement reliability.Particle Image Velocimetry (PIV) is often used to measure flow velocity but its resolution may be hardware restricted and limited to specific planes of interest (2D planes in a full 3D field), whilst pressure measurements are difficult to obtain in the bulk flow without being intrusive.
In conjunction with the large growth in the use of data-based methodologies to tackle a wide range of fluid mechanics problems [5,6], data assimilation methods have enabled the augmentation of low-fidelity numerical simulations with experimental data, to reconstruct full and accurate mean flow descriptions from the latter that obey the RANS equations.Such methods thus allow, among others, to extract information not present in the original dataset, interpolate missing data by super-resolving the flow fields (such as coarse PIV fields), and infer missing fields (such as pressure and Reynolds stresses).Different methods exist for tackling such a reconstruction problem: variational methods [7,8], Bayesian inference [9,10], real-state observer [11,12] and, amongst various machine-learning-based approaches, Physics-Informed Neural Networks (PINNs), which have recently been used in conjunction with RANS [13][14][15][16].As further detailed in the following, both variational approaches and PINNs formulate data assimilation as an inverse problem, and more specifically as an optimisation one, where one aims to minimise the discrepancies between data and the reconstructed mean flow.However, variational approaches and PINNs vastly differ in the solution method of such an optimisation problem, such that one might expect significant differences between their respective results.
As above-mentioned, the variational approach formulates data assimilation as an optimisation problem, where the discrepancies between the reconstructed flow and data are minimised.This minimisation is performed in conjunction with the strong imposition of the flow governing equations, here the RANS equations, through the adjoint method.Both the RANS equations and their adjoint counterpart are solved based on classical Computational Fluid Dynamics (CFD) discretisation methods, such as the finite-element method that is employed here.Foures et al. [17] notably applied variational data assimilation method to super-resolve the mean flow for a laminar 2D circular cylinder.By providing sparse, synthetic experimental data (data that has been generated using high-fidelity computational methods but treated as experimentally generated) on a rectangular grid, the mean flow was successfully super-resolved.Notably, a Helmholtz decomposition is applied to the divergence of the Reynolds stress tensor, separating it into a potential component and a divergence-free component, the latter forming the control vector in the optimisation procedure.This work was expanded by Symon et al. [18] using the same variational data assimilation method but applied to a turbulent mean flow around an aerofoil at Re = 13500 using experimental two-component planar PIV data.The work by Franceschini et al. [8] extends this work further.Instead of the Helmholtz decomposition, the Reynolds stress is broken down into a eddy viscosity term and a non eddy viscosity forcing term.The eddy viscosity term is solved by including the SA model to the governing equations, whilst the variational data assimilation is used to find the non eddy viscosity forcing and also a corrective forcing within the SA equation itself.This is intended to reduce the corrections required in regions where the LEVM is insufficient.Whilst the present study focuses on flow reconstruction, it may still be worth mentioning that the above-mentioned corrections that are identified by variational data assimilation may then form a dataset to build machine-learning-based turbulence models, as performed in the field-inversion machine-learning approach [19][20][21].
As an alternative to variational data assimilation, we here also consider the PINN approach first presented by Raissi et al. [22].In conventional neural network approaches, models are optimised by minimising a cost function, which is purely driven by error to data.The PINN framework enforces weakly the governing equations, by adding the residual PDE errors to the cost function.This is similar to the aforementioned variational approaches, although PINNs only softly constrain the governing equations.Constraining the optimisation with fully or partially known governing equations (i.e.Navier Stokes or RANS), in addition to measurement data, ensures that solutions are physical, improving model performance, reduces the possible solution space during optimisation and the dependence on high volumes of data to train models.This approach has already shown great promise across a range of inverse flow problems (water wave problem in Jagtap et al. [23] and high speed and supersonic flows in Mao et al. [24] and Jagtap et al. [25]) and can even be extended to convolutional neural networks, as demonstrated in Gao et al. [26] and Kelshaw et al. [27].Eivazi et al. [14] and Hasanuzzman et al. [28] used PINNs that are constrained by the RANS equations to reconstruct the flow field within a subdomain, by providing data at the domain boundaries for all flow fields (from experimental PIV data).Sliwinski and Rigas [29] approached the same problem as Foures et al. [17], albeit using the PINN framework instead.Firstly, the Helmholtz decomposition is applied to the unknown divergence of the Reynolds stress tensor (forcing vector), as used in Foures et al. [17], to reconstruct the flow from mean velocity data (first order statistics) alone.Secondly, by considering the full Reynolds stress tensor and providing both sparse first-and second-order velocity statistical data, the pressure field was also successfully extracted.A similar approach was taken by von Saldern et al. [15], using PINNs in combination with velocity data, to reconstruct the mean flow of a swirling turbulent jet.Pioch et al. [16] also uses PINNs, constrained by the RANS equations with the addition of turbulence models, to reconstruct the mean flow over a backward facing step.Additional research by Molnar et al. [30] propose a new PINN-based process for background-oriented Schlieren imaging (BOS), to reconstruct the density field reconstruction from experimental data, with greater accuracy than traditional BOS image processing algorithms.
One avenue to increase the capability of PINNs is through the development of novel PINN architectures.The PINN framework has been expanded with Conservative PINNs (Jagtap et al. [31]) and Extended PINNs (Jagtap and Karniadakis [32]), with further development by Shukla et al. [33].These methods divide the domain into subdomains, each with its own PINN, increasing the representative capacity of the framework, as discussed in [34].Additionally, Mishra and Molinaro [35] and De Rycket al. [36] complete error analysis on the use of (extended) PINNs.Alternatively, one can apply advancements in other data-assimilation techniques to the existing PINN framework, in order to improve on previous studies, as performed here.This study thus aims to demonstrate that similar strategies as proposed in variational approaches may be developed in the PINN framework, to achieve the reconstruction of turbulent mean flows from sparse mean-velocity measurements, extending the laminar PINN framework of Sliwinski and Rigas [29].
We here introduce a turbulence model-augmented PINN (RANS with SA) following an approach that is inspired by the above-mentioned variational data assimilation study by Franceschini et al. [8]: the Reynolds-stress term in the RANS equations is separated into a modelled component and a corrective one.The modelled component involves an eddy viscosity that is assumed to verify the SA equation, whilst the corrective component forms one of the PINN's output, along with the mean velocity, pressure and eddy viscosity fields.Compared to previous studies investigating PINNs in the context of RANS [13][14][15][16], the present approach exploits traditional turbulence modelling, adding physical constraints on the reconstructed mean flow and facilitating the handling of sparse data as demonstrated in the following, whilst at the same time not being restricted by modelling deficiencies (e.g.non-compliance to the Boussinesq hypothesis) through the introduction of the correction term in the RANS equations.Moreover, the similarities between the present PINN-based methodology and the variational framework of Franceschini et al. [8] allow to perform here detailed comparisons between PINN-based and variational data assimilation, which, to the extent of the authors' knowledge, has never yet been performed in the context of RANS and mean flow reconstruction, while only very few such comparisons for unsteady flow reconstruction may be found in the literature [37].
The following sections will be set out as follows.Section II will begin with presentation of the mean flow RANS equations and detail the use of the SA turbulence model, which will be used in this work.This is followed by a description of the data assimilation procedure and its implementation using PINNs or a variational approach.Section III introduces the test case and its numerical set up: a turbulent periodic hill at Re = 5600.This is followed by a description of the PINN architecture and training procedure in addition to the numerical setup for the variational data assimilation.Section IV will present the data assimilation results of two parts.The first segment presents results of the PINN mean flow reconstruction, without use of a turbulence model.This will demonstrate the basic PINN's capabilities and its limitations.The second part will present results, showing the use of turbulence model augmentation (SA) with PINNs.These results will be compared first to the baseline PINN (without turbulence model) and then compared with results from the equivalent variational-based method.Section V contains concluding remarks.

II. FLOW EQUATIONS AND DATA ASSIMILATION METHODOLOGY
In this article, the inverse problem is solved, where a mean (time-averaged) flow field is reconstructed from a set of partial observations.The observations are N m sparse high-fidelity mean velocity measurements at points x m = (x m , y m , z m ), i.e.U m,i , where U i is the i-th component of the mean velocity and the subscript m denotes the m-th sparse measurement.The reconstructed flow field is sought by solving the constrained optimisation problem, such that error to measurements, J, is minimised, where while simultaneously the reconstructed field satisfies the under-determined mean flow governing equations.Here, Ûi (x m ) is the reconstructed mean velocity at x m .Whilst the equations in this section are presented for general flows, the results demonstrate use for a two dimensional flow.Firstly, the governing mean flow equations without (baseline) and with (SA) turbulence model are presented in Section II A. Secondly, in Section II B, the details of two data assimilation approaches used to solve the optimisation are described -an adjoint-based variational approach and a PINN approach.

A. Formulation of the Reynolds stress closure
By applying the Reynolds decomposition, an unsteady flow field can be separated into a mean and fluctuating component.Through time-averaging of the Navier-Stokes equations, one can obtain the steady RANS equations.The data-driven techniques used in this study will be applied to mean flow data, which satisfy these RANS equations where P is the mean pressure field and u i ′ u j ′ are the Reynolds stress tensor terms.S ij is the mean strain rate tensor equivalent to 1 2 ∂Ui ∂xj + ∂Uj ∂xi .For a three-dimensional flow (3D), this results in a system of 4 equations with 10 unknowns (assuming no spanwise mean flow, this can be reduced down to 3 equations and 6 unknowns).Given the dependence of the mean flow solution, (U i , P ), on the terms of the unknown Reynolds stress tensor, u ′ i u ′ j , determination of these terms is critical to the performance and capability of data-assimilation methods for reconstructing the mean flow field.The remainder of Section II A will present formulations of the mean flow equations, such that data can be used to infer the Reynolds stress (closure) terms.

Baseline formulation of mean flow equations
To solve the inverse problem, an additional constraint is placed on the closure term to improve its inference.Firstly, the divergence of the Reynolds stress tensor (as it appears in the RANS equations (2b)) is considered as a forcing vector (Reynolds forcing vector).This reduces the closure term from 6 individual Reynolds stresses to 3 individual forcing terms (or 2 forcing terms from 3 stresses with no spanwise mean flow).Subsequently, a Helmholtz decomposition is applied, as in Foures et al. [17] and Sliwinski and Rigas [29], decomposing the forcing into a scalar, potential part (ϕ) and a divergence-free, solenoidal, vectorial component (f s,i ).This latter condition (divergence-free) provides an additional equation.Introducing decomposition (3) into the RANS equations (2), gives In addition to the new equation, the Reynolds forcing vector now consists of 4 terms (3 in two dimensions), in the form of a solenoidal forcing vector and a scalar potential field.In this system, there are now 5 equations and 8 unknowns (4 equations and 6 unknowns in 2D).By combining the pressure and potential terms into a single scalar field (P − ϕ), another unknown can be removed.The resultant flow field, (U i , P − ϕ, f s,i ) is composed of the mean velocity components, combined pressure and potential forcing term and the solenoidal forcing components.Providing mean velocity data alone, at discrete locations, is sufficient to close the system at those points.However, the P − ϕ term can not be separated into its constituent parts and the pressure field cannot be determined.

Turbulence-model augmentation of mean flow equations
Whilst the above approach combines data with a formulation of the governing equations to infer the unclosed divergence of the Reynolds stress tensor, the closure problem still poses a challenge away from data points, since an infinite number of potentially physically unrealisable solutions satisfy the under-determined governing equations.The following methodology is proposed to augment the aforementioned approach, by adding additional physical constraints.
When performing numerical simulations of the RANS equations, turbulence modelling is utilised to approximate the Reynolds stresses and close the equations.Typically, the Boussinesq approximation is used, where the Reynolds stresses are assumed to be isotropic and depend on the eddy viscosity ν t .However, neglecting the anisotropic component, as well as other modelling assumptions made in the calculation of ν t , introduce errors which propagate to the mean flow solution.Instead, as applied in Franceschini et al. [8], one can decompose the Reynolds forcing into a modelled, isotropic part and a corrective component which can include the anisotropic part, in addition to any isotropic component not captured by the modelled part, such that The first term represents the modelled forcing component and the second term the corrective forcing component, such that the mean flow solution matches the data (high-fidelity measurements).In this work, the Helmholtz decomposition is also applied to f c,i .By augmenting the system with a turbulence model, one can solve for ν t .This work will use the one-equation Spalart-Allmaras turbulence model, which provides a governing equation for a pseudo eddyviscosity variable ν.The actual eddy viscosity ν t may then be obtained from ν through an algebraic expression.
The SA equations are given in Appendix A 1. The decomposition (5) has a closed (modelled) part, ν t , and unclosed (corrective) part, f c,i .Equation ( 5) can then be substituted into Equation (2b) to derive the final SA augmented RANS equations, Applying the turbulence model allows an isotropic component of Reynolds stresses to be approximated using the SA equation, leaving the anisotropic (and potentially an isotropic) component underdetermined.As a large portion of the Reynolds forcing is now modelled, data-assimilation techniques only need to apply smaller corrections, compared with Equation (4), to close the solution.

B. Data assimilation techniques
Using data assimilation techniques, to solve the constrained optimisation, a solution Ûi (x), P (x)− φ(x), fs,i (x), ν(x) to the mean flow equations ( 4) or ( 6) is sought, such that J, the measurement error (1), is minimised.In the specific case where turbulence model augmentation is not applied, as in Equation ( 4), the eddy viscosity field ν(x) = 0. Two techniques are implemented to solve the constrained optimisation problem: a neural network approach, specifically the PINN approach from Sliwinski and Rigas [29] and an adjoint-based variational approach, used by Foures et al. [17] and Franceschini et al. [8].

Physics informed neural network
For a PINN, the solution ( Ûi , P − φ, fs,i , ν)(x; θ), is defined by weights and biases, θ, which define a neural network.Figure 1 shows this PINN schematic.For any neural network, the optimal set of weights and biases, θ, is sought by minimising some loss function C, In the case of PINNs, the loss function, C, is used to enforce the governing equations and the boundary conditions.Whilst these two terms are sufficient to define an optimisation for the forward problem, when using data assimilation (the inverse problem) the cost function must also include the measurement error term, J.To distinguish the use of PINNs for the inverse problem (as in this paper) compared to the forward problem, it will be referred to as PINN-DA hereon-in to emphasise the use of data.For the mean flow problem described in Section II A, the loss function can be defined as where the governing laws are evaluated at N c collocation points (subscripted by k), and boundary conditions are imposed at N b points along the boundary (subscripted by b).The data loss is subscripted by m.RAN S is the residual of the mean flow equations (( 4) or ( 6)) and BC is the error to the imposed boundary conditions.To evaluate the governing laws, the gradients in the RANS equations are calculated accurately to machine-precision using automatic differentiation (leveraging the chain rule to trace the derivatives of all the constituent operations in the mapping N (x, θ) from Figure 1).The data term, physics term and boundary condition terms are weighted by factors λ D , λ P and λ B , respectively.Selection of these weights is discussed in Section IV A.
The decomposition of the forcing into a modelled and corrective component ( 5) is not unique, as discussed in Foures et al. [17].As a result, additional regularisation is applied to ensure uniqueness of the distribution between modelled and corrective forcing during optimisation.For the PINN-DA, an L 2 regularisation is used for its effectiveness and ease of implementation.This regularisation term (also evaluated at collocation points), weighted by λ R , controls the magnitude of the corrective forcing field and subsequently the decomposition between modelled and corrective forcing.By increasing the penalisation of forcing magnitude, the importance of the modelled forcing component will grow and the governing equations will tend towards the RANS-SA equations (no corrective forcing).The effect of regularisation and importance of λ R is demonstrated in Appendix B 1. In the case without turbulence-model augmentation (effectively ν(x) = 0), λ R = 0 and the SA equation can be removed from the loss function.The total loss C function ( 7) is minimised using backpropagation by iteratively adjusting the weights and biases based on the gradient of loss with respect to the network parameters, ∂C/∂θ.

Variational assimilation
In the variational approach, the corrective forcing f s,i is adjusted such that the measurement error, J, is minimised, constrained by the RANS equations (6).The solution to this minimisation problem is given by the stationary points of the Lagrangian The inner product is defined as where a, b are two spatial fields and Ω is the volume over which the field is defined.U * i , P * , ν * are the Lagrange multipliers (or adjoint variables) enforcing the constraints.By setting the variation of the Lagrangian with respect to the direct flow variables U i , P − ϕ, ν, one can derive the adjoint equations Furthermore, by calculating the gradient of Lagrangian, L, with respect to changes in forcing, ∂L ∂fs,i , one finds ∂L ∂fs,i = dJ dfs,i = U * i .Specifically, this means the gradient dJ dfs,i used to set the minimisation direction for optimisation is the solution to the adjoint equations (11).
As with the PINN-DA-SA, regularisation is required to ensure uniqueness of the corrective forcing.Furthermore, the measurements are sparse, inducing a pointwise forcing of the adjoint momentum equations (11b) according to ), which may lead to discontinues in the adjoint field and in the gradient dJ dfs,i , and thus ultimately in the reconstructed forcing and corresponding mean-flow solution.This can be circumvented through H 1 -like regularisation of the gradient, which consists in getting a smoothed gradient dJ R dfs from the original gradient dJ dfs through the inversion of the following system where β = l 2 r and l r may be interpreted as a filter length.In the following, l r is chosen as the spacing between measurement locations.Eq. ( 12) also involves a pseudo-pressure field π which is introduced to preserve the divergencefree character of the gradient through the regularisation procedure.A complete derivation and details of the variational approach can be found in [17], [8] and [38].
It has to be noted that for high Reynolds numbers turbulent cases, the variational-based data assimilation requires the use of the SA turbulence model augmentation in order to numerically solve the discretised equations, unlike the PINNs which work (albeit with different accuracy) with and without the SA model.Previously, the laminar studies in [17] have not required this augmentation.To distinguish this approach it shall be referred to as variational-DA-SA.For comparison, the DNS dividing streamline is shown in grey.

Comparison: PINNs vs variational assimilation
Equations ( 8) and ( 9) highlight the similarities as well as key differences between the variational-DA-SA and PINN-DA approaches.Firstly, whilst the variational-DA-SA directly adjusts the corrective forcing, through the gradient dJ/df s,i , the PINN-DA indirectly adjusts this forcing through adjustment of the network parameters, θ.Secondly, by the nature of PINNs, a continuous mapping of the flow solution is obtained, which can be queried at any spatial location.Additionally, the governing equations are evaluated and enforced at discrete locations within the domain (collocation points) and the spatial gradients in this mesh-free approach are determined through automatic differentiation.In the variational-DA-SA approach, solution of the direct and adjoint equations introduces mesh dependency, where the conservation laws are enforced across the entire domain, whilst spatial gradients are calculated through discretisation (here Finite Element Method).As a result, discretisation errors are introduced, which are not present in PINNs.PINNs softly constrain the conservation laws in their continuous form whilst these are hard constraints in their discrete form for the variational-DA-SA.Lastly, the forcing regularisation method applied to the two approaches are different.The PINN-DA uses an L 2 -based regularisation which was found to be sufficient to produce a smooth forcing solution.However, for the variational-DA-SA, an L 2 approach with pointwise data enforcement leads to unsmooth, pointwise forcing, thus the use of a H 1 -based method, as discussed in Section II B 2. Appendix B contains a detailed breakdown on the effect of regularisation on PINNs and the effect of different regularisers for the variational-DA-SA approach.

III. NUMERICAL SIMULATIONS AND MODEL SETUP
The PINN and variational approaches described above are applied to a canonical turbulent test case, the flow over a periodic hill.This section will describe the details of the flow case and the numerical setup of the two methods employed.

A. The periodic hill set-up
The high-fidelity data for the periodic hill is obtained from the DNS database found in [39] using the Incompact3D solver [40] and the simulation setup is detailed in Xiao et al. [41].The DNS mean flow (time-averaged) solution, for streamwise velocity U is shown in Figure 2(a).The streamline (solution to the streamfunction) dividing the recirculation region from the rest of the flow is also shown.The periodic hill is simulated at Re = U b H/ν = 5600, where H is the hill height, U b is the bulk velocity over the hill apex and ν the kinematic viscosity.Out of the several geometry configurations in the database, this work uses the canonical configuration with the domain length given by L x /H = 3.858α + 5.142 for which the hill stretch factor is α = 1.0.At this Reynolds number, separation, recirculation and reattachment are present.The no-slip condition is enforced along the upper and lower walls, whilst periodicity is enforced in the streamwise direction.In order to maintain the required massflow and thus Reynolds number, body forcing is applied in the periodic streamwise direction.This is determined at each iteration, such that the Reynolds number matches the preset Reynolds number.
To analyse the effect of data resolution on mean flow reconstruction accuracy, the high-fidelity DNS data is spatially

B. PINN set-up
The PINNs are implemented using the DeepXDE Python library by Lu et al. [42], as in [29].This toolbox simplifies the process of building PINNs (with the Raissi et al. [22] paradigm).The code and data to run the PINN cases from this paper are available at https://github.com/RigasLab/PINN_SA.
All PINNs use the same architecture.The neural network was constructed as a multilayer perceptron with 7 fully connected hidden layers, each with 50 nodes.The input layer consists of two nodes, representing spatial coordinates, x/H and y/H.The output layer consists of varying number of nodes depending on the test case.The PINN initially contains five output nodes with two mean velocity components, the combined mean pressure and potential forcing field, P − ϕ, and two solenoidal forcing nodes.When the SA turbulence model augmentation was used, an additional output node is required for ν.Each node used a tanh activation function for its smoothness and its second order differentiability.This architecture was deemed to perform well across a range of data resolutions.A grid search method was employed to find optimal hyper-parameters (number of nodes and layers, activation function and learning rate).The effect of hyper-parameters on reconstruction accuracy are included in Appendix C. For PINNs, the weights and biases defining the network are randomly initialised using the Glorot uniform algorithm.This contrasts the variational-DA-SA approach, in which an "initial guess" (the RANS-SA solution) was used, around which the system is initially linearised.
For PINN optimisation, collocation points must be specified, in order to define where the governing laws are evaluated.In this work, the PINNs use 10000 collocations points which are distributed using a Hammersely distribution, such that there are sufficient points across all regions of the domain.This number of collocation points for the PINNs was selected after the reconstruction performance showed little sensitivity to further increases of the number of points.An additional 1000 points were distributed to the walls and periodic boundaries to enforce the boundary conditions.The collocation points (and domain) are shown in Figure 3.As in [29], a two step optimisation was used for the PINNs.Initially an ADAM optimisation is applied, followed by an L-BFGS-B phase.L-BFGS-B is a quasi-Newton method and thus its convergence depends on the initial guess.By using an ADAM optimiser first, one can reach an acceptable loss level to provide a good initial guess for the L-BFGS-B optimiser.Details on PINN convergence are discussed further in Section IV A. The PINNs were trained on one NVIDIA RTX 6000 GPU.

C. Modification to PINNs for Spalart-Allmaras augmentation
Using PINNs in combination with the SA turbulence model introduces additional complexity, which require modifications to the neural network architecture and the problem formulation to achieve convergence.
Firstly, the SA transport equation contains 1 d 2 terms in the production and destruction terms, where d is the true wall distance, resulting in singularities at the wall, where d = 0.This is not a problem for the finite element method (RANS and variational) where the residual is evaluated at a cell centre, which is not on the wall.For the PINN paradigm used, however, the governing equations are evaluated exactly at the wall boundaries and so this leads to singularities and subsequently failure to convergence.Consequently, the equation was reformulated by multiplying the SA equation with d 2 to eliminate the singularity.This new formulation can be found in Appendix A 2. Problems arising from singularities when used with PINNs have also been mentioned in [30].
The second change is to eliminate the effect of clipping negative values of ν.The negative SA equation formulation used in the variational approach has distinct equations for positive and negative ν (detailed in Appendix A 1) to ensure the transport equation is continuous across both positive and negative ν.The final eddy viscosity, ν t , is then clipped at negative values of ν and set to zero, resulting in only positive values of ν t .Convergence of the PINNs is improved when the complexity caused by the discontinuity from clipping ν t is removed.Thus, to enforce positive ν, the output of the neural network is defined as √ ν.A transformation is then applied to the output, squaring the value, to produce a positive ν.This is shown in Figure 1.

D. Variational assimilation set-up
The solutions of the RANS-SA system (6) and its adjoint (11) are obtained through a finite-element method (FEM) spatial discretisation implemented in FreeFEM++ [43].Contrary to the presentation in Section II B 2, the adjoint model is obtained following a discrete adjoint approach, where the Jacobian matrix that is associated to the RANS-SA equations is transposed.This Jacobian matrix is also used to invert the nonlinear RANS-SA equations based on a Newton method.Piecewise-linear functions that are enriched by bubble functions are used for velocity and pseudo-turbulent viscosity variables, while piecewise-linear functions are used for pressure.The FEM is known to be unstable at high Reynolds numbers.Accordingly, both streamline-upwind Petrov-Galerkin (SUPG) [8,44] and graddiv [45] formulations are here employed to stabilise the method.A low-memory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS) [46] algorithm is used to exploit the gradient dJ R dfs from Eq. ( 12) in order to perform the minimisation of the cost function J in Eq. (1).As mentioned above, the first-guess for the optimisation procedure is f s = 0, i.e. the baseline solution of the RANS-SA equations.The full optimisation procedure was carried out on 28 cores Intel Xeon Broadwell E5-2680v4, 2.4 GHz.The total number of degrees of freedom in the RANS and adjoint simulations is around 2 × 10 5 .Further details on the implementation of the variational-DA-SA can be found in [8,38].

IV. RESULTS
We compare the mean flow reconstruction from sparse mean velocity measurements for the turbulent periodic hill case using three approaches: PINN-DA-Baseline (without turbulence model), PINN-DA-SA (with SA) and variational-DA-SA (with SA).Firstly, an evaluation of the PINN-DA-Baseline approach is presented.Secondly, we will show that significant reduction in reconstruction error is achieved by augmenting the PINN-DA-Baseline approach with the SA turbulence model (PINN-DA-SA).Thirdly, we compare it directly to the variational-DA-SA showing the increased accuracy of the PINN-DA-SA relative to the variational-DA-SA.
To evaluate the mean flow reconstruction accuracy, a volume weighted L 2 error will be used and is calculated as where Ω is the volume of the domain, δΩ n is the cell volume of each point used to determine the error (as defined by the DNS grid) and N is the number of points across the domain.U n,i is the high fidelity, DNS mean velocity at x n and Ûi (x n ) is the corresponding reconstructed mean velocity.The variational solution is interpolated onto the DNS grid using a linear interpolation method.Additionally, the employed finite-element mesh corresponds to typical mesh sizes that are close to those in DNS.

A. Optimisation
Figure 4 shows the optimisation convergence of the PINN-DAs (a,b) and the variational-DA-SA approach (c).The PINN-DA loss convergence curve is decomposed into three loss components -PDE residual, data error and boundary condition error.For PINN-DA-SA the additional loss from the SA transport equation contributes also to the PDE loss.Each individual loss component appearing in Eq. ( 8) was given initial weight values λ based on the relative magnitudes of the flow variables.Thereafter, the weights were manually adjusted to find the combination with the lowest loss, as compiled in Table I.The weight sensitivities are presented in Appendix C.During the first phase of the training, ADAM optimisation is performed for 150000 iterations.The second phase applies L-BFGS-B optimisation for 300000 iterations or until the model has sufficiently converged to the floating point tolerance.Figures 4(a,b) shows that during optimisation (for both PINN-DA methodologies), the ADAM phase reduces loss from the high value caused by the initial randomised state to approximately 10 −3 .During the ADAM optimisation phase, one can observe frequent jumps in the loss.This is an effect of the adaptive step size used in ADAM.The step size is inversely proportional to the exponential average of squared (loss) gradient (as described in [47]).As the loss flattens out, the adaptive step size becomes very large, causing the model to "jump" too far to a higher loss level.The L-BFGS-B phase is then used to fine tune the optimisation, reducing total loss to below 10 −5 .
The convergence of the variational-DA-SA, seen in Figure 4(c) shows the data loss, corresponding to J/N m , where J is defined in (1).This data loss is equivalent to the corresponding curves in Figures 4(a,b).Whilst convergent, the data loss is an order of magnitude higher for the variational-DA-SA compared with the PINN-DAs.Whilst there is a vast difference in number of iterations between the PINN and variational approaches, an iteration is not equivalent across both methods.For PINNs, a single iteration involves evaluating the PINN at all collocation and data points, calculating the cost function and back-propagating the loss to update the model weights (and thus indirectly the flow solution).For the variational approach, in each iteration the direct and adjoint equations are solved after which the forcing solution is directly updated.An average PINN iteration takes approximately 0.84 s/iteration, whilst for the variational approach this is approximately 55 s/iteration.This difference in computational cost can be attributed to computational complexity, more mesh points in the variational method (than PINN collocation points), reduced number of degrees of freedom in the neural network (updating only weights and biases instead of the actual solution across the domain) and hardware.

B. PINN-DA-Baseline
Here, the PINN-DA-Baseline approach is applied to reconstruct the mean flow velocity field when coarse mean flow velocity data is available with spacing ∆L = 0.5.This is compared against the DNS and RANS-SA velocity profiles TABLE I. Loss function weights used in PINN optimisation.For PINN-DA-Baseline, the ν and λ R weights can be ignored.Additionally, this regularisation weight is only valid for the ∆L = 0.5 data resolution.For other resolutions, see Figure 13. in Figure 5.The data assimilation of the high fidelity measurements has improved the prediction of the velocity field compared to RANS-SA.The total mean velocity L 2 error is 3.60 × 10 −2 instead of the 8.99 × 10 −2 error from RANS-SA.High accuracy of reconstruction is observed in the bulk flow domain.However, high mean velocity reconstruction error is observed in areas with high velocity gradients, such as near the walls and at separation.This is most apparent at the upper wall boundary layer.This near wall error accounts for over 40% of the total L 2 error as the PINN-DA-Baseline fails to accurately capture high-gradient flow features.These trends are highlighted in the streamwise velocity contours in Figure 6(a,d).At data points, the PINN-DA-Baseline reproduces the measurements accurately.Furthermore, in spite of the higher near-wall reconstruction errors, Figure 6(a) shows improved reconstruction of the recirculation bubble, even with limited data in this area.This is most noticeable when comparing the PINN-DA-Baseline bubble from Figure 6(a) with the RANS bubble in Figure 2(b).The provision of sparse high-fidelity datapoints in the PINN-DA-Baseline approach has enabled better prediction of the shape of the recirculation region.However, the error field makes the limitation of the PINN-DA-Baseline approach more apparent, with the near wall errors particularly noticeable.
The results of PINN-DA-Baseline show that the methodology used in [29] for laminar flow can also be used for turbulent flows.In [29], the base formulation of the RANS equations (as in ( 2)) was also tested.It was shown that by providing data for both first order (mean velocity) and second order (Reynolds stresses) statistics, the PINN-DA-Baseline was able to reconstruct the pressure field for the laminar case.This approach has been demonstrated in Appendix D. With similar mean velocity reconstruction accuracy, this second approach was also able to infer the pressure field for this turbulent case.

C. PINN-DA-SA
This section contains results from the augmentation of the SA model to PINNs.First, a comparison between the PINN-DA-Baseline and PINN-DA-SA results will be presented, followed by a discussion on the improved reconstruction when the SA model is used.Finally, the results of the variational data assimilation using SA (variational-DA-SA) are also presented, followed by a parametric analysis of the data resolution.with the SA turbulence model.The high mean flow error along the upper wall has largely reduced indicating that use of the SA model has enabled a better reconstruction of the near-wall region.There are also similarly significant improvements along the lower wall both at separation but also after reattachment.The improvements of PINN-DA-SA compared to PINN-DA-Baseline corroborate well with the original design purpose and performance of the SA turbulence model, which shows improved performance for turbulent boundary layers in adverse pressure gradients.

Why PINN-DA-SA?
In order to understand the difference in performance between PINN-DA-Baseline and PINN-DA-SA, we examine the individual terms of the cost function, C, for each approach.This consists of the measurement error term (data loss) and the enforcement of the governing equations (PDE loss).
Firstly, we examine the data loss component.Inspection of Figure 6(d) reveals that the PINN-DA-Baseline reconstructs the mean velocity measurements accurately.This is comparable to the performance of the PINN-DA-SA, with the convergence plots in Figure 4(a,b) showing similar converged data loss values.Away from the data points, the reconstruction error increases substantially for the PINN-DA-Baseline.This is most apparent along the upper wall where a high error region exists above the data at y/H = 2.5.Given that the reconstruction error is low at the data points for PINN-DA-Baseline (and also similar to PINN-DA-SA), further decrease of the error to measurements, J, does not reduce the mean velocity error away from datapoints, where most of the reconstruction error lies.
Secondly, one can consider the other component of the cost function, the residual PDE error.Figure 7(a) shows the PINN-DA-Baseline momentum residual error (for x-and y-momentum equations).The PDE residual error contour shows that there is a limited correlation between mean velocity reconstruction error and residual error from the conservation laws -most apparent along the upper and lower walls.A similar analysis for the PINN-DA-SA corroborates this trend (Figure 7(b)).The total residual PDE loss from both PINN-DA-Baseline and PINN-DA-SA (seen in Figure 4) are of similar magnitudes (approx.10 −5 ) and the residual PDE fields are comparable, despite the significant differences in reconstruction accuracy between the two PINN-DA approaches.This indicates that further reduction in the residual PDE error does not lead to increase in reduction in reconstruction error and thus cannot be attributed as the cause of the difference between PINN-DA-Baseline and PINN-DA-SA.
The difference between PINN-DA-Baseline and PINN-DA-SA can be consequently attributed to the effect of the turbulence model augmentation.Due to the unclosed (underdetermined) nature of the RANS equations, there are infinite solutions to the forcing (divergence of Reynolds stress tensor) in the mean flow equations.Providing highfidelity data measurements "anchors" the PINN-DA solution at these discrete data points.However, away from them, there are still many candidate solutions to the mean flow, as the governing equations are still under-determined.These solutions all satisfy the governing equations and thus all have low residual PDE error.There are consequently no terms within the optimised cost function (Eq.( 8)) that can be minimised further, to reduce the mean velocity reconstruction error across the domain for PINN-DA-Baseline.The current optimisation for PINN-DA-Baseline is converged towards a solution that obeys the under-determined RANS equations.This is also briefly mentioned in Foures et al. [17].The SA model adds physics constraints that limit further the infinite number of admissible solutions far from the measurement points and as a direct consequence lead to improved reconstruction.

Variational-DA-SA results and comparison with PINNs
The PINN-DA-SA can be directly compared to the variational-DA-SA, which also uses SA augmentation, since both methods use the same data and physics constraints.Figures 6(c,f) show the reconstructed mean streamwise velocity and the respective absolute error for the variational-DA-SA.As with all the data assimilation results presented, the key flow features (separation, recirculation bubble and reattachment) are well reconstructed.As with both PINN-DA methodologies, the variational-DA-SA reconstructs the shape of the recirculation bubble more accurately than the RANS-SA, despite the limited measurements in this region.The reconstruction error matches the features seen in the PINN-DA results as well as in [17], with high error between data points, near walls and at the initial separation.
The key feature seen in the comparison between variational-DA-SA and PINN-DA-SA is that the latter has demonstrably lower error across all regions of flow.The highest reduction in error in the PINN-DA-SA is observed in the region between hills (y/H ≤ 1.0), specifically the recirculation and separation.Whilst all approaches have topologically similar error fields (such as the initial separation), the PINN-DA-SA demonstrates a vastly reduced error.At data points, the reconstruction error is much lower for PINN-DA approaches than the variational-DA-SA.Figure 4 shows that the PINN-DA's have measurement error of the order of 10 −7 whilst the variational-DA-SA is 10 −5 .The relative performance of the variational-DA-SA and PINN-DA-SA compared with both RANS-SA and the DNS flow can be summarised in the velocity profiles in Figure 8. Whilst both data assimilation methods represent a significant improvement over a RANS-SA simulation, the PINN-DA-SA is almost identical to the DNS results, whilst small discrepancies are still more apparent for the variational-DA-SA.
Comparing both PINN-DA-SA and variational-DA-SA to PINN-DA-Baseline show that application of SA has been successfully used to reduce error for near-wall gradients, as a result of the additional physical description provided by the turbulence model.The specific choice of turbulence model (SA) is also a contributing factor, given the aforementioned good performance in boundary layers with adverse pressure gradients.
The lower reconstruction error of the PINN-DA-SA approach compared to the variational-DA-SA approach has been confirmed also for the laminar case, in the absence of the SA model.The test case for this was the unsteady cylinder flow at low Reynolds numbers.Appendix E includes details on the laminar flow case, followed by a comparison of results.

Forcing fields
The inferred forcing given by Eq. ( 5) consists of solenoidal and modelled eddy viscosity (SA contribution) components.For PINN-DA-Baseline, the modelled eddy viscosity component is zero.A further analysis of the forcing terms allows isolation of the effect of SA in the accurate reconstruction of the mean flow field.
The PINN-DA-Baseline and PINN-DA-SA forcing fields shown in Figure 9 are qualitatively similar in the uniform flow region between y/H = 1.5 and y/H = 2.5, where turbulent production is low.However, the regions of the flow where PINN-DA-Baseline shows high levels of error (as seen in Figure 6(d)) result in larger forcing changes in the PINN-DA-SA solution.Introduction of a modelled component (via SA) has led to significant changes in the corrective forcing, most notably at the separation but also near the walls.These changes underpin the differences in reconstruction error between PINN-DA-Baseline and PINN-DA-SA.Constraining the solution with additional physical definition (via SA) has provided additional modelled forcing in regions where data-assimilation methods with the under-determined RANS had the highest errors.Whilst the PINN-DA-Baseline had converged to a poor solution, the PINN-DA-SA has reduced the "amount" of under-determined forcing resulting in lower reconstruction errors.However, the addition of the corrective term, which is still unclosed, has allowed the data-assimilation approach to correct the SA model.
The PINN-DA-SA forcing is different to the variational-DA-SA forcing, in accordance with the differences in mean flow reconstruction.In general, the PINN-DA-SA forcing is closer to the variational-DA-SA forcing than the PINN-DA-Baseline forcing.However, there are several differences in forcing, appearing in regions with the largest differences in reconstruction.This is most apparent in the significantly less smooth forcing around the separation in the variational-DA-SA forcing compared to PINN-DA-SA.The variational-DA-SA also exhibits a high forcing concentration in f 2 forcing along the rear hill, centred around a datapoint.This datapoint (in the ∆L = 0.5 case) has shown to be particularly sensitive.Whilst the gradient regularisation procedure employed in the variational approach ensures that each update of forcing is smooth, it does not guarantee that the accumulation of updates and thus the final forcing will also be smooth.Whilst other regularisation strategies were investigated (see Appendix B 2), they led to poorer reconstructed mean velocity fields.
To evaluate the accuracy of the data assimilation methods, the forcing is also compared with the true DNS forcing.As the potential forcing field is absorbed into the pressure term and thus is inseparable, we instead compare the curl of forcing, to remove the contribution of the potential forcing (since ∇ × ∇ϕ) = 0).
A direct comparison can be performed then against the curl of the DNS forcing The curl of the inferred forcing for all three reconstruction methods employed here is shown in Figure 10.Best agreement with the DNS is achieved for PINN-DA-SA, which aligns with the most accurate mean flow reconstruction.The PINN-DA-Baseline differs from DNS in several key regions -the initial separation, the recirculation, the rearward hill apex and along the walls.Differences between the variational-DA-SA and PINN-DA-SA forcing (and curl) correlate with differences in the error fields.These are most apparent at the separation region and along the upper and lower walls.

Data Assimilation with varying data resolution
Here, we evaluate the dependence of the reconstruction error on the data resolution of the available high-fidelity measurements.In Figure 11, we compare variational-DA-SA, PINN-DA-Baseline and PINN-DA-SA approaches, for a range of data-resolutions, ranging from ∆L = 0.05 to ∆L = 1.0.The reconstruction error increases monotonically as the (equidistant/square) spacing between data points increases.All data assimilation methods improve the RANS-SA solution, irrespective of resolution.For data spacing ∆L < 1, the turbulence model augmented PINN-DA-SA method achieves the lowest reconstruction error against all the other methods.Additionally, whilst the ∆L = 0.8 has not converged as well for the variational-DA-SA approach, the reconstruction error remains close to the neighbouring points.
Comparing the PINN-DA-Baseline to the PINN-DA-SA results shows the significant improvement achieved through the turbulence model augmentation.The gap between the two can mostly be attributed to the error near the walls and in the separation.Additionally, the rate of error increase with respect to data resolution is steeper for PINN-DA-Baseline.This suggests adding a modelled (closed) forcing component (i.e.PINN-DA-SA) reduces the effect coarsening data resolution has on reconstruction accuracy compared to PINN-DA-Baseline.As data resolution becomes finer, the two PINN-DA approaches converge.This is expected -as more data is provided, there are more points at which the PINN-DA-Baseline can close the solution and thus the advantage the turbulence model augmentation provides is  13)) versus data resolution for the three reconstruction methodologies (solid lines with symbols).All three approaches improve the RANS-SA predictions (dashed line). diminished.
Comparing the PINN-DA-SA versus the variational-DA-SA solution shows that, for most of the resolution range, the PINN-DA-SA also demonstrates consistently lower error at equivalent resolutions.In fact, the PINN-DA-SA at ∆L = 0.8 matches the reconstruction accuracy of the variational-DA-SA at ∆L = 0.3.However, at the extreme end of coarse data (∆L = 1.0), whilst the PINN-DA-SA is still a significant improvement over the PINN-DA-Baseline (and also RANS-SA), the variational-DA-SA approach has a lower reconstruction error.This suggests that the PINNs have a greater dependence on data for good performance.As data resolution reduces (spacing increases), the optimisation problem converges to solving the direct RANS-SA problem, for which the numerical framework within the variational-DA-SA method is more robust i.e.PINNs are less effective forward solvers.Algorithmic improvements for PINNs, including adaptive weighting algorithms, can potentially provide further improvements of the PINN approach (as demonstrated by [48]), which is currently under investigation.On the opposite end, for fine data resolutions, the variational-DA-SA is less accurate than the PINN-DA-Baseline.
One may attribute the difference between the reconstruction accuracy of the PINN-DA-SA and variational-DA-SA to the differences in the way the RANS equations are enforced in the data assimilation procedure.Concerning the variational-DA-SA approach, it may be first noticed that the latter solves for only the weak form of the RANS equations, following the finite-element approach.In addition, as mentioned in Section III D, stabilisation terms are added in the governing equations.Then there is the process of spatial discretisation through, among others, the choice of shape functions, which will determine the order of spatial accuracy.In contrast, PINNs enforce the RANS equations in their original, strong formulation, and in a continuous way.Moreover, the effect of applying discrete pointwise forcing at measurement locations poses a greater challenge for the variational approach than PINNs and thus the method of regularisation has a greater impact on reconstruction performance.Appendix B 2 contains analysis of regularisation choice.By avoiding both these sources of error, PINN-DA-SA achieves a lower reconstruction error than the variational-DA-SA.This is evidenced by the reconstruction error for finer data resolutions.Even as data resolution increases, the aforementioned discretisation error remains.Consequently, the reconstruction error of the variational-DA-SA reduces at a much slower rate than both PINN-DAs.

V. CONCLUSION
Three approaches to turbulent mean flow reconstruction from sparse data measurements have been presented and compared.Firstly, it was shown that for all methods, inferring the closure of the underdetermined RANS equations or the corrective closure of the RANS-SA equations by assimilating high-fidelity sparse data improved the mean flow reconstruction compared to RANS-SA predictions.
Secondly, the use of Spalart-Allmaras turbulence model, as used in [8] (in the variational-DA-SA approach) was proposed for use in PINNs, particularly in the context of turbulent flows.The SA turbulence model augmented PINN (PINN-DA-SA) proved to be the most accurate flow reconstruction method followed by the equivalent approach albeit with a variational formulation.Comparing the PINN-DA-Baseline (PINN without turbulence SA model) with PINN-DA-SA showed the isolated effect of turbulence model augmentation.The PINN-DA-SA performs significantly better in regions with complex flow features, such as near walls and at separation.Injection of additional physical constraints in the form of a turbulence model provides an additional modelled component to the underdetermined Reynolds forcing, reducing the complexity of the optimisation and narrowing the solution space for the optimisation.A corrective forcing term is then tuned to correct the rest of the closure.
Thirdly, the variational-DA-SA approach, whilst an improvement over the PINN-DA-Baseline, has higher errors compared to the PINN-DA-SA, most notably off the separation.These observations for mean flow reconstruction accuracy are also associated with equivalent error in Reynolds forcing.The difference in variational-DA-SA and PINN-DA-SA can be attributed to the effect of discretisation, which is avoided by PINNs.These observations are consistent across different data resolutions.
Given that PINNs is a growing field of research, there are potentially many avenues to improve the reconstruction capability further.The most apparent example is hyper-parameter tuning, such as loss weighting.In this study, a static weighting method was used with manual iterations to tune them.This is a laborious exercise, particularly as systems become more complex with more boundary conditions and less uniform topologies.Adaptive weighting algorithms [48][49][50] or novel architectures, such as Competitive PINNs [51], would both simplify this problem and also generalise the process to selecting weights.This has been demonstrated in Appendix C.

Importance of regularisation for variational DA
For the variational DA, the choice of regularisation (of solenoidal forcing f s,i ) was more important.Three approaches were trialled: no regularisation; L 2 regularisation; H 1 regularisation.
Firstly, without regularisation, whilst data is well enforced (as seen in Figure 14(a)) even compared with the final H 1 based approach, it results in a non-physical solenoidal forcing that appears strongly pointwise, clearly seen in Figure 14(d), as large corrections are applied at the data locations.The final reconstruction is thus less physical.Additionally, regularisation is needed to enforce convergence to a unique decomposition.
For the L 2 regularisation approach, error is higher than either H 1 or no regularisation methods.This is due to the inefficiency of the approach in the variational framework.To achieve sufficient smoothing of the solution, avoiding the problems of the none regularised approach, very large penalisation weighting is needing.However, this comes at the cost of poor data enforcement, as seen in Figure 14(b).
For the variational-DA-SA, H 1 regularisation was selected to smooth out the solution.As demonstrated, for the variational approach, H 1 regularisation is a much more suitable approach.It allows much better matching of the data (as compared with L 2 regularisation), whilst also giving a more physical solution.As a result this was selected as the final regularisation method.
For PINN-DA-SA, good data enforcement and a smooth non-pointwise solenoidal forcing was always ensured.Regularisation is exclusively used to control the distribution between modelled and corrective forcing.Given its computational simplicity, L 2 was selected for PINN-DA-SA.This was not the case for variational-DA-SA and regularisation was needed to smooth the solution in addition to the distribution of modelled and corrective forcing.results, the architecture, as discussed in Section III B, was selected.A similar analysis was completed for the loss weighting.A summary of the key results are shown in Table IV.Repeats of the final model were found to converge within 3.0% of the loss presented in Section IV A (after 3 repetitions).A tanh activation was selected due to the lowest value of cost function.Whilst an initial learning rate of 1 × 10 −3 demonstrated a slightly lower final loss value, this model was less repeatable than the chosen learning rate (5 × 10 −4 ).Similarly, the model architecture (nodes and layers) was selected for demonstrating the lowest final loss value.Whilst adding more layers and nodes performs similarly, the additional model complexity increases computational cost.Although both approaches achieve similar reconstruction error, this second approach requires second order statistics (Reynolds stress measurements) along with first order ones (mean flow field).However, this formulation can be used to reconstruct the pressure field, as shown before for laminar cases [29].This is demonstrated in the turbulent case, as seen in Figure 16(b,c).This has a high potential for practical applications in order to obtain the pressure field without performing intrusive measurements.Examination of the RANS equations (2) shows how provision of mean velocity and Reynolds stress data at discrete points defines all terms in the equations except for the pressure term.At these measurement locations, the equations are closed and the resultant neural network can infer a pressure field, both at data points but also away from them as well.This resultant pressure field has highest reconstruction error in the regions with highest pressure gradients in the flow, matching the trend seen with mean velocity error.The pressure reconstruction error is less correlated with the location of data points, than the mean velocity.It should be highlighted that the neural network pressure output must also include an integration constant (as the RANS equations only contain derivatives of pressure).To determine this constant and recover the actual, unshifted pressure field, one needs to fix the pressure at one point.This can be done during training (by providing pressure data at a single point) or as in this case, the pressure field is shifted in post processing.Taking the pressure data at a singular point, one can find a constant, such that by adding it to the model pressure at that point, the predicted pressure exactly matches the measurement.This constant can then be applied to the full solution.

FIG. 1 .
FIG. 1. Structure of the Physics-Informed Neural Network.A continuous mapping between spatial coordinates and flow variables is approximated by a deep neural network, denoted by N , defined by a set of weights and biases, θ.High-fidelity coarse data measurements and the (turbulence model augmented) RANS equations constrain the training of the network.k, m subscript the physics and data loss, respectively.i represents the i-th component of a vector.

FIG. 2 .
FIG. 2. High-fidelity (DNS) solution of the flow over a periodic hill compared with the RANS solution using an SA model.The dividing streamline between the recirculation region and rest of the flow is shown.(a) Streamwise mean velocity -DNS; (b) streamwise mean velocity -RANS; (c) absolute error of streamwise velocity relative to DNS solution (Uerr = UDNS − URANS).For comparison, the DNS dividing streamline is shown in grey.

FIG. 3 .
FIG.3.Turbulent periodic hill domain used for data assimilation, with ∆L = 0.5 data resolution (green crosses) used to compare the variational and PINN approaches.The collocation points, where the PINN evaluates residual error, are in red and were sampled across the domain using a Hammersley distribution.

Figure 6 (FIG. 6 .
Figure 6(a,d)  shows the PINN-DA-Baseline and Figure6(b,e) the PINN-DA-SA reconstruction results.Streamwise velocity and error contours relative to the DNS solution are shown, with similar trends holding for the vertical velocity component (not shown here).Whilst the error field topologies are qualitatively similar, the magnitude error is significantly reduced.For PINN-DA-SA, the recirculation bubble is almost indistinguishable from the DNS one.As a consequence of introducing the Spalart-Allmaras model, reconstruction error has reduced by 63% (ε 2 = 1.35×10 −2 ).For the bulk flow region above the hills and away from the walls (between y/H = 1.5 and y/H = 2.5) and at datapoints, the PINN-DA-SA and PINN-DA-Baseline show similarly low reconstruction errors.This is expected in this region given the lack of mean shear and production of turbulent kinetic energy.Comparing the PINN-DA-Baseline error against the PINN-DA-SA highlights the effect of augmenting the PINN

FIG. 7 .
FIG. 7. PDE residual error (x-momentum and y-momentum).(a) PINN-DA-Baseline and (b) PINN-DA-SA.The residual error between approaches are comparable despite the considerable improvements in mean velocity reconstruction error demonstrated by the PINN-DA-SA.

FIG. 11 .
FIG.11.Volume-weighted L2 mean velocity reconstruction error (as defined in Eq. (13)) versus data resolution for the three reconstruction methodologies (solid lines with symbols).All three approaches improve the RANS-SA predictions (dashed line).

FIG. 16 .
FIG. 16.PINN-DA-Baseline mean flow reconstruction.(a) Comparison between the DNS, PINN-DA-Baseline and PINN-DA-Baseline (with Reynolds stress data) velocity profiles for the U component.(b) Pressure field inferred from the PINN-DA-Baseline approach from sparse mean velocity and Reynolds stress data.(c) Absolute mean pressure error from PINN-DA-Baseline pressure reconstruction.

TABLE II .
Sensitivity of final model loss function, C, to activation functions and learning rate on a model with 7 layers and 50 nodes/layer relative to final choice (tanh and 5 × 10 −4 learning rate)

TABLE III .
Sensitivity on final model loss function, C, to architecture parameters (number of layers and nodes) relative to final architecture (7 layers, 50 nodes).

TABLE IV .
Sensitivity on final model loss function, C, to perturbations on weights from TableI

TABLE V .
Weights used in PINN optimisation when both mean velocity and Reynolds stress data are provided