CFD Online Discussion Forums - OF 2.2.0: “fixedMean” bug

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Bugs (https://www.cfd-online.com/Forums/openfoam-bugs/)

- - OF 2.2.0: “fixedMean” bug (https://www.cfd-online.com/Forums/openfoam-bugs/122224-2-2-0-fixedmean-bug.html)

OF 2.2.0: “fixedMean” bug

Hello,

there seems to be a bug in the "fixedMean" condition in OF-2.2.0. Small values for the pressure lead to wrong results. So using "0" isn't possible any more compared to OF-1.6-ext.

Here is the description of my case, a sketch is attached:

Differences in total pressure between In- and Outlet (deltaPtot) for a static (non-rotating) Pipe without any GGI/AMI after 500 iterations using simpleFoam and different OF-versions.

deltaPtot (OF-1.6-ext): meanValue = 0 -> deltaPtot = -136.18; meanValue = 1e-30 -> deltaPtot = -136.18; meanValue = 1e-10 -> deltaPtot = -136.178

deltaPtot (OF- 2.2.0): meanValue = 0 -> deltaPtot = -386.64; meanValue = 1e-30 -> deltaPtot = -386.64; meanValue = 1e-10 -> deltaPtot = -118.662

Has anyone a shared experience with this BC combined with Of 2.2.0???

//Markus

Did you use the same settings for both cases?
Did you check the difference in source code between the simpleFoam solvers of the different distributions? And of the boundary conditions?

settings are the same!
I have not checked the code yet. I am not well versed in this;)

//Markus

Hi Markus,

I might be a bit to late into the discussion, but could you please share with us the correct answer. I suppose that you have an analytical solution, which you are trying to compare with the numerical result?

Secondly, could you try to undertake the same exercise with a square cross section of the pipe, such that you are not having any non-orthogonalities in the mesh, i.e. a simple blockMesh and perhaps with a bit of stretching away from the boundaries if you fancy.

Thirdly. Are you using any turbulence models in your setup? If the answer is no, could you please try to run the circular cross section with any turbulence model of your liking and report your findings.

Fourthly, are the models converged after 500 iterations, or is this number arbitrarily chosen based on prior experience?

This was a lot of questions, once we have an answer to these, we might help you forward, but I would emphasise the importance of point 2 and 3.

Kind regards

Niels

The mesh was created with ICEM, the result of checkMesh can be seen below:

Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
points: 100542
faces: 294476
internal faces: 287644
cells: 97020
boundary patches: 4
point zones: 0
face zones: 3
cell zones: 1

Overall number of cells of each type:
hexahedra: 97020
prisms: 0
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 0
polyhedra: 0

Checking topology...
Boundary definition OK.
Point usage OK.
Upper triangular ordering OK.
Face vertices OK.
Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces ...
Patch Faces Points Surface topology
INLET 1260 1289 ok (non-closed singly connected)
OUTLET 1260 1289 ok (non-closed singly connected)
WALL_STAT 2128 2240 ok (non-closed singly connected)
WALL_ROT 2184 2240 ok (non-closed singly connected)

Checking geometry...
This is a 3-D mesh
Overall domain bounding box (0 -1 -1) (15 1 1)
Mesh (non-empty, non-wedge) directions (1 1 1)
Mesh (non-empty) directions (1 1 1)
Mesh (non-empty, non-wedge) dimensions 3
Boundary openness (2.60553e-17 2.62945e-16 2.45788e-17) Threshold = 1e-06 OK.
Max cell openness = 2.57603e-15 OK.
Max aspect ratio = 237.872 OK.
Minumum face area = 9.99999e-05. Maximum face area = 0.0346467. Face area magnitudes OK.
Min volume = 1.14835e-05. Max volume = 0.00160775. Total volume = 46.9561. Cell volumes OK.
Mesh non-orthogonality Max: 26.6077 average: 4.2815 Threshold = 70
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 0.618008 OK.

Mesh OK.

End

So there are no problems with non-orthogonal cells, the quality of the whole mesh is good. The cells are also expanded with approaching the In- and Outlet. The turbulence model is the kOmegaSST.

I also did some calculations with >2000 iterations but there is no difference in totalpressure from about 200 iterations -> 500 iterations are enough.

The analytical calculation leads to a result of deltaPtot = 104Pa, so the CFD-results of deltaPtot ~120-130Pa are quite proper but a deltaPtot of ~386Pa isn't adequate!! -> Problems with the fixedMean condition (OF 2.2.0) in combination with a very small pressure value (tested: 1e-10, 1e-30 and 0)!!

Unfortunately, I can't upload the case because of its size of 1.8MB!

//Markus

Hi Markus,

The reason I am asking about the non-orthogonalities is that the treatment of the non-orthogonal corrections have changed between the two versions [1], which could be the reason, why your see the large difference in pressure drop between 1.6-ext and 2.2.0. As you will see in the checkMesh results, your mesh is non-orthogonal (not severely), but non-orthogonal corrections will be finite.

Still, I agree with you that such a large impact from a small change in a boundary condition is worrying. Have you tried looking at the pressure at the inlet and outlets to see, whether the pressures are reasonably uniformly distributed?

Kind regards

Niels

[1]: I have been told that the change in delta coefficients used in the discretisation of the Laplacian operators has had an impact on the non-orthogonal correction as well.

Hy Niels,

there is no non-orthogonal correction in use for these calculations!

The differences in static pressure between the two OF-versions can be located both at the In- and Outlet!

//Markus

I had a quick read of the discussion thus far....have you tried to use the fixedMeanValue BC in 2.2.x directly? There was one provided

http://www.cfd-online.com/Forums/ope...tml#post418371 Post 40 for a 2.2.x compatible version.

that you can compile. I have not tested that version extensively in 2.2.x, but it will compile. Then you may be able to compare and eliminate some of the reasons or causes for the difference.

Good evening,

I went into the code, as it was a bit puzzling, and my previous thought are obviously redundant, so here it goes:

Code:

File: OpenFOAM-1.6-ext/src/finiteVolume/fields/fvPatchFields/derived/fixedMeanValue/fixedMeanValueFvPatchField.C



// Update the coefficients associated with the patch field

template<class Type>

void fixedMeanValueFvPatchField<Type>::updateCoeffs()

{

    if (this->updated())

    {

        return;

    }



    if (curTimeIndex_ != this->db().time().timeIndex())

    {

        Field<Type>& patchField = *this;



        Field<Type> pif = this->patchInternalField();



        patchField = meanValue_ - gAverage(pif) + pif;



        curTimeIndex_ = this->db().time().timeIndex();

    }



    fixedValueFvPatchField<Type>::updateCoeffs();

}

where the version 2.2.0 reads:

Code:

File: OpenFOAM-2.2.0/src/finiteVolume/fields/fvPatchFields/derived/fixedMean/fixedMeanFvPatchField.C



template<class Type>

void fixedMeanFvPatchField<Type>::updateCoeffs()

{

    if (this->updated())

    {

        return;

    }



    Field<Type> newValues(this->patchInternalField());



    Type meanValuePsi =

        gSum(this->patch().magSf()*newValues)

       /gSum(this->patch().magSf());



    if (mag(meanValue_) > SMALL && mag(meanValuePsi)/mag(meanValue_) > 0.5)

    {

        newValues *= mag(meanValue_)/mag(meanValuePsi);

    }

    else

    {

        newValues += (meanValue_ - meanValuePsi);

    }



    this->operator==(newValues);



    fixedValueFvPatchField<Type>::updateCoeffs();

}

The interesting part is that the two versions are supposed to work identical for small meanValue, as the 2.2.0 version then enters the else-statement. But, here the things become a bit tricky, because the function call

Code:

gAverage(pif)

does exactly this (i.e. simple average over all processors), whereas the OF-2.2.0 version computes an area weighted average. Therefore, the two methods does not give the same pressure at the boundary for small meanValue. For 'large' mean values, you can potentially enter the if-statement in OF2.2.0, but it depends on the ratio of the the area weight average to the meanValue.

I suppose that you use fixedMean in order to robustly handle such things as vortices approaching the outlet? This probably will not happen in your case, so could you try to use fixedValue on the pressure instead to check, whether it is that particular boundary condition or something else, where the answer is clouded by the behaviour of the boundary condition.

Kind regards

Niels

Ok, with fixedValue p = 0 I got deltaPtot = 119Pa - the same value as with fixedMean p = 1e-10 which leads also to the correct value of deltaPtot!

So how can I adptat the fixedMean condition in order to get a proper result with p = 0?? Because I'll need the fixedMean condition for my futher investigations.

//Markus

Quote:

I had a quick read of the discussion thus far....have you tried to use the fixedMeanValue BC in 2.2.x directly? There was one provided

Pressure outlet boundary condition Post 40 for a 2.2.x compatible version.

that you can compile. I have not tested that version extensively in 2.2.x, but it will compile. Then you may be able to compare and eliminate some of the reasons or causes for the difference.

__________________
Dan

on twitter @dancombest

I'll check this and give you reply about the results!!

//Markus

Ok, with this adapted "fixedMeanValue" the result is deltaPtot = 120Pa...so this is similar to the 119Pa I have reached with p=1e-10 in combination with fixedMean!!

Looks quite good!

Has anyone shared experiences with this adpated "fixedMeanValue" for OF 2.2.x???

//Markus

OK, so as I read it, two values of the flow resistance are suggested as the true answer. Namely original 1.6-ext computations and the recent fixedValue results. As neither compare well with the analytical result, I would suggest that you try the directMapped value on the inlet for U. This removes the development of the boundary layer through the pipe, so comparison with uniform pipeflow becomes easier.

After that you will have to decide on a formulation, which is consistent across an average value of 0. The implementation part will probably be the easy part.

Kind regards,

Niels