ENH: New schemes for non-orthogonal meshes

Summary

Aim:

Improvement of numerical stability of non-orthogonal mesh simulations only via changes in the diffusion operator handling.​

(Theory) What is face non-orthogonality:

• "The angle between the centroid vector and the unit normal vector is the orthogonality angle." (Wimhurst, 2019)

(Theory) Why is face non-orthogonality important:

• Can affect the discrete diffusion operator:​
  • Numerical stability (hence cost)​
  • Level of errors​
• Numerical stability​
  • Larger explicit treatment -> less stability​
• Level of errors​
  • Larger corrections -> larger truncation errors​
  • Cannot be reduced by mesh refinement

(Theory) What is the technical challenge in face non-orthogonality:

• "The difficulty is to evaluate the dot product of the normal vector and the velocity gradient on the cell face." (Wimhurst, 2019)

(Theory) Treatments in OpenFOAM:

• Explicit correction by using an over-relaxed-type unit-normal vector decomposition

Resolved bugs

N/A

Methodology

Problem definition:

• Convergence rate​
  • Explicit component can become limiting factor for numerical stability​
  • During initial iterations, especially for steady-state runs or when starting from uniform fields in transient runs​
  • Yet its removal reduces accuracy
• Convergence order​
  • Non-orthogonality treatments are neglected at boundaries (except cyclic conditions)​
  • Even minor non-orthogonality on patches reduce the convergence order for cases where Dirichlet/Neumann boundary conditions are applied. (Noriega et al., 2018)

A solution: Field relaxation

• Field relaxation for explicit non-orthogonality components​
  • Under-relaxation for numerical stability​
  • Over-relaxation for numerical cost
• Implementation: a Laplacian scheme​
  • Based on an idea from nextFoam​
  • Named relaxedNonOrthoLaplacian​
• Implementation: a surface-normal gradient scheme:​
  • Named relaxedSnGrad

Test case:

• The DNS study of a smooth-wall plane channel flow by (Moser et al., 1999) where the friction Reynolds number of the flow is ReTau=395.​

Metrics:

• Numerical stability
  • Simulation crashes or not
• Prediction fidelity
  • Flow field predictions vs DNS statistics

Control variable:

Results

Discussion

What happened?

• New schemes​
  • Did not dampen initial residuals unlike nextFoam observations (at least for this particular case)​
  • Improved prediction fidelity for above 50 degrees​
  • Prevented simulations crashing after a certain non-orthogonality angle (above 50)​
• Numerical cost​
  • For below-50-degree non-orthogonality cases, the runtime remained similar to that of the base​
  • For above-50-degree cases, the numerical cost seems to be reduced

What do the results mean in practise? (How will it change in how we do things?)​

• New schemes may be used​
  • to avoid any simulation crashes due to non-orthogonality treatments​
    • without affecting the fidelity of numerical predictions​
    • without affecting numerical cost estimations​
  • more tests are needed

Risks

• No changes in existing output.
• No changes in existing user input.

Constraints

• No finite-area numerics​
• Needs extensive tests to assume safe for single-phase, multiphase and overset finite-volume applications​
• No boundary treatments for Dirichlet and Neumann conditions

Tests

• Compilation (incl. submodules):
  • linux64ClangDPInt32Opt (clang11)
  • linux64GccDPInt32Opt
  • linux64GccSPDPInt64Debug
• Alltest: No change in output with respect to the develop HEAD + no error

References

• Wimhurst, A. (2019). Mesh Non-Orthogonality. https://tinyurl.com/49b55tzw​
• Wimhurst, A. (2019). Mesh Non-Orthogonality 2: The Over-Relaxed Approach​ https://tinyurl.com/234h7hbt​
• Noriega et al. (2018). A case-study in open-source CFD code verification, Part I: Convergence rate loss diagnosis. DOI: 10.1016/j.matcom.2017.12.002

src/finiteVolume/Make/files

@@ -456,6 +456,7 @@ $(snGradSchemes)/orthogonalSnGrad/orthogonalSnGrads.C
 $(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGrads.C
 $(snGradSchemes)/linearFitSnGrad/linearFitSnGrads.C
 $(snGradSchemes)/skewCorrectedSnGrad/skewCorrectedSnGrads.C
+$(snGradSchemes)/relaxedSnGrad/relaxedSnGrads.C
 
 convectionSchemes = finiteVolume/convectionSchemes
 $(convectionSchemes)/convectionScheme/convectionSchemes.C
@@ -466,6 +467,7 @@ $(convectionSchemes)/boundedConvectionScheme/boundedConvectionSchemes.C
 laplacianSchemes = finiteVolume/laplacianSchemes
 $(laplacianSchemes)/laplacianScheme/laplacianSchemes.C
 $(laplacianSchemes)/gaussLaplacianScheme/gaussLaplacianSchemes.C
+$(laplacianSchemes)/relaxedNonOrthoGaussLaplacianScheme/relaxedNonOrthoGaussLaplacianSchemes.C
 
 finiteVolume/fvc/fvcFlux.C
 finiteVolume/fvc/fvcMeshPhi.C