One-equation (simplified) LM transition model for SST and SA by rois1995 · Pull Request #1901 · su2code/SU2

rois1995 · 2023-01-27T18:26:29Z

Proposed Changes

One-equation (simplified) version of the LM transition model: only the intermittency γ is transported, and the transition onset is computed from local correlations instead of the Re_θt equation. It is enabled with KIND_TRANS_MODEL= LM and LM_OPTIONS= (SLM, MENTER_SLM), optionally with CROSSFLOW.

Supported combinations. Only the combinations published in the literature are allowed; any other one stops with a config error.

Turbulence model	Cross-flow	Reference
SST	off	Menter, Smirnov, Liu, Avancha, A One-Equation Local Correlation-Based Transition Model, Flow Turbul. Combust. 95, 2015 (DOI)
SST	on	Vallinayagam Pillai & Lardeau, Accounting crossflow effects in one-equation local correlation-based transition model, AIAA 2017-3159
SA	off / on	Lee & Baeder, Prediction and validation of laminar-turbulent transition using SA-γ transition model, AIAA 2021-1532 (including their Langtry and Menter-Smirnov C1 cross-flow models)

The equations are implemented as in these papers, with equation numbers in the code comments. Notes:

SST: Menter et al. (2015) use the Kato-Launder production in the k equation without the production limiter. SU2 prints a warning if KATO-LAUNDER is not in SST_OPTIONS (the SST production limiter stays active in any case).
SA: the intermittency is rescaled (γ_s, Lee & Baeder eq. 18) so that it goes to zero in the laminar boundary layer, and it multiplies the SA production and destruction (eq. 17). The constants of Re_θc are blended with the turbulence intensity (eq. 13).
CODER_SLM and MOD_EPPLER_SLM are rejected: Coder & Maughmer (AIAA 2012-672) use these correlations with the two-equation Langtry-Menter model, not with the one-equation model.
Roughness (cross-flow): the roughness height is limited to at least 1e-8 (it enters log(h/θ_t) and h/h0; the papers give no calibration limit), and the value used is printed at startup. h0 = 0.25 µm assumes a mesh in meters.

Fixes to the existing LM model (found while bringing this branch up to date):

Restart: INTERMITTENCY_SEP and INTERMITTENCY_EFF were read as solution variables, but they are not in compact restart files. Restarts are now consistent (tested round trips for the two-equation and one-equation models, compact and full).
Free-stream Re_θt: the low-Tu branch used 0.27 instead of 0.027, and the Tu ≥ 0.027 % limit was applied to only one term of the correlation.
AD: γ and γ_eff were not registered as preaccumulation inputs of the SA source.
Discrete adjoint: the adjoint solvers create no transition solver, so KIND_TRANS_MODEL= LM with a discrete adjoint crashed (null pointer). It now stops with a config error. Adjoint support for transition models would be a separate PR.

Tests: new regression tests slm_spheroid_* (serial, parallel, hybrid) for SST and SA, with and without cross-flow, on a coarse grid of the 6:1 prolate spheroid at 15°. They restart from converged transitional solutions (su2code/TestCases#207), so the transition terms are active. Temporary: .github/workflows/regression.yml uses the feature_Trans_SLM branch of TestCases until that PR is merged.

Validation (prolate spheroid, airfoils) is in progress.

Related Work

Two-equation LM model: Langtry and Menter transition model #1751, LM2015 for SA and cross-flow: Implementation of LM model for SA and cross-flow corrections (LM2015) #1810.
Grid and restart solutions of the new tests: Prolate spheroid grid and solutions for the one-equation LM transition tests TestCases#207.

PR Checklist

I am submitting my contribution to the develop branch.
My contribution generates no new compiler warnings (try with --warnlevel=3 when using meson).
My contribution is commented and consistent with SU2 style (https://su2code.github.io/docs_v7/Style-Guide/).
I used the pre-commit hook to prevent dirty commits and used pre-commit run --all to format old commits.
I have added a test case that demonstrates my contribution, if necessary.
I have updated appropriate documentation (Tutorials, Docs Page, config_template.cpp), if necessary.

🤖 Generated with Claude Code

rois1995 · 2023-01-27T18:28:50Z

Technically, this model has only one equation, the one for intermittency. Most of the functions are the same as the original LM model. Should I consider the Simplified model (SLM) as an option for the LM model to avoid duplicates? This is how I started, but I am open to discussions and suggestions.

pcarruscag

If the changes to the solver don't become too intrusive this approach sounds good to me

- Added correlations for Simplified LM. - There is a bug on the computation of grad(n*U)*n

rois1995 · 2023-01-31T17:28:58Z

I have a question: for each point in the mesh I am trying to compute the dot product between the velocity vector and the normal to the wall of the nearest point on the wall. How do I access such information? I've found that in the CPoint class I have the ClosestWall_Elem variable which stores the index of the closest element on a wall. However, when I try to assess the information with a number of cores greater than 2, it crashes. Moreover, to recover the normal of the element I perform a mean of the normals on the nodes of that element. Is there a structure that has the normals saved for each element of the primal grid? The part that I am referring to is from line 208 in CTransLMSolver.cpp .

- Modified LM_OPTIONS to include cross-flow effects: from LM2015 to CROSSFLOW

pcarruscag · 2023-01-31T19:25:16Z

See what is done at the bottom of CGeometry.cpp in CGeometry::ComputeWallDistance.
You need to communicate the normals you need first.

github-advanced-security

CodeQL found more than 10 potential problems in the proposed changes. Check the Files changed tab for more details.

WallyMaier · 2023-02-01T05:47:54Z

+   * \brief Get the index of the closest wall element.
+   * \param[in] iPoint - Index of the point.
+   */
+  inline unsigned long GetClosestWall_Elem(unsigned long iPoint) {return ClosestWall_Elem(iPoint);}


Can you add the param[out] for these please

Sure, I'll do it now.

WallyMaier · 2023-02-01T05:50:10Z

    SU2_MPI::Error("Two correlations selected for LM_OPTIONS. Please choose only one.", CURRENT_FUNCTION);
  }
+  if (NFoundCorrelations_SLM > 1) {
+    SU2_MPI::Error("Two correlations selected for Simplified model into LM_OPTIONS. Please choose only one.", CURRENT_FUNCTION);


Suggested change

SU2_MPI::Error("Two correlations selected for Simplified model into LM_OPTIONS. Please choose only one.", CURRENT_FUNCTION);

SU2_MPI::Error("Two correlations selected for simplified LM_OPTIONS. Please choose only one.", CURRENT_FUNCTION);

I'd stick with my change since the options are for the simplified model. They are not simplified options.

WallyMaier · 2023-02-01T05:50:50Z


-    /*--- Check if problem is 2D and LM2015 has been selected ---*/
-    if (lmParsedOptions.LM2015 && val_nDim == 2) {
-      SU2_MPI::Error("LM2015 is available only for 3D problems", CURRENT_FUNCTION);


is LM2015 gone? or is crossflow the same?

I have changed the option to CROSSFLOW, since I will use it also for the Simplified model.

WallyMaier · 2023-02-01T05:54:18Z

+        // This is not reported in the paper
+        //FPG = max(FPG, 0.0);
+


Suggested change

// This is not reported in the paper

//FPG = max(FPG, 0.0);

I have left it there since I do not know if I have to keep it or not.

bigfooted · 2023-02-10T16:02:42Z

I have a question: for each point in the mesh I am trying to compute the dot product between the velocity vector and the normal to the wall of the nearest point on the wall. How do I access such information? I've found that in the CPoint class I have the ClosestWall_Elem variable which stores the index of the closest element on a wall. However, when I try to assess the information with a number of cores greater than 2, it crashes. Moreover, to recover the normal of the element I perform a mean of the normals on the nodes of that element. Is there a structure that has the normals saved for each element of the primal grid? The part that I am referring to is from line 208 in CTransLMSolver.cpp .

There is the geometry toolbox for dot product and normal:
https://github.com/su2code/SU2/blob/master/Common/include/toolboxes/geometry_toolbox.hpp
just look in the code for examples, in vscode you can search for GeometryToolbox::DotProduct or GeometryToolbox::Norm

rois1995 · 2023-02-14T11:50:53Z

I have a question: for each point in the mesh I am trying to compute the dot product between the velocity vector and the normal to the wall of the nearest point on the wall. How do I access such information? I've found that in the CPoint class I have the ClosestWall_Elem variable which stores the index of the closest element on a wall. However, when I try to assess the information with a number of cores greater than 2, it crashes. Moreover, to recover the normal of the element I perform a mean of the normals on the nodes of that element. Is there a structure that has the normals saved for each element of the primal grid? The part that I am referring to is from line 208 in CTransLMSolver.cpp .

There is the geometry toolbox for dot product and normal: https://github.com/su2code/SU2/blob/master/Common/include/toolboxes/geometry_toolbox.hpp just look in the code for examples, in vscode you can search for GeometryToolbox::DotProduct or GeometryToolbox::Norm

The problem is more related to the finding of the wall-normal for a point within the volume mesh, not to the computations that it will be involved in.

pcarruscag · 2023-02-14T15:05:26Z

The solution I suggested didn’t work?

rois1995 · 2023-02-14T15:46:22Z

The solution I suggested didn’t work?

I still have to check if the implementation is correct but with more than two cores the code breaks. However, I found out that the wall-normal of a volume point can be computed as the normalized gradient of the wall-distance. Does this sound correct to you?

However, there is a problem: I am using the aux variables to compute these gradients, but to compute dot(n, U) I first need n, thus I cannot compute them simultaneously. Since these computations are performed in the Preprocessing of the solvers, I was thinking to compute the normal within the FLOW_SOL preprocessing and the dot(n, U) in the TRANS_SOL preprocessing since the flow solver comes before the trans solver. Is this right?

pcarruscag · 2023-02-14T16:09:19Z

It's not correct. You need to follow the pattern from CGeometry::ComputeWallDistance
First you collect all the normals that are local to the rank
Then you collect the normals across all ranks (using the NDFlatner), then you access this information by (ClosestRank, ClosestMarker, ClosestVertex)
You cannot access the local information directly because ClosestMarker and ClosestVertex may refer to a different rank.

- Added variables only for debug

rois1995 · 2023-02-24T13:15:26Z

I managed to implement the computation of grad(n*U)*n, but at the moment it is located into CTransLMSolver::PreProcessing. It seems to work with a structured mesh on a flat plate. Currently, I am testing with a 2D profile too. However, being into the PreProcessing of the transition solver, the normals are computed at each iteration, thus it is not computationally efficient if non-deforming meshes are used. I am looking into where to put it such that it is computed just if the mesh is updated (and at the first iteration of course).

Plus, I have added a whole lot of variables to the output, but they will be removed in the final version. They are just used as debug.

pcarruscag · 2023-02-24T15:11:00Z

Do what we do with wall roughness and do it in the same place (computewalldistance and store in CPoint)

The transition numerics of develop are a single edge flux kernel, CScalarFlux_TransLM: it now also handles the simplified model (one equation, intermittency only) and the solver dispatches it for 2 or 1 equations. The SLM convection and diffusion classes and their setup in the driver are removed with the old numerics. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

…ource AuxVar (wall-normal derivative of the wall-normal velocity) is only set for the Menter correlation; with the Coder and modified Eppler correlations lambda_theta, which the cross-flow term and the output use, was computed from an uninitialized value. Initialize it to zero. The preaccumulation registered only nVar = 1 turbulence variables, so omega was missing for SST, and AuxVar was not registered at all. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

In 2D the vertex normals have two components, reading three went past the end of the array. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

…production Menter et al. (2015, Eq. 24) compute the k production of SST with the Kato-Launder form mu_t*S*Omega and without the production limiter. Print a warning when MENTER_SLM is used with SST, suggesting KATO-LAUNDER when it is not selected, and noting that the SST production limiter stays active. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

LoadRestart read the separation and effective intermittencies at index + 2 and + 3, but they are PRIMITIVE outputs: in compact restarts these are the values of the next point, in full restarts other fields. They are not solution variables and are recomputed by Postprocessing at the end of LoadRestart, so they are no longer read. With the simplified model (one solution variable) RE_THETA_T and TU were in the SOLUTION group and therefore in the restart files; they are now PRIMITIVE outputs. The duplicated definitions of INTERMITTENCY_SEP and INTERMITTENCY_EFF in the solution fields are removed, they are defined with the primitive fields and written for all LM variants. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

…ta_t correlation The limit Tu >= 0.027 % was applied only to the 1/Tu^2 term, the linear term still used the unlimited value. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

The one-equation model is available with MENTER_SLM and SST (Menter et al. 2015, cross-flow of Vallinayagam Pillai and Lardeau, AIAA 2017-3159) or SA (Lee and Baeder, AIAA 2021-1532). The CODER_SLM and MOD_EPPLER_SLM correlations of Coder and Maughmer (AIAA 2012-672) belong to the Langtry-Menter intermittency equation, not to the one-equation model implemented here, so they now stop with an error. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

…1-1532) - Onset functions, F_turb and R_T = mu_t/mu as in Eqs. 4-7 (the same F_onset2 and F_onset3 as with SST). - Re_theta_c with C_TU1 and C_TU2 blended with the freestream Tu between the constants of Colonia et al. and the original ones (Eqs. 10-13). - SA equation: production times gamma_s, destruction times max(gamma_s, 0.1), with the scaled intermittency of Eq. 18 (Eq. 17). The two-equation LM coupling is unchanged. - Positivity of the implicit operator of the intermittency source (Eq. 24). - Cross-flow: the Langtry et al. stationary cross-flow criterion (Eqs. 36-43) and the Menter-Smirnov C1 criterion (Eqs. 25-35), with the cross-flow strength from the gradient of the vorticity direction. The CODER_SLM and MOD_EPPLER_SLM branch of the source, rejected at configuration, is removed. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

The SA source uses the intermittency of the LM models but did not declare it as a preaccumulation input, so its derivative was lost in AD. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

The stationary cross-flow Reynolds number of Langtry et al. contains log(h / theta_t), which diverges for HROUGHNESS= 0. Limit h from below to 0.25 micrometers, the smallest roughness for which the correlation was validated by Lee and Baeder (AIAA 2021-1532) and the reference height h0 of Vallinayagam Pillai and Lardeau (AIAA 2017-3159). Applied to the SA cross-flow of the simplified model and to the LM2015 option of the two-equation model. The roughness factor C_r of the SST cross-flow model is finite for h = 0 and is unchanged. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

…ST cross-flow model The implementation matches Eqs. 2-10 and 17-18 of AIAA 2017-3159. The sign of Eq. 2 differs from the printed one, which gives a negative critical Reynolds number; the comment now explains it. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

The papers give no calibration limit for small roughness heights, so the limit of the previous commit (0.25 micrometers, the smallest validated height) is replaced by h >= 1e-8, in the units of HROUGHNESS. It applies to log(h/theta_t) of the Langtry et al. correlation (SA cross-flow of the simplified model, and the cross-flow option of the two-equation model) and to h/h0 of the SST cross-flow model of Vallinayagam Pillai and Lardeau. The limit and whether HROUGHNESS is below it are printed at startup. h and theta_t are both in mesh length units (reference length 1), so log(h/theta_t) does not depend on REF_DIMENSIONALIZATION. Also add the missing end of line after the name of the simplified model without cross-flow. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

The two branches of the conditional had an AD expression and a double. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

The discrete adjoint solvers do not create the transition solver, but the SA and SST sources read its solution when KIND_TRANS_MODEL is set: with MATH_PROBLEM= DISCRETE_ADJOINT and LM (two-equation or simplified) SU2_CFD_AD crashed with a segmentation fault (null pointer) on the E387 case. Stop with a clear error instead. This also applies to the two-equation LM model of develop. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

Four tests on the 6:1 prolate spheroid at 15 degrees (INC_RANS, coarse unstructured grid of 100958 points): SST and SA with MENTER_SLM, with and without cross-flow. Each runs 5 iterations from a converged SST or SA transitional solution; the cross-flow tests restart from the solution of their base model. The grid and the solutions are in su2code/TestCases, branch feature_Trans_SLM. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

Temporary, for the prolate spheroid data of su2code/TestCases#207. To be reverted to develop before merging. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

rois1995 · 2026-09-24T09:56:24Z

Sorry for the long pause on this PR. I have now brought it up to date and reviewed it against the literature:

Merged develop. The one-equation model now uses the new edge-flux numerics of develop (CScalarFlux_TransLM with one variable), which also makes the convective flux of γ conservative (it was missing the density).
Checked the implementation against Menter et al. (2015), Vallinayagam Pillai & Lardeau (2017) and, for SA, Lee & Baeder (2021). The SA coupling previously reused the constants of the two-equation SA-LM coupling (Medida & Baeder); it now follows Lee & Baeder, who coupled exactly this one-equation model to SA.
Combinations not found in the literature (e.g. CODER_SLM, MOD_EPPLER_SLM with the one-equation model) now stop with a config error.
Fixed some bugs, also in the existing LM model: restart consistency, the free-stream Re_θt limit, missing AD inputs, and a crash with the discrete adjoint (now a clear error).
Added regression tests on the 6:1 prolate spheroid (TestCases#207).

I updated the description with the details. Validation runs are next.

bigfooted · 2026-09-26T20:45:21Z

Great, thanks! I made some changes to multigrid, so cases that use it will converge differently, hopefully in a positive way.

…onset with the cross-flow terms - LM_OPTIONS= LM2015 (the name in develop) is accepted again, as CROSSFLOW. - PRODLIM was parsed but never used: P_k^lim is part of the one-equation model of Menter et al. (2015, Eq. 25) and is always applied with it. - The F_ONSET output of the one-equation model now includes the cross-flow corrections, as the value used by the source term. - config_template.cfg documents the LM_OPTIONS of both models. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

…RB_INDEX - BC_Inlet of the transition solver scaled gamma and Re_theta_t as k and omega when an inlet profile file was used, although they are dimensionless and are not read from the file, and read a second variable that the one-equation model does not have. The nVar free-stream values are now used as they are. - With SA the turbulence intensity of the two-equation correlations is limited to 0.027 %, as with k-omega, since the correlations divide by it. - The AD preaccumulation registered two turbulence variables also for SA. - TURB_INDEX: the friction velocity is sqrt(nu |Omega|) at the wall in 2D and 3D (it used mu in 3D and the skin friction coefficient in 2D). - Debug comments replaced. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

EvertBunschoten

Thanks for the work, I think there are a couple of areas in the code that should be addressed to avoid code duplication and improve readability.

EvertBunschoten · 2026-09-30T14:56:06Z

+   */
+  template <typename Normals_type>
+  inline void SetNormal(unsigned long iPoint, Normals_type const& normal) {
+    for (unsigned long iDim = 0; iDim < nDim; iDim++) Normals(iPoint, iDim) = normal[iDim];


Use auto for loop index variables.

Done in afe6d25.

EvertBunschoten · 2026-09-30T14:59:46Z

+      su2vector<su2vector<su2matrix<su2double>>> WallNormal_container;
+      WallNormal_container.resize(nZone) = su2vector<su2matrix<su2double>>();
+      for (int iZone = 0; iZone < nZone; iZone++) {
+        const CConfig* config = config_container[iZone];
+        const CGeometry* geometry = geometry_container[iZone][iInst][MESH_0];
+        WallNormal_container[iZone].resize(geometry->GetnMarker());
+        for (auto iMarker = 0; iMarker < geometry->GetnMarker(); iMarker++) {
+          if (config->GetViscous_Wall(iMarker)) {
+            WallNormal_container[iZone][iMarker].resize(geometry->GetnElem_Bound(iMarker), 3);
+
+            for (auto iElem = 0u; iElem < geometry->GetnElem_Bound(iMarker); iElem++) {
+              su2vector<su2double> NormalHere;
+              NormalHere.resize(3) = su2double(0.0);
+
+              for (unsigned short iNode = 0; iNode < geometry->bound[iMarker][iElem]->GetnNodes(); iNode++) {
+                // Extract global coordinate of the node
+                unsigned long iPointHere = geometry->bound[iMarker][iElem]->GetNode(iNode);
+                long iVertexHere = geometry->nodes->GetVertex(iPointHere, iMarker);
+                for (auto iDim = 0u; iDim < geometry->GetnDim(); iDim++)
+                  NormalHere[iDim] += geometry->vertex[iMarker][iVertexHere]->GetNormal(iDim);
+              }
+
+              for (auto iDim = 0u; iDim < 3; iDim++) NormalHere[iDim] /= geometry->bound[iMarker][iElem]->GetnNodes();
+
+              su2double NormalMag = 0.0;
+              for (auto iDim = 0u; iDim < 3; iDim++) NormalMag += NormalHere[iDim] * NormalHere[iDim];
+              NormalMag = sqrt(NormalMag);
+
+              for (auto iDim = 0u; iDim < 3; iDim++) NormalHere[iDim] /= NormalMag;
+
+              for (auto iDim = 0u; iDim < 3; iDim++)
+                WallNormal_container[iZone][iMarker](iElem, iDim) = NormalHere[iDim];
+            }
+          } else {
+            WallNormal_container[iZone][iMarker].resize(1, 3) = su2double(0.0);
+          }
+        }
+      }


Please refactor this so it is easier to read.

Split in afe6d25 into WallElementNormal and CollectWallNormals.

EvertBunschoten · 2026-09-30T15:02:10Z

+                unsigned long iPointHere = geometry->bound[iMarker][iElem]->GetNode(iNode);
+                long iVertexHere = geometry->nodes->GetVertex(iPointHere, iMarker);


Be consistent with the point and vertex data type or just use auto.

Done in afe6d25, for the point, vertex and element indices.

EvertBunschoten · 2026-09-30T15:02:57Z

+              for (auto iDim = 0u; iDim < 3; iDim++) NormalHere[iDim] /= geometry->bound[iMarker][iElem]->GetnNodes();
+
+              su2double NormalMag = 0.0;
+              for (auto iDim = 0u; iDim < 3; iDim++) NormalMag += NormalHere[iDim] * NormalHere[iDim];


use the GeometryToolbox when possible.

Done: GeometryToolbox::Norm3 (afe6d25).

EvertBunschoten · 2026-09-30T15:15:49Z

+      auto normal_i =
+          make_pair(nZone, [config_container, geometry_container, iInst, WallNormal_container](unsigned long iZone) {
+            const CConfig* config = config_container[iZone];
+            const CGeometry* geometry = geometry_container[iZone][iInst][MESH_0];
+            const auto nMarker = geometry->GetnMarker();
+            const auto WallNormal = WallNormal_container[iZone];
+
+            return make_pair(nMarker, [config, geometry, WallNormal](unsigned long iMarker) {
+              auto nElem_Bou = geometry->GetnElem_Bound(iMarker);
+              if (!config->GetViscous_Wall(iMarker)) nElem_Bou = 1;
+
+              return make_pair(nElem_Bou, [WallNormal, iMarker](unsigned long iElem) {
+                const auto dimensions = 3;
+
+                return make_pair(dimensions, [WallNormal, iMarker, iElem](unsigned short iDim) {
+                  return WallNormal[iMarker](iElem, iDim);
+                });
+              });
+            });
+          });


Also refactor this

Done in afe6d25: separate marker/matrix helpers, and the lambdas read the collected normals by reference.

EvertBunschoten · 2026-10-05T15:58:58Z

+  VectorType normal_x;
+  VectorType normal_y;
+  VectorType normal_z;


Done in afe6d25: one MatrixType with three columns.

EvertBunschoten · 2026-10-05T16:04:36Z

+  /*--- Langtry-Menter correlation, the turbulence intensity (in percent) is limited to 0.027 to avoid the singularity. ---*/
+  const su2double Intensity = max(config->GetTurbulenceIntensity_FreeStream()*100.0, 0.027);
+  if (Intensity <= 1.3) {
+    Re_ThetaT_FreeStream = 1173.51 - 589.428*Intensity + 0.2196/(Intensity*Intensity);


Are the coefficients here constants from literature? If so, replace them with variables with names that can be easily retrieved from the reference.

Named the Langtry–Menter coefficients and linked the farfield correlation in afe6d25.

EvertBunschoten · 2026-10-05T16:06:44Z

  if (Intensity <= 1.3) {
-    if(Intensity >=0.027) {
-      ReThetaT_Inf = (1173.51-589.428*Intensity+0.2196/(Intensity*Intensity));
-    }
-    else {
-      ReThetaT_Inf = (1173.51-589.428*Intensity+0.2196/(0.27*0.27));
-    }
+    ReThetaT_Inf = (1173.51-589.428*Intensity+0.2196/(pow(max(Intensity, 0.027), 2.0)));
  }
  else if(Intensity>1.3) {
    ReThetaT_Inf = 331.5*pow(Intensity-0.5658,-0.671);
  }


This sequence is the same as used in CEulerSolver. Consider replacing it with a static method to avoid code duplication.

Done: both solvers use TransLMCorrelations::FreestreamReThetaT (afe6d25), including the same low-intensity limit.

EvertBunschoten · 2026-10-05T16:10:22Z

+
+    auto* flowNodes = su2staticcast_p<CFlowVariable*>(solver_container[FLOW_SOL]->GetNodes());
+    auto* turbNodes = su2staticcast_p<CTurbVariable*>(solver_container[TURB_SOL]->GetNodes());
+
+    SU2_OMP_FOR_STAT(omp_chunk_size)
+    for (unsigned long iPoint = 0; iPoint < nPoint; iPoint ++) {
+
+      // Here the nodes already have the new solution, thus I have to compute everything from scratch
+
+      const su2double rho = flowNodes->GetDensity(iPoint);
+      const su2double mu  = flowNodes->GetLaminarViscosity(iPoint);
+      const su2double muT = turbNodes->GetmuT(iPoint);
+      const su2double dist = geometry->nodes->GetWall_Distance(iPoint);
+      su2double VorticityMag = GeometryToolbox::Norm(3, flowNodes->GetVorticity(iPoint));
+      su2double StrainMag =flowNodes->GetStrainMag(iPoint);
+      VorticityMag = max(VorticityMag, 1e-12);
+      StrainMag = max(StrainMag, 1e-12); // safety against division by zero
+      const su2double Intermittency = nodes->GetSolution(iPoint,0);
+      const su2double Re_v = rho*dist*dist*StrainMag/mu;
+      const su2double vel_u = flowNodes->GetVelocity(iPoint, 0);
+      const su2double vel_v = flowNodes->GetVelocity(iPoint, 1);
+      const su2double vel_w = (nDim ==3) ? flowNodes->GetVelocity(iPoint, 2) : 0.0;
+      const su2double VelocityMag = max(sqrt(pow(vel_u, 2) + pow(vel_v, 2) + pow(vel_w, 2)), EPS);
+      su2double omega = 0.0;
+      su2double k = 0.0;
+      if(TurbFamily == TURB_FAMILY::KW){
+        omega = turbNodes->GetSolution(iPoint,1);
+        k = turbNodes->GetSolution(iPoint,0);
+      }
+
+      su2double Re_t = 0.0;
+      su2double Corr_Rec = 0.0;
+
+      if (options.SLM) {
+        Re_t = nodes->GetRe_t(iPoint);
+        Corr_Rec = nodes->GetCorr_Rec(iPoint);
+      } else {
+        Re_t = nodes->GetSolution(iPoint,1);
+        su2double Tu = 1.0;
+        if(TurbFamily == TURB_FAMILY::KW)
+          Tu = max(100.0*sqrt( 2.0 * k / 3.0 ) / VelocityMag,0.027);
+        if(TurbFamily == TURB_FAMILY::SA)
+          Tu = config->GetTurbulenceIntensity_FreeStream()*100;
+
+        Corr_Rec = TransCorrelations.ReThetaC_Correlations(Tu, Re_t);
+      }
+
+      su2double R_t = 1.0;
+      if(TurbFamily == TURB_FAMILY::KW)
+        R_t = rho*k/ mu/ omega;
+      if(TurbFamily == TURB_FAMILY::SA)
+        R_t = muT/ mu;
+
+      const su2double f_reattach = exp(-pow(R_t/20,4));
+
+      su2double f_wake = 0.0;
+      if(TurbFamily == TURB_FAMILY::KW){
+        const su2double re_omega = rho*omega*dist*dist/mu;
+        f_wake = exp(-pow(re_omega/(1.0e+05),2));
+      }
+      if(TurbFamily == TURB_FAMILY::SA)
+        f_wake = 1.0;
+
+      const su2double theta_bl   = Re_t*mu / rho /VelocityMag;
+      const su2double delta_bl   = 7.5*theta_bl;
+      const su2double delta      = 50.0*VorticityMag*dist/VelocityMag*delta_bl + 1e-20;
+      const su2double var1 = (Intermittency-1.0/50.0)/(1.0-1.0/50.0);
+      const su2double var2 = 1.0 - pow(var1,2.0);
+      const su2double f_theta = min(max(f_wake*exp(-pow(dist/delta, 4)), var2), 1.0);
+      su2double Intermittency_Sep = 2.0*max(0.0, Re_v/(3.235*Corr_Rec)-1.0)*f_reattach;
+      Intermittency_Sep = min(Intermittency_Sep,2.0)*f_theta;
+      Intermittency_Sep = min(max(0.0, Intermittency_Sep), 2.0);
+      nodes->SetIntermittencySep(iPoint, Intermittency_Sep);
+      nodes->SetIntermittencyEff(iPoint, Intermittency_Sep);


Move to separate function

Done: SetSeparationIntermittency (afe6d25).

EvertBunschoten · 2026-10-05T16:14:41Z

+  /*!
+   * \brief Set Value of Transition Momentum Thickness Reynolds number from correlations.
+   */
+  inline virtual void SetCorr_Rec(unsigned long iPoint, su2double val_Corr_Rec) {};


It's a lot of setters and getters. Would it be possible to replace these with a single accessor pair that accesses a struct with all the relevant data? Something like the parsed options used for some of the config options.

Done in afe6d25: a mutable/const accessor pair for TransitionLMData, shared by the source and output.

rois1995 · 2026-10-06T14:56:50Z

Addressed the review comments in afe6d25.

The MPI build and all 41 unit cases pass (74372 assertions). All eleven two-rank MPI runs finish successfully. The ten regular cases preserve the printed histories and volume fields; binary restart differences are at most 2.3e-13. All four existing SLM iteration-5 references pass at 1e-5.

The shared freestream correlation also makes the LM initialization use the same 0.027% intensity floor as the flow solver. The zero-intensity SA-LM check confirms this change; no reference values were updated.

All five forward/reverse CoDiPack interface syntax checks pass. These are compilation checks; sensitivities were not tested.

rois1995 added 3 commits January 27, 2023 18:21

- Just copied the CTransLMSolver header and source to commit

a7136d3

- Start of implementation of Simplified LM transition model

e9226db

- Modify number of variables if SLM option is used

d03742a

pr-triage Bot added the PR: unreviewed label Jan 27, 2023

pcarruscag changed the title ~~Implementation of Simplified LM Transition model [WIP]~~ [WIP] Implementation of Simplified LM Transition model Jan 27, 2023

pr-triage Bot removed the PR: unreviewed label Jan 27, 2023

pcarruscag reviewed Jan 27, 2023

View reviewed changes

github-advanced-security AI found potential problems Jan 27, 2023

View reviewed changes

Comment thread SU2_CFD/src/solvers/CTransLMSolver.cpp Fixed

- Updates on the transition model.

e17872d

- Added correlations for Simplified LM. - There is a bug on the computation of grad(n*U)*n

github-advanced-security AI found potential problems Jan 31, 2023

View reviewed changes

- Added Cross flow effects for SLM

4b9c615

- Modified LM_OPTIONS to include cross-flow effects: from LM2015 to CROSSFLOW

github-advanced-security AI found potential problems Jan 31, 2023

View reviewed changes

WallyMaier reviewed Feb 1, 2023

View reviewed changes

- Fixed Intermittency_Sep computation for SLM.

04cd56d

- Added computation for normals of grid points

e29b4df

- Added variables only for debug

Fixed normal extraction from structure

0319d62

Added computation of wall normals to the CPoint structure

b1edf61

rois1995 and others added 17 commits September 23, 2026 17:00

Use the codi and medi submodules of develop

2c0255e

Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

Remove debug and whitespace changes unrelated to the simplified LM model

4d4b0af

Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

Fix the AD build of the cross-flow vorticity direction

aefa3c9

The two branches of the conditional had an AD expression and a double. Co-Authored-By: Claude Opus 5.5 (1M context) <[email protected]> Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6

rois1995 mentioned this pull request Sep 24, 2026

Prolate spheroid grid and solutions for the one-equation LM transition tests su2code/TestCases#207

Open

rois1995 changed the title ~~Implementation of Simplified LM Transition model~~ One-equation (simplified) LM transition model for SST and SA Sep 24, 2026

github-advanced-security AI found potential problems Sep 24, 2026

View reviewed changes

Comment thread Common/src/CConfig.cpp Fixed

rois1995 and others added 2 commits September 28, 2026 15:59

rois1995 mentioned this pull request Sep 30, 2026

Align the SST model with the NASA TMR definitions and fix the boundary energy (#1851) #2329

Closed

5 of 6 tasks

Merge branch 'develop' into feature_Trans_SLM

2a187cb

EvertBunschoten requested changes Oct 5, 2026

View reviewed changes

Refactor LM transition diagnostics and correlations for review

afe6d25

	SU2_MPI::Error("Two correlations selected for Simplified model into LM_OPTIONS. Please choose only one.", CURRENT_FUNCTION);
	SU2_MPI::Error("Two correlations selected for simplified LM_OPTIONS. Please choose only one.", CURRENT_FUNCTION);

		unsigned long iPointHere = geometry->bound[iMarker][iElem]->GetNode(iNode);
		long iVertexHere = geometry->nodes->GetVertex(iPointHere, iMarker);

Conversation

rois1995 commented Jan 27, 2023 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Proposed Changes

Related Work

PR Checklist

Uh oh!

rois1995 commented Jan 27, 2023

Uh oh!

pcarruscag left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

rois1995 commented Jan 31, 2023

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

pcarruscag commented Jan 31, 2023

Uh oh!

github-advanced-security AI left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

bigfooted commented Feb 10, 2023

Uh oh!

rois1995 commented Feb 14, 2023

Uh oh!

pcarruscag commented Feb 14, 2023

Uh oh!

rois1995 commented Feb 14, 2023

Uh oh!

pcarruscag commented Feb 14, 2023

Uh oh!

rois1995 commented Feb 24, 2023

Uh oh!

pcarruscag commented Feb 24, 2023

Uh oh!

rois1995 commented Sep 24, 2026

Uh oh!

Uh oh!

bigfooted commented Sep 26, 2026

Uh oh!

EvertBunschoten left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

rois1995 commented Jan 27, 2023 •

edited

Loading