Physically constrained eigenspace perturbation for turbulence model uncertainty estimation

Aerospace design is increasingly incorporating Design Under Uncertainty based approaches to lead to more robust and reliable optimal designs. These approaches require dependable estimates of uncertainty in simulations for their success. The key contributor of predictive uncertainty in Computational Fluid Dynamics (CFD) simulations of turbulent flows are the structural limitations of Reynolds-averaged Navier-Stokes models, termed model-form uncertainty. Currently, the common procedure to estimate turbulence model-form uncertainty is the Eigenspace Perturbation Framework (EPF), involving perturbations to the modeled Reynolds stress tensor within physical limits. The EPF has been applied with success in design and analysis tasks in numerous prior works from the industry and academia. Owing to its rapid success and adoption in several commercial and open-source CFD solvers, in depth Verification and Validation of the EPF is critical. In this work, we show that under certain conditions, the perturbations in the EPF can lead to Reynolds stress dynamics that are not physically realizable. This analysis enables us to propose a set of necessary physics-based constraints, leading to a realizable EPF. We apply this constrained procedure to the illustrative test case of a converging-diverging channel, and we demonstrate that these constraints limit physically implausible dynamics of the Reynolds stress tensor, while enhancing the accuracy and stability of the uncertainty estimation procedure.

As there is the need for Verification and Validation (V&V) of novel CFD methods, this paper addresses the underlying modeling rationale of this framework.Recently, we have already proposed a novel advancement in the context of the EPF, that focuses on ensuring realizable Reynolds stresses and consistency between the envisioned conceptual and the implemented computational model 12 .While the theoretical modeling structure and limitations of the eigenvalue perturbation have been exhaustively discussed 11 , the analysis of the eigenvector perturbation remains incomplete so far.
In this article, we undertake a detailed examination of the foundation and ramifications of the eigenvector perturbations.This thorough analysis of the Reynolds stress tensor's eigenvector perturbation in the context of RANS equations, enables us to show that eigenvector perturbation, as they are currently implemented, may lead to non-realizable Reynolds stress tensor dynamics.Moreover, we highlight numerical stability issues that may arise as a consequence, potentially preventing broader application of this approach.Therefore, we derive and propose a novel idea to prevent implausible Reynolds stress tensor dynamics in the current paper.

II. ACCOUNTING FOR TURBULENCE MODELING UNCERTAINTY USING THE TENSOR PERTURBATION FRAMEWORK
Despite the ongoing increase in computational resources, solving the set of Navier-Stokes equations for turbulent flows by scale-resolving simulations in the design phase for industrially relevant devices operating at high Reynolds numbers cannot be expected in the near future.As engineers and system designers are rather interested in rapid iteration cycles, the ability to make decisions based on statistical consideration of the mean flow is still industrial practice.Hence, all flow quantities can be split into a mean and a fluctuating part, according to φ = φ + φ ′ .To accommodate this Physically constrained uncertainty estimation need for compressible flows, a density weighted average (Favre-average) is performed, whereby φ = φ + φ ′′ and ρ φ = ρφ (1)   holds for all instantaneous quantities except density ρ and pressure p.In the scope of this paper, we will use the term RANS for the favre-averaged Navier-Stokes equations, although Reynoldsaveraging was initially developed for incompressible flows.The statistically Favre-averaged momentum equation following Einstein's notation convention describes the change of the mean momentum in both time and space, attributed to acting mean forces such as pressure gradients and divergence of viscous stresses (for the sake of simplicity, gravitational forces and forces due to rotating frames of reference are neglected).
Note: To shorten and simplify the notation, we denote the mean velocities by a capital letter u i → U i and omit the overline for density ρ → ρ and pressure p → p. Additionally we use x, y, z for x 1 , x 2 , x 3 in the following.
Based on Stokes' hypothesis the mean viscous stresses σ depend on the strain-rate tensor S i j = 1 2 ∂ x i and kinematic viscosity denoted as ν: In addition to these stresses, the right hand side of the equation contains unknown correlations of fluctuating velocities τ i j = u ′′ i u ′′ j , called the turbulent stresses or Reynolds stresses 33 .To close the set of equations and facilitate computational simulations, there exist numerous approximation methods.A widely used modeling assumption is the representation of Reynolds stresses as an isotropic function of the scalar eddy viscosity ν t and the mean rate of strain tensor, drawing an analogy to the representation of viscous stresses The equation mentioned above, also known as the Boussinesq approximation, ensures, that the trace of the resulting tensor is twice the turbulent kinetic energy k = 1 2 τ kk .State-of-the-art twoequation turbulence models, such as Menter's SST k −ω model 34 , typically solve additional partial differential transport equations for the turbulent kinetic energy and the turbulent dissipation rate and reconstruct the eddy viscosity afterwards to close the set of equations.The assumed linear Physically constrained uncertainty estimation relationship between Reynolds stresses and strain-rate tensor, however, is not universally valid, as already discussed in Section I. Consequently, any simulation using the Boussinesq assumption contains inherent epistemic uncertainty.The perturbation of Reynolds stress tensor's eigenspace 6 is the method of choice in order to account for turbulence modeling uncertainty on QoI.The underlying methodology is described in the following.
The symmetric, positive semi-definite Reynolds stress tensor τ i j can be decomposed into an anisotropy tensor a i j and an isotropic component Eddy viscosity based turbulence models assume that the tensorial characteristics of the anisotropy tensor are solely dictated by the mean rate of strain tensor (see Eq. ( 4)) The epistemic discrepancy in the evaluation of Reynolds stresses can be represented by the tensor Q i j , such that the true Reynolds stresses are Building upon the concept of the eigenspace perturbation approach, the structural uncertainty of the Reynolds stress tensor can be split into contributors of shape, alignment and amplitude of the tensor.Therefore, the anisotropy tensor can be represented by an eigenspace decomposition The orthonormal eigenvectors form the matrix v in while the traceless diagonal matrix Λ nl contains the corresponding ordered eigenvalues λ k .When Boussinesq approximation is used, the eigenvectors of the anisotropy tensor coincide with those of the strain-rate tensor, while the eigenvalues λ k are solely dependent on the strain-rate tensor's eigenvalues γ k and its trace Evidently, the Reynolds stress tensor features identical eigenvectors as well, however the eigenvalues of the Reynolds stress tensor are Inserting Eq. ( 8) into Eq.( 7) leads to: Because of the tensorial properties, the tensor Q i j can be decomposed into whereby ∆ describe the error terms for turbulent kinetic energy (amplitude), alignment (eigenvectors) and shape (eigenvalues).
As precisely quantifying the uncertainty of the turbulence model in representing the modeled Reynolds stress tensor is a challenging task, the developers and founders of the methodology rather try to estimate the uncertainty by sampling from possible solution space.Hence, it is not the aim to apply a correct Reynolds stress tensor τ true i j but a perturbed, physicallyrealizable one, which is called τ * i j 30 .Following the line of argument above, the EPF, considered in this work, creates a perturbed state of the Reynolds stress tensor defined as where a * i j indicates the perturbed anisotropy tensor, Λ * nl represents its perturbed eigenvalue matrix and v * in is the perturbed eigenvector matrix.Adhering to the procedure established in the majority of previously published works, there is no explicit modification of the turbulent kinetic k energy.
Instead the level of turbulence is manipulated indirectly by altering the production of turbulence due to affirmative perturbations of eigenvalues and eigenvectors, as will be clarified in subsequent sections.

A. Eigenvalue perturbation
As the components of the symmetric anisotropy tensor are bounded according to the realizability constraints 35 , the respective eigenvalues can be transformed into barycentric coordinates 36 .
By defining the vertices x 1C , x 2C , x 3C of an equilateral triangle, representing the componentiality of turbulence (three-component, isotropic limit (3C), two-component axisymmetric limit (2C) and the one-component limit (1C)) 37 , the mapping from anisotropy eigenvalues to barycentric coordinates is defined as The envisioned perturbation of the eigenvalues of the anisotropy tensor within physically permissible limits is grounded on shifting the barycentric state within the borders of the barycentric triangle 6 .Using the pseudoinverse of B, any perturbed eigenvalues are expressed through remapping where the relocated position x * results from linear interpolation between starting point x and target The relative distance ∆ B ∈ [0, 1] controls the magnitude of eigenvalue perturbation as presented in Fig. 1.Traditional eddy viscosity-based turbulence models assume that this eddy viscosity is a scalar, known as the isotropic eddy viscosity.Thus, turbulence behaves as an isotropic medium.
The eigenvalue perturbation modulates this to an orthotropic medium, where turbulence behaves differently along each eigen-direction 11 , accounting for the sensitivity of the model with respect to the anisotropic characteristics of turbulence.

B. Eigenvector perturbation
Given that LEVM rely on the Boussinesq approximation in Eq. ( 4), the Reynolds stress, the anisotropy and the strain-rate tensor share identical eigendirections, as already discussed in Section II.However, this relationship results in inaccuracies in predicting certain flows, e.g.involving flow separation and reattachment.Nevertheless, even for simple turbulent boundary layer flow there is significant misalignment between scale-resolving simulation (such as Direct Numerical Simulation (DNS)) and RANS model predicted eigenvectors of the Reynolds stress tensor 38 .
Hence, the eigenspace perturbation idea adds a perturbation to the eigenvectors.In contrast to the eigenvalues, there are no actual bounds for the orientation of the Reynolds stress tensor ellipsoid.
To address this issue, Iaccarino et al. 10 suggest to make use of the boundedness of the Frobenius inner product of the Reynolds stress and the strain-rate tensor, called the turbulent production P k Physically constrained uncertainty estimation of the turbulent kinetic energy transport equation.Based on the relationship of the strain-rate and Reynolds stress tensor for LEVM (see Eq. ( 4)), the bounds of the turbulent production term can be written in terms of their eigenvalues ψ i and γ i 39 : As the Reynolds stress and the strain-rate tensor share the same eigenvectors for LEVM, the lower bound of the turbulent production term can be obtained by commuting the first and third eigenvector of the Reynolds stress tensor that manipulates the relationship between eigenvalues and respective eigendirections.The permutation of first and third eigenvector results in a reconstructed Reynolds stress tensor based on Eq. ( 13), which is equivalent to the one obtained by rotating the eigenvector matrix v around the second eigenvector by π/2, see Appendix A. Whereas keeping the ordering of eigenvectors in case of LEVM, evidently leads to the upper limit of the turbulent production.In the subsequent section, we outline, why the eigenvector perturbation can lead to implausible dynamics of the Reynolds stress tensor combined with an unrealistically derived turbulent production term.The significant advantage of the EPF lies in its ability to generate a perturbed and realizable Reynolds stress tensor from an unperturbed one.This is accomplished by ensuring that the realizability condition, saying that the Reynolds stress tensor must be positive semi-definite, is met 35 .To illustrate, when perturbing the eigenvalues of the modeled Reynolds stress tensor, choosing ∆ ≤ 1 (see Fig. 1) inevitably leads to fulfillment of the realizability condition as the perturbed Reynolds stress anisotropy eigenvalues remain inside the barycentric triangle.Recently, we addressed an appropriate way to incorporate eigenvector perturbations in a self-consistent manner in order to obtain the desired realizable Reynolds stresses 12 .However, while the current formulation of the realizability principle is valuable, it is not comprehensive or adequate in ensuring that the evolution of the Reynolds stress, from one physically permissible state to another, remains physically plausible.Indeed, under certain conditions, the realizable Reynolds stress tensor, obtained through eigenspace perturbation, may become physically implausible leading to turbulent stress dynamics, which are rather not realizable.
An exemplary case to illustrate these conditions is the turbulent boundary layer, whereby we consider the flow to be steady, 1D and fully developed.This is equivalent to analyzing half of a symmetric infinite channel flow, as sketched in Fig. 2. Hence, by setting The diffusion based on viscous stresses has to be balanced by a source term, associated with a streamwise pressure gradient.Applying the isotropic eddy viscosity assumption (see Eq. ( 4)), the Reynolds stress tensor for 1D boundary layer flow becomes The eigenvalues ∂ y come along with the respective eigen- . By means of the eigenspace decomposition τ i j = v in Ψ nl v jl and employing the eigenvector matrix v = (v 1 , v 2 , v 3 ), the shear stress component of the Reynolds stress tensor can be reformulated Inserting the unperturbed eigenvectors and eigenvalues of the Reynolds stress tensor, as outlined above, into Eq.( 20) results in a strictly negative Reynolds shear stress component τ 12 = 1 2 (ψ 3 − ψ 1 ) as ψ 1 ≥ ψ 3 .However, when perturbing the eigenspace orientation according to the approach of Iaccarino et al. 10 , we obtain τ * 12 = 1 2 (ψ 1 − ψ 3 ) ≥ 0 as ψ 1 ≥ ψ 3 .In conclusion, the simple permutation of first and third eigenvector leads to different sign of the relevant shear stress shaping the boundary layer profile.At first glance, this already seems to violate obvious flow physics.
Nevertheless, we aim to present a conceptual explanation for this phenomenon.In order to qualitatively assess the physical relationships related to a change in sign of the Reynolds shear stress, we insert Stokes' hypothesis (see Eq. ( 3)) and the eddy viscosity hypothesis (see Eq. ( 4)) into Eq.( 18) Consequently, a change in sign of Reynolds stress component of higher momentum is transported to regions of lower momentum.This implies a countergradient transport that is physically implausible.Additionally, this phenomenon is associated with positively correlated Reynolds stresses and mean velocity gradients.Following the example of steady, fully developed 1D boundary layer, the turbulent production term P k = −τ i j will be negative, if τ 12 becomes positive because of eigenvector permutation.Such negative turbulent production denotes transferring energy from the turbulent scales to the mean kinetic energy, a process that is deemed physically implausible as well 40 .
In contrast to the eigenvector perturbation, a pure eigenvalue perturbation is incapable of inducing a change in the sign of τ 12 for fully developed boundary layer flow.Indeed, Fig. 3b additionally serves to illustrate the observation, that applying eigenvector perturbation lead to negative turbulent production and, consequently, negative effective eddy viscosity for any eigenvalue perturbation that falls within the bounds of the barycentric triangle.As depicted in Fig. 3, the absolute value of the turbulent production reaches its maximum at the one-component limiting state of turbulence, whereas it becomes zero for an isotropic Reynolds stress tensor, which is in accordance to the finding of Gorlé et al. 16 .This illustrative example demonstrates that, in the context of wall-bounded, boundary layer like flows, the suggested eigenvector permutation of first and third eigenvector can give rise to non-realizable Reynolds stress tensor dynamics in the set of RANS equations.Therefore, there is the need for a physics-based constraint that ensures not only realizable Reynolds stresses but also plausible Reynolds stress tensor dynamics.Subsequently, we derive this constraint, verify its validity and suggest its future usage within the EPF.
Physically constrained uncertainty estimation

B. Simplified derivation of realizable eigenvector perturbation dynamics for wall-bounded flows
As the second eigenvector of the Reynolds stress tensor in Eq. ( 19), is v 2 = (0, 0, 1) T , the rotation matrix for any rotation around this eigenvector simplifies to (choosing α = π/2 results in Iaccarino's permutation of first and third eigenvector 10 see Appendix A).The general rotation of the Reynolds stress tensor ellipsoid around its second eigenvector is sketched in Fig. 4a.
The objective is to derive a condition that evidently causes a change of sign for the shear Reynolds stress component τ 12 , ultimately resulting in non-realizable Reynolds stress tensor dynamics.Therefore, we formulate the rotated eigenvector matrix based on the unperturbed eigen- Hence, the resulting Reynolds shear stress based on Eq. ( 20) becomes: Consequently, Eq. ( 24) holds true for isotropic turbulence, as ψ 1 = ψ 3 and any rotation angle The relationship of the rotation angle and the shear stress component is verified by a step-by-step analysis presented in Fig. 4b.The resulting dependency of the Reynolds shear stress component is exactly the analytically derived one in Eq. (24).
Note: The rotation of the Reynolds stress tensor is symmetric to π/2, which means that any rotation around π/2 − δ results in the same tensor as any rotation around π/2 + δ .Therefore, α = π/4 is the appropriate choice as the smallest angle at which a sign change occurs.
The mean of the cosine in Fig. 4b has to be zero in order to obtain zero crossing of the turbulent production at exactly α = π/4.In other words, it is required that the maximum and the minimum presented in Eq. ( 19) are rotated by α.The resulting τ 12 and P k (see Eq. ( 17)) are evaluated subsequently.
value of the turbulent production have equal absolute magnitude but opposite signs.Equating the lower and upper bound of the inner Frobenius product Eq.( 17) leads to Thus, rotating the orthogonal eigenvectors around the second eigenvector by an angle of π/4, results in zero turbulent production if γ 1 = −γ 3 and γ 2 = 0.This conditions always holds true for fully-developed 1D boundary layers, as there is only a single velocity gradient present in the flow.However, any 2D flow featuring vanishing divergence of the velocity field does also satisfy Eq. ( 25).This means, that ∂U 1 ∂ x = − ∂U 2 ∂ y , given that U 3 is the vanishing velocity component.

C. A posteriori validation of the suggested constraint on the eigenvector perturbation
To corroborate our findings, we analyse another generic flow scenario, which is the 2D converging-diverging channel flow 41 .This test case consists of viscous walls with and without curvature as sketched in Fig. 12a.Our numerical setup uses inlet boundary conditions extracted forced using a boundary controller, which adjusts the static pressure at the outlet of the computational domain.We conducted a RANS grid independence study using a low-Reynolds resolution (y + ≤ 1) at solid walls and by applying Menter SST k − ω LEVM 34 .Based on this, we conduct analysis in post processing for the finest mesh featuring a resolution of 242x242x1 grid points (see Fig. 12b).The resulting velocity gradients, the eddy viscosity and the turbulent kinetic energy are used to determine the Reynolds stress tensor following Boussinesq's approximation (see Eq. ( 4)).According to our derivation above, the eigenvectors of these Reynolds stress tensors are rotated around the second eigenvector by α = π/4 in the entire domain as a first step.The rotated Reynolds stress tensors are composed using Eq. ( 13).Subsequently, we can compare the resulting turbulent production term (see Eq. ( 17)) after rotating the eigenvectors with the one based on the initial Reynolds stress tensor.The comparison, presented in Fig. 5, reveals a reduction in the effective turbulent production due to the rotation as expected.This observation confirms the exemplarily derived relationship of the turbulent production term with respect to eigenvector rotation of the Reynolds stress tensor.As a second step, we further validate the derivations by solving an optimization problem for achieving zero turbulent production by an eigenvector rotation of the Reynolds stress tensor given the velocity gradients of the previously performed RANS simulation of the 2D converging-diverging channel.Figure 6 shows the appropriate rotation angle α 0 relative to π/4 that would lead to zero turbulent production term.The deviations from to the derived α under idealised conditions of a 1D boundary layer flow, can be ascribed to the fact, that the flow is not fully divergence-free in the outer parts of the boundary layers.However, as the optimized α 0 differs by only 10% at maximum from π/4, we believe, that restricting the eigenvector rotation of the Reynolds stress tensor to π/4 is a reasonable choice also for more complex flows.
To sum up, we have shown through mathematical analysis, that a simple eigenvector perturbation involving permuting the first and third eigenvector, may lead to implausible, not realizable Reynolds stress tensor dynamics.Based on that, a constraint that facilitates physically meaningful Reynolds stress tensor perturbations with respect to rotation of the eigenspace has been derived for wall-bounded, boundary layer like flows.Additionally, we have substantiated the derivations by presenting illustrative proofs.In the next section, we apply the proposed eigenvector rotations in the EPF .

IV. APPLICATION OF PHYSICALLY CONSTRAINED EIGENVECTOR PERTURBATION
The idea of the EPF is to sample from possible solution space for certain QoI attributed to perturbing the Reynolds stress tensor within the discussed physical bounds.The entire framework was implemented in DLR's in-house CFD solver TRACE 42 .TRACE is a parallel Navier-Stokes flow solver that has been developed at DLR's Institute of Propulsion Technology in close cooperation with MTU Aero Engines AG.In the present work, we use the finite-volume method to discretize the compressible RANS equations.The EPF and can be subdivided in several steps within each pseudo-time step of steady simulations: 1. determine anisotropy tensor (see Equation Eq. ( 5)) 2. perturb eigenvalues of anisotropy tensor by choosing the relative perturbation strength ∆ B (cf. Eq. ( 15) and Eq. ( 16)) 3. perturb eigenvectors of anisotropy (Reynolds stress) tensor by choosing the rotation angle α.The rotation is done by rotating the eigenvector matrix around the second eigenvector Physically constrained uncertainty estimation v * = Rv with R according to Eq. (A1).
4. reconstruct perturbed Reynolds stress tensor τ * i j according to Eq. ( 13) 5. update of the viscous fluxes using τ * i j 6. update of the turbulent production term The simulations of the flow within a converging-diverging channel serve to exemplify the application of the EPF and further validate the proposed constraint on the eigenvector perturbation.We compare against DNS data by Laval et al. 41 .The two-equation, Menter SST k − ω 34 LEVM is chosen to be the baseline model in the present investigation.Hence, the uncertainty estimates presented subsequently based on the EPF can be attributed to the structural uncertainties within this particular turbulence model.As the amount of considered structural uncertainty increases with increasing eigenvalue perturbation, the most conservative estimation of the modelling uncertainty is obtained by choosing ∆ B = 1.0.Nevertheless, according to latest publications 30 , intense Reynolds stress tensor perturbations may cause numerical convergence issues.
Following the approach proposed in our previous work 12 , the relative perturbation magnitude with respect to the relative shift in barycentric coordinates ∆ B has to be adjusted as a consequence of the convergence issues.In the present study, we seek to apply a ∆ B as large as possible by steps of 0.1.Consequently, while the full Reynolds stress tensor perturbation could be used for the 2C and 1C corners, the perturbation towards the isotropic corner had to be adjusted by ∆ B < 1, as approaching the isotropic state results in a reduction of turbulent kinetic energy production.
Although, we have just derived that the maximum eigenvector rotation angle has to be α ≤ π/4, this constraint is necessary but not sufficient for practical applications.The eigenvector modification by applying α ≤ π/4 may result in states of the Reynolds stress tensor that are indeed realizable and physically plausible but still lead to numerical stability issues.Therefore, we iteratively decrease the rotation angle by fractions of 10% with respect to the maximum value of π/4.Besides examining the overall residuals and convergence of the static outlet pressure (controlled to maintain the prescribed mass flow outlet boundary condition) of each simulation, we evaluate the evolution of the the streamwise velocity.Therefore, we record iterative data at x/H ∈ [0.5, 1, every 1000 iterations and evaluate the relative error (standard deviation divided by the mean) over the last 100 snapshots.In order to distinguish between an unacceptable unstable and an acceptable converged solution, we use a maximum tol-  7.Estimated turbulence model uncertainty for the streamwise velocity inside of the convergingdiverging channel based on the EPF.U 1 0,max is the maximum streamwise velocity of the baseline simulation at x/H=0.The settings for every eigenspace perturbation of the Reynolds stress tensor can be found in Table I. erable relative error of 1.5% in each considered location.The numerically achievable perturbations leading to converged RANS results for this study of the convergence-divergence channel flow are summarized in Table I.In order to verify, that the eigenvector perturbation proposed by Iaccarino et al. 10 leads to unstable CFD simulations as a result of non-realizable Reynolds stress tensor dynamics, we have conducted one exemplary simulation, presented in Appendix C, applying eigenvector permutation without any eigenvalue perturbation.
In the subsequent section, we discuss the resulting estimated uncertainty intervals based on the eigenspace perturbation.The analysis refers to the presented QoI in Fig. 7 to Fig. 10.The estimated uncertainty for the streamwise velocity field is shown in Fig. 7. Perturbing the eigenspace of the Reynolds stress tensor has minor effect upstream of the diverging section (x/H ≈ 5), where indicate, that they are extracted at x/H = 0.The settings for every eigenspace perturbation of the Reynolds stress tensor can be found in Table I.
the baseline RANS simulation closely aligns with the DNS data.Due to the increased turbulent production at the one-and two-component limiting state of turbulence (as can be observed in the turbulent kinetic energy distributions in Fig. 10), the velocity profiles become sharper with an increased gradient at the wall.This is also reflected in higher friction coefficients in Fig. 8 at the bottom and top wall compared with the baseline simulation.Conversely, the simulations featuring more isotropic Reynolds stresses (simulation #5 and #6 with ∆ B = 0.2), lead to rounder velocity profiles and reduced friction coefficients.Larger deviations between RANS and DNS results arise, when the flow experiences the adverse pressure gradient in the diffusor section.This is further reflected in an increased sensitivity of the velocity field with respect to the shape and orientation of the Reynolds stress tensor.Due to the indirectly manipulated turbulent production behaviour, the turbulent kinetic energy evolves differently in the simulation domain for every perturbation (see Fig. 10).This is also in accordance with the described dependency of the turbulent production term on the eigenvalues in the front section of the diffusor (see also Section II A).The turbulent kinetic energy production significantly affects the size of the separation bubble due to the adverse pressure gradient.While the simulations aiming for the one-and two-component turbulence  I.
state completely suppress the separation zone at the lower wall, the simulations #5 and #6, featuring more isotropic turbulence both separate early and over-predict the reattachment length (see Fig. 8).While the static wall pressure reduction in the converging section is not affected by the eigenspace perturbation, the pressure recovery in the diffusor section shows growing turbulence model uncertainty (see Fig. 9).These uncertainty intervals on the pressure coefficient, underline the increased model-form uncertainty when it comes to adverse pressure gradient flows.The reduced turbulent production of the simulations #5 and #6 in the converging section x ≤ 5 results in increased turbulent kinetic energy in the massively separated section (see Fig. 10).As already described in the theoretical parts of this paper (see Section II B and Section III), the rotation of the eigenvector matrix mainly affects the turbulent production process indirectly.This is especially highlighted by decreased turbulent kinetic energy patterns in the upstream section of the domain.
It is noticeable, that the simulations aiming for the two-component limiting state of the Reynolds stress tensor, show the best agreement with the DNS data in the diffusor section for the considered QoI.Additionally, the DNS data are included in the turbulence model uncertainty estimates in most of the plots, which also validates the EPF to a certain extent, although the au-  I).
thors are aware of the fact, that this is not the main goal of the perturbation methodology.The interested reader is referred to the discussion on the restrictions and capabilities of the framework in Matha et al. 30 .Regarding the potentially large uncertainty intervals concerning the considered QoI, we must note that the eigenspace perturbations of the Reynolds stress tensor were chosen to be sufficiently large, allowing the CFD solver to handle it just well enough to produce a convergent This need is addressed by our investigation.
In this work, our primary focus centers on the eigenvector perturbation of the Reynolds stress tensor that has received limited attention in the literature.We systematically analyze that the eigenvector perturbation may violate Reynolds stress tensor dynamics under specific conditions.
The present study derives and introduces physics-based constraints for eigenvector perturbations, adhering to the realizability and stability of the uncertainty estimation procedure.The application of these constraints to the flow within a converging-diverging channel illustrates improved stability and accuracy in capturing the turbulent behavior.The flow characteristics of this case encompass turbulent boundary layer, separation and reattachment regions, revealing deviations of RANS and DNS results.Applying the EPF unveils significant sensitivity of the considered QoI based on Reynolds stress anisotropy and its eigenvector alignment with the strain-rate tensor.Based on the present paper and our previous research 12 , we have successfully identified challenges in the application and interpretability and proposed potential solutions.Our future investigations will focus on implications of the physics-constrained Reynolds stress tensor perturbation method for more complex engineering flows, such as those encountered in turbomachinery components.This will provide valuable insights into the method's practical utility.
Physically constrained uncertainty estimation tating v around v 2 with α = π/2 leads to: As it can be seen in Eq. (A3), the rotation result in permuting v * 1 = −v 3 and v * 3 = v 1 .This is illustrated in Fig. 11 and compared with the formerly proposed eigenvector permutation.Note: In contrast to rotating the eigenvectors, simple permutation of v 1 and v 3 evidently lead to a lefthanded oriented eigenvector matrix.The resulting Reynolds stress tensor, when reconstructed based on Eq. ( 13), is identical due to the characteristic of the spectral decomposition.In contrast, the streamwise velocity snapshots of the perturbed simulations, used in Section IV, converge over runtime of the simulations (Fig. 15 presents one example).FIG. 15.Streamwise velocity inside of the converging-diverging channel for simulation #2 (see Table I) using moderated eigenvector perturbation and eigenvalue modification towards the one-component turbulent limiting state.The snapshots are taken every 1000 iteration, while the mean U 1 and the standard deviation std(U 1 ) are determined between iterations 400.000 to 500.000.
FIG. 1. Systematic representation of the eigenvalue perturbation within the barycentric triangle and its effect on the shape of the Reynolds stress tensor ellipsoid.

1 ∂ y ≥ 0 σ 12 ≥ 0 τ 12 ≤ 0 FIG. 2 .
FIG. 2. Schematics of steady, fully developed 1D boundary layer flow FIG. 3. Comparison of the effect of eigenspace perturbation on the turbulent production term P k in case of fully developed boundary layer flow.Effect of pure eigenvalue perturbation is shown in (a), while (b) presents the effect, when combining the permutation of the eigenvectors v 1 and v 3 and eigenvalue perturbation within the barycentric triangle.

FIG. 4 .
FIG. 4. Rotation of the eigenvector matrix of the Reynolds stress tensor around second eigenvector v 2 by α.The schemtical impact of the rotation on the Reynolds stress tensor ellipsoid is shown in (a).(b) shows the effect of eigenvector rotation on the Reynolds shear stress component and the turbulent production.This plot is created based on assuming 1D boundary layer flow, as sketched in Fig. 2. The eigenvectors of τ i j

FIG. 5 .FIG. 6 .
FIG. 5. Distribution of turbulent production term P k rotated , when rotating the eigenvectors of the Reynolds stress tensor around second eigenvector by α = π/4.For better interpretability the resulting production is scaled by the unperturbed turbulent production P k .The magenta line indicates U 1 = 0

FIG. 8 .
FIG. 8. Turbulence model uncertainty based on the EPF for the friction coefficient c f = τ w / 1 2 ρ 0 U 2 1 0,max at upper wall (a) and bottom wall (b) of the converging-diverging channel.The quantities with subscript 0

FIG. 9 .
FIG. 9.Turbulence model uncertainty based on the EPF for the pressure coefficient c p = (p − p 0 ) / 1 2 ρ 0 U 2 1 0,max at upper wall (a) and bottom wall (b) of the converging-diverging channel.The quantities with subscript 0 indicate, that they are extracted at x/H = 0.The settings for every eigenspace perturbation of the Reynolds stress tensor can be found in TableI.

FIG. 14 .
FIG.14.Streamwise velocity inside of the converging-diverging channel based on pure eigenvector permutation without any eigenvalue perturbation.The snapshots are taken every 1000 iteration, while the mean U 1 and the standard deviation std(U 1 ) are determined between iterations 400.000 to 500.000.

TABLE I .
Selected turbulent target state (componentiality), ∆ B for eigenvalue and α for eigenvector perturbation of Reynolds stress tensor perturbation of flow within converging-diverging channel.
solution.On the one hand, this enables exploring the capabilities of the considered EPF, and on Uncertainty estimation in the context of RANS turbulence modeling is crucial in industrial applications as it provides a quantifiable measure of the reliability of CFD simulations.Accurate assessment of physically plausible uncertainties ensures the credibility of simulation results, allowing engineers and designers to make informed decisions under uncertainty.The Eigenspace Perturbation Framework has established itself as the only physics-based uncertainty quantification approach for turbulence model uncertainty.It has been applied to problems in Aerospace, Civil, Environmental, Mechanical Engineering.It is widely used in leading CFD software.This underscores the need to ensure Verification & Validation of this framework.However due to its newness, in depth verification of the rationale and application of this framework have not been conducted.
the other hand, it represents the most conservative estimate of turbulence model uncertainty.For upcoming design decisions considering this framework, design engineers would likely aim for a non-overly conservative estimate of turbulence model uncertainty, as they may introduce expert knowledge into the analysis.V. CONCLUSION & OUTLOOK 21 v 21 v 22 − v 23 v 21 v 23 + v 22 v 21 v 22 + v 23 v 2 22 v 22 v 23 − v 21 v 21 v 23 − v 22 v 22 v 23 + v 21 v 11 v 21 v 31 v 12 v 22 v 32 v 13 v 23 v 33