[Samuel Twain]

Dear all:
I'm doing a simulation about wind aroud a flat roof. The model can be seen below.
[image]https://www.cfd-china.com/assets/uploads/files/1611984747439-wholemodel.png[/img]

In the reference paper, the author applied 'second order upwind' for space discretization in Fluent. And I have used the same model parameters (eg. domain size, minimum element size, boundary condition, turbulence model). simpleFOAM and realizable k-epsilon model are used. However I find that if I adopt different divergence schemes (namely, gauss upwind, linearUpwind and limitedLinear), the results will be different from each other and only the result of 'gauss upwind' agrees with that of reference paper.

From what I gather, in openfoam 'gauss linearUpwind' represents 'second order upwind' which was adopted in the reference paper. But when I use 'gauss linearUpwind' for the divergence schemes in OpenFOAM, the simulated result (namely, friction velocity on roof) differs greatly from that of reference paper. In comparison, the performance improves if I adopted 'gauss upwind' which represented first order upwind. The results of different div schemes can be seen below.
[image]https://www.cfd-china.com/assets/uploads/files/1611984885997-0f0345fa-d59b-4030-9578-0dcb60ab0d6f-image.png[/img]

To describe my question detailedly, some details of my fvSchems and fvSolutions are also posted.

ddtSchemes
{
    default         steadyState;
}

gradSchemes
{
    default         Gauss linear;

    limited         cellLimited leastSquares 1;
    //grad(U)         $limited;
    grad(k)         $limited;
    grad(epsilon)     $limited;
}

divSchemes
{
    default         none;

    div(phi,U)      Gauss limitedLinear 1;//bounded Gauss linearUpwind limited, Gauss upwind;

    turbulence      bounded Gauss linearUpwind limited;
    div(phi,k)      $turbulence;
    div(phi,epsilon) $turbulence;

    div((nuEff*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}

wallDist
{
    method meshWave;
}

solvers
{
    p
    {
        //solver          PCG;
        //preconditioner  DIC;
        solver          GAMG;
        tolerance       1e-7;
        relTol          0.1;
        smoother        DICGaussSeidel;
    }

    "(U|k|epsilon|omega)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        //solver          smoothSolver;
        //smoother        symGaussSeidel;
        tolerance       1e-7;
        relTol          0.001;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 2;
    consistent      yes;

    residualControl
    {
        p               1e-7;
        U               1e-7;
        "(k|epsilon|omega|f|v2)" 1e-7;
    }
}

relaxationFactors
{
    fields
    {
        p               0.3;
    }
    equations
    {
        U               0.7;
        "(k|omega|epsilon).*" 0.7;
    }
}

Actually, in many papers, first order upwind scheme are not encouraged to use. But in my simulation, it seems that I have to choose 'gauss upwind' to get a 'right' result compared with reference papar of which the author has adopted second order upwind schemes. And higher-order div schemes lead to something wrong.
That makes me confused. Could anyone please give me some suggestions?
Thank you
samuel

[hjasak]

You are using the wrong kind of schemes on bounded scalars - this will destabilize the solver.

If numerical diffusion is causing trouble, you can simply refine the mesh.

Hrv

[Samuel Twain]

Originally Posted by hjasak

You are using the wrong kind of schemes on bounded scalars - this will destabilize the solver.

If numerical diffusion is causing trouble, you can simply refine the mesh.

Hrv

Thank you, Hrv. I will remove "bounded" , refine the mesh and find out whether it will cause peoblems.
Before that, I have found that my “k” boundary type seems to be unreasonable. Originally, I used kLowReWallFunction, after changing it to kqrWallFunction, the result becomes better.
This also makes me confused. According to description, kLowREWallFunction is suitable for both low and high RE condition, while kqrWallFunction is only for high RE condtion. Is kLowREWallFunction unsuitable for realizable k-epsilon turbulence model?

[hbulus]

Originally Posted by hjasak

You are using the wrong kind of schemes on bounded scalars - this will destabilize the solver.

If numerical diffusion is causing trouble, you can simply refine the mesh.

Hrv

You mean limitedLinear would be better for scalars, and linearUpwind for U, right ?

[piu58]

You may use Gauss upwind if you are looking for robust schemata. More accurate is linear upwind, but more oscillatory too.