View Issue Details
ID | Project | Category | View Status | Date Submitted | Last Update |
---|---|---|---|---|---|
0003238 | OpenFOAM | Bug | public | 2019-05-16 13:38 | 2019-05-28 00:31 |
Reporter | RicardoFE95 | Assigned To | henry | ||
Priority | high | Severity | major | Reproducibility | always |
Status | closed | Resolution | unable to reproduce | ||
Platform | GNU/Linux | OS | Other | OS Version | (please specify) |
Summary | 0003238: The Reynolds normal stresses are being computed as negative in a recirculating flow by LienCubicKE model | ||||
Description | When running a simulation for an axisymmetric recirculating flow using the LienCubicKE model the Reynolds normal stresses come out as negative with a strange ripple pattern, particularly near the axis. This ripple pattern can be seen in other variables like p and nut. This physically impossible results have only been found while using the cubic model. The ShihQuadraticKE model does not show this behaviour when being ran on the same case. The standard kEpsilon model has also been used and physically possible results have been obtained as well. It seems that only LienCubicKE has this problem Several grid resolutions have been tried. 2D wedge meshes, 3D cyclic meshes and a full 360ยบ mesh have been tried as well and the same behaviour is observed. In every case the simpleFoam solver is used. OpenFOAM is being ran on a Ubuntu 16.04 virtual machine. | ||||
Steps To Reproduce | 1. Setup a case 2. Generate mesh 3. Initialize with potentialFoam or a converged kEpsilon solution 4. Run the simulation with LienCubicKE | ||||
Additional Information | I have found some posts in the cfd-online forums that describe problems with the LienCubicKE model, in particular this thread started by cfdmarkus in 2010 https://www.cfd-online.com/Forums/openfoam-solving/72527-numerical-problems-non-linear-ras-models.html I noticed that in the implementation of the cubic model Cmu is declared as a model coefficient at a constant value of 0.09. However, it is also calculated as a variable in the same implementation. I am not sure if this creates conflict with each other. In the paper referenced for the implementation of the ShihQuadraticKE model, it is said that if Cmu is constant it can cause negative stress values due to realizability. Perhaps this is the cause of the bug as in the quadratic model Cmu is not declared as a model coefficient. I have included the case I have been testing in the uploaded files. I had to delete the polyMesh directory in order to make it as small as possible. The mesh can be generated just by running blockMesh. | ||||
Tags | No tags attached. | ||||
|
|
|
Can you provide a patch to fix the issue with the LienCubicKE model? Have you tested using the field form of Cmu in the rest of the model? If there is no clear way to correct this model I can remove it. |
|
I have tried to compile a modified version of the LienCubicKE model by commenting out the declaration of Cmu as a model coefficient but have been running into trouble compiling mostly due to inexperience. I'm not sure if I'm following correctly, I'm assuming the field form of Cmu is what should be calculated rather than the constant value. That is implemented in the model but only inside the function correctNonlinearStress if I'm not mistaken. Do you mean to remove the model entirely? I will upload my attempt at correcting the model in case it is useful to you. |
|
As far as I am aware the current implementation is correct according to the original documents: Description Lien cubic non-linear low-Reynolds k-epsilon turbulence models for incompressible flows. This turbulence model is described in: \verbatim Lien, F.S., Chen, W.L. & Leschziner, M.A. (1996). Low-Reynolds-number eddy-viscosity modeling based on non-linear stress-strain/vorticity relations. Engineering Turbulence Modelling and Experiments 3, 91-100. \endverbatim Implemented according to the specification in: <a href= "http://personalpages.manchester.ac.uk/staff/david.d.apsley/specturb.pdf" >Apsley: Turbulence Models 2002</a> In addition to the low-Reynolds number damping functions support for wall-functions is also included to allow for low- and high-Reynolds number operation. If you can show that it is not implemented correctly and provide a fix that would be great. If the model is implemented correctly but does not operate in a useful way it can be removed from OpenFOAM. |
|
The current location for David Apsley's report on turbulence models is: https://personalpages.manchester.ac.uk/staff/david.d.apsley/turbmod.pdf |
|
I think there is something I am missing, that document states in section 2.7 that Cmu should not be a constant value but rather a function of both strain and vorticity. So why is it defined as a model coefficient? |
|
Cmu is not constant: volScalarField Cmu((2.0/3.0)/(Cmu1_ + sBar + Cmu2_*wBar)); volScalarField fMu(this->fMu()); nut_ = Cmu*fMu*sqr(k_)/epsilon_; nut_.correctBoundaryConditions(); nonlinearStress_ = fMu*k_ *( // Quadratic terms sqr(k_/epsilon_)/(Cbeta_ + pow3(sBar)) *( Cbeta1_*dev(innerSqr(S)) + Cbeta2_*twoSymm(S&W) + Cbeta3_*dev(symm(W&W)) ) // Cubic terms - pow3(Cmu*k_/epsilon_) *( (Cgamma1_*magSqr(S) - Cgamma2_*magSqr(W))*S + Cgamma4_*twoSymm((innerSqr(S)&W)) ) ); There is also a constant value used in the wall functions. |
|
If you believe the model is incorrect I would be happy to delete it. |
|
I don't think that the model is incorrect, just that there might be a problem somewhere in the implementation since I am not the only one to run into issues. The lines of code which are highlighted are from the function correctNonlinearStress which, if I am interpreting the code correctly, is only used at the last step. My worry behind this was that since Cmu is also being defined in lines 263 to 271 as: Cmu_ ( dimensioned<scalar>::lookupOrAddToDict ( "Cmu", coeffDict_, 0.09 ) ), That this might be causing issues, since in the other NLEVM implemented, which doesn't have this issue, this definition is absent. I do not mean for the model to be deleted as I believe it can be quite a useful alternative to more commonly used linear two equation models. |
|
As I said before Cmu_ is used for the wall damping functions which are defined with a constant Cmu. If you want to use the varying Cmu also for the wall damping functions feel free to do so but this is not consistent with their definition. |
|
As far as I am aware the LienCubicKE model implementation is correct. If you still believe there is an error in the implementation please provide details of the difference between the implementation and the definition of the model in the papers. |
Date Modified | Username | Field | Change |
---|---|---|---|
2019-05-16 13:38 | RicardoFE95 | New Issue | |
2019-05-16 13:38 | RicardoFE95 | File Added: LienCubicBugCase.zip | |
2019-05-16 14:41 | henry | Note Added: 0010466 | |
2019-05-16 15:06 | RicardoFE95 | File Added: myLienCubicKE.zip | |
2019-05-16 15:06 | RicardoFE95 | Note Added: 0010467 | |
2019-05-16 15:17 | henry | Note Added: 0010469 | |
2019-05-16 15:25 | henry | Note Added: 0010470 | |
2019-05-16 15:45 | RicardoFE95 | Note Added: 0010471 | |
2019-05-16 15:56 | henry | Note Added: 0010472 | |
2019-05-16 15:57 | henry | Note Added: 0010473 | |
2019-05-16 16:43 | RicardoFE95 | Note Added: 0010474 | |
2019-05-16 16:51 | henry | Note Added: 0010475 | |
2019-05-28 00:31 | henry | Assigned To | => henry |
2019-05-28 00:31 | henry | Status | new => closed |
2019-05-28 00:31 | henry | Resolution | open => unable to reproduce |
2019-05-28 00:31 | henry | Note Added: 0010496 |