[avila.vc]

Dears,

I am running a case with twoPhaseEulerFoam solver (OF 30x). It is diverging (the error message is at the end of this post). I am trying to simulate a membrane bioreactor (membranes as porous medium), but it is diverging even with porosity model off (like a simple box with water and air). Could you, please, give me some light?

My case is attached, as well as a screenshot of the mesh. The results from checkMesh are below.

My case was based on the bubbleColumn tutorial.

It is a 3D "box" (0.5x0.5x1.1m) filled with water until 1m, with a baffle in the center (4 connected plates with about 2mm thick). Inside the baffle there are 4 boxes representing membrane modules and right below these there are 3 cilinders representing the aerators.

What I consider to be the main aspects that are different from the bubbleColumn tutorial are the following:
- time step ajustable with maxCo = 0.2 and initial time step = 0.00005;
- laplacian and snGrad schemes both corrected;
- solver and smoother for pressure: GAMG and GaussSeidel;
- tolerances: 1e-6 for pressure and 1e-5 for the others;
- nNonOrthogonalCorrectors 3;
- Max aspect ratio = 7.73712 OK.
- Mesh non-orthogonality Max: 68.481 average: 3.12309 Non-orthogonality check OK.
- Max skewness = 2.88954 OK.

Please, could you take a look and give me some hint?
Thank you.

/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.4.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
Build  : 2.4.x-51457432eeff
Exec   : checkMesh
Date   : Dec 08 2015
Time   : 23:44:24
Host   : "vini-VirtualBox"
PID    : 2253
Case   : /media/vini/vini1/OpenFOAM/vini-2.4.x/run/liu2
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           1345829
    faces:            4068617
    internal faces:   3760225
    cells:            1304807
    faces per cell:   6
    boundary patches: 7
    point zones:      0
    face zones:       0
    cell zones:       2

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

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
   *Number of regions: 9
    The mesh has multiple regions which are not connected by any face.
  <<Writing region information to "0/cellToRegion"
  <<Writing region 0 with 1093724 cells to cellSet region0
  <<Writing region 1 with 30752 cells to cellSet region1
  <<Writing region 2 with 43870 cells to cellSet region2
  <<Writing region 3 with 43870 cells to cellSet region3
  <<Writing region 4 with 43870 cells to cellSet region4
  <<Writing region 5 with 43870 cells to cellSet region5
  <<Writing region 6 with 1617 cells to cellSet region6
  <<Writing region 7 with 1617 cells to cellSet region7
  <<Writing region 8 with 1617 cells to cellSet region8

Checking patch topology for multiply connected surfaces...
    Patch               Faces    Points   Surface topology                  
    IN                  544      580      ok (non-closed singly connected)  
    OUT                 15       24       ok (non-closed singly connected)  
    BAFFLE              124000   62000    multiply connected (shared edge)  
    CILINDROS           9240     4680     multiply connected (shared edge)  
    MEMBRANA            93872    46944    multiply connected (shared edge)  
    OUTAR               6240     6405     ok (non-closed singly connected)  
    WALLS               74481    74687    ok (non-closed singly connected)  
  <<Writing 113504 conflicting points to set nonManifoldPoints

Checking geometry...
    Overall domain bounding box (-2.77556e-17 -1.38778e-17 -0.01) (0.5 1.1 0.51)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (6.52007e-15 5.58754e-16 3.41726e-17) OK.
    Max cell openness = 3.23933e-16 OK.
    Max aspect ratio = 7.73712 OK.
    Minimum face area = 6.08885e-06. Maximum face area = 0.000170173.  Face area magnitudes OK.
    Min volume = 1.97146e-08. Max volume = 1.1764e-06.  Total volume = 0.275101.  Cell volumes OK.
    Mesh non-orthogonality Max: 68.481 average: 3.12309
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 2.88954 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

Courant Number mean: 8.27337e-13 max: 8.35536e-07
Max Ur Courant Number = 0.248167
deltaT = 1.16476e-107
Time = 0.074858

PIMPLE: iteration 1
MULES: Solving for alpha.air
MULES: Solving for alpha.air
alpha.air volume fraction = 0.0918828  Min(alpha.air) = -0.016039  Max(alpha.air) = 1
Constructing face momentum equations
smoothSolver:  Solving for e.air, Initial residual = 1, Final residual = 1, No Iterations 0
smoothSolver:  Solving for e.water, Initial residual = 1, Final residual = 1, No Iterations 0
min T.air 298
min T.water 297.979
GAMG:  Solving for p_rgh, Initial residual = 2.08589e-10, Final residual = 2.08589e-10, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 2.08589e-10, Final residual = 2.08589e-10, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 2.08589e-10, Final residual = 2.08589e-10, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 2.08589e-10, Final residual = 2.08589e-10, No Iterations 0
PIMPLE: iteration 2
MULES: Solving for alpha.air
MULES: Solving for alpha.air
alpha.air volume fraction = 0.0918828  Min(alpha.air) = -0.016039  Max(alpha.air) = 1
Constructing face momentum equations
smoothSolver:  Solving for e.air, Initial residual = 1, Final residual = 1, No Iterations 0
smoothSolver:  Solving for e.water, Initial residual = 1, Final residual = 1, No Iterations 0
min T.air 298
min T.water 297.979
GAMG:  Solving for p_rgh, Initial residual = 2.15761e-10, Final residual = 2.15761e-10, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 2.15761e-10, Final residual = 2.15761e-10, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 2.15761e-10, Final residual = 2.15761e-10, No Iterations 0
GAMG:  Solving for p_rgh, Initial residual = 2.15761e-10, Final residual = 2.15761e-10, No Iterations 0
PIMPLE: iteration 3
MULES: Solving for alpha.air
MULES: Solving for alpha.air
alpha.air volume fraction = 0.0918828  Min(alpha.air) = -0.016039  Max(alpha.air) = 1
Constructing face momentum equations
[1] #0  Foam::error::printStack(Foam::Ostream&) at ??:?
[1] #1  Foam::sigFpe::sigHandler(int) at ??:?
[1] #2  ? in "/lib64/libc.so.6"
[1] #3  Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> > Foam::operator*<double, Foam::fvPatchField, Foam::volMesh>(Foam::dimensioned<double> const&, Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> > const&) at ??:?
[1] #4  Foam::fv::EulerDdtScheme<double>::fvcDdt(Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&, Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> const&) at ??:?
[1] #5  ? at ??:?
[1] #6  ? at ??:?
[1] #7  __libc_start_main in "/lib64/libc.so.6"
[1] #8  ? at /usr/src/packages/BUILD/glibc-2.11.1/csu/../sysdeps/x86_64/elf/start.S:116
[r1i0n14:05650] *** Process received signal ***
[r1i0n14:05650] Signal: Floating point exception (8)
[r1i0n14:05650] Signal code:  (-6)
[r1i0n14:05650] Failing at address: 0xc2f00001612
[r1i0n14:05650] [ 0] /lib64/libc.so.6(+0x329e0) [0x2b70a31ab9e0]
[r1i0n14:05650] [ 1] /lib64/libc.so.6(gsignal+0x35) [0x2b70a31ab945]
[r1i0n14:05650] [ 2] /lib64/libc.so.6(+0x329e0) [0x2b70a31ab9e0]
[r1i0n14:05650] [ 3] /home/u/avilavc/OpenFOAM/OpenFOAM-3.0.x/platforms/linux64GccDPInt32Opt/bin/mytwoPhaseEulerFoam(_ZN4FoammlIdNS_12fvPatchFieldENS_7volMeshEEENS_3tmpINS_14GeometricFieldIT_T0_T1_EEEERKNS_11dimensionedIdEERKS9_+0x692) [0x497e92]
[r1i0n14:05650] [ 4] /home/u/avilavc/OpenFOAM/OpenFOAM-3.0.x/platforms/linux64GccDPInt32Opt/lib/libfiniteVolume.so(_ZN4Foam2fv14EulerDdtSchemeIdE6fvcDdtERKNS_14GeometricFieldIdNS_12fvPatchFieldENS_7volMeshEEES8_S8_+0x765) [0x2b70a06bebc5]
[r1i0n14:05650] [ 5] /home/u/avilavc/OpenFOAM/OpenFOAM-3.0.x/platforms/linux64GccDPInt32Opt/bin/mytwoPhaseEulerFoam() [0x480fda]
[r1i0n14:05650] [ 6] /home/u/avilavc/OpenFOAM/OpenFOAM-3.0.x/platforms/linux64GccDPInt32Opt/bin/mytwoPhaseEulerFoam() [0x42bb98]
[r1i0n14:05650] [ 7] /lib64/libc.so.6(__libc_start_main+0xe6) [0x2b70a3197bc6]
[r1i0n14:05650] [ 8] /home/u/avilavc/OpenFOAM/OpenFOAM-3.0.x/platforms/linux64GccDPInt32Opt/bin/mytwoPhaseEulerFoam() [0x43f961]
[r1i0n14:05650] *** End of error message ***
--------------------------------------------------------------------------
mpiexec noticed that process rank 1 with PID 5650 on node r1i0n14 exited on signal 8 (Floating point exception).

[BlnPhoenix]

Hey, couple of thougts:

I think you may have an issue with your mesh. It says in CheckMesh:
"The mesh has multiple regions which are not connected by any face."

So you should check online this message, i think flow can not be calculated between these regions when not connected. I think this is most important.

Than, are you solving for Energy also? The solver seems not to iterate properly, on the other hand when switched off the temperatures should be mainly the same. So you can check "diameterModel" in phaseProperties and set both phases to "constant" and check the ThermophysicalProperties and set EquationOfState to "rhoConst". Furthermore you can set "minIter" for the Energy equations in fvSolution to "0". This should exclude any effect of energ transfer in your simulation if wish to do that. I experienced it can be the cause of divergence sometimes in twoPhaseEulerFoam.

Lastly: You should check the results from first couple of iterations (when you save every 0.01 s time step) and look for maybe unsuitable boundary conditions.
But work on the first thing from above first.

Just my two cents. Good luck.

[avila.vc]

BlnPhoenix, thank you for the time and tips you gave. I tried what you suggested, but it didnt work :/ . The boundary conditions seem OK to me, most of them are the same as in the bubbleColumn tutorial and the initial (first iterations) flow behaviour (with paraview) seems correct.

Well, it is still diverging. Any other ideas?

[BlnPhoenix]

Hm, have you been able to get rid of the message "The mesh has multiple regions which are not connected by any face." when you run checkMesh command?

[avila.vc]

BlnPhoenix, no. But I runned checkMesh on another mesh that worked with this case/solver and the same message was shown ("The mesh has multiple regions which are not connected by any face."). Therefore I thought it shouldnt be it the problem.

The mesh that converged OK had just a cube inside that big box and the checkMesh output showed there were 2 regions not connected by any face. I think it is just the regions delimited by the geometries. In this case here (that diverged) there are 9 regions (4 membranes modules, 3 aerators cilinders, 1 baffle and the external region). However, I am not sure, it is just a guess. What do you think? have you experienced something like this?

Thanks

[BlnPhoenix]

Hm, ok. I'm not sure about this error either but let's assume it's not ther root problem:

From your running log it's evident that you get a very small time step, so you can try an figure out which cells are causing this tiny time step, check the ones where the Courant Number gets large. You can do this in paraview http://www.cfd-online.com/Forums/ope...nt-number.html

Try to localize regions with high courant number and refine the mesh there. Also: Courant Number of 0.2 seems a bit low. Is there a reason for this? I usually go with something like 0.75 in twophaseeulerfoam.

Maybe also try to run the case in laminar first.

[avila.vc]

Already runned in laminar.

Unfortunately, the Courant number that gets higher than the maximum I estipulated is the Ur Co, (related to the relative velocity between phases, I guess). And the Co utility only writes the "normal" Courant number (the one that stays below 0.2)

However, I believe Co and Ur Co usually have the same behaviour, so I visualized the cells with high Co (0.19 - 0.20). It turns out that those cells are located at my outlet (my less refined region), which makes sense now that you said to refine the mesh region with high Co.

I am going to try this and return to you.

Thank you very much!!!

[avila.vc]

I forgot: actually there is no reason for Max Co = 0.2 yet. I was running like this before because I thought my problem got something to do with high Co number. I'll try higher Co also. Thanks

[avila.vc]

Originally Posted by BlnPhoenix

Hm, ok. I'm not sure about this error either but let's assume it's not ther root problem:

From your running log it's evident that you get a very small time step, so you can try an figure out which cells are causing this tiny time step, check the ones where the Courant Number gets large. You can do this in paraview http://www.cfd-online.com/Forums/ope...nt-number.html

Try to localize regions with high courant number and refine the mesh there. Also: Courant Number of 0.2 seems a bit low. Is there a reason for this? I usually go with something like 0.75 in twophaseeulerfoam.

Maybe also try to run the case in laminar first.

Didn't work :/ Anyone knows if twoPhaseEulerFoam solver has any problems running with big meshes (over 1 million)? When I run with pure water or pure air it works... maybe it got something to do with drag model? any other ideas? Thanks

[BlnPhoenix]

Hi,

shame it did not work. As far as i know twophaseeulerFoam is not limited to a number of cells ~1.3 Mio. cells should be perfectly ok, as long as your PC is ok with it. So the problem must be sth else.

TwoPhaseEulerFoam is very mesh sensitive from what i know and at your checkmesh log i see now: " mesh non-orthogonality Max: 68.481. I read somewhere that you should use Non-ortho corrector loops for a max non ortho value higher than 50 and completly re-mesh for a value higher than 60, could be wrong though, maybe you can research on that. So may advise would be re-mesh and try to get better non-ortho values and try again. And i'm still non perfectly convinced about this non connected areas in your domain..^^

Good luck!

[avila.vc]

My mistake: I had set a wrong boundary condition for the U.air on the aerators (air inlet cilynders). The correct one is surfaceNormalFixedValue

It seems to be working so far.

[jason]

Hi Vinícius,

Did you ever get twoPhaseEulerFoam working with porous regions?

I have tried to specify porosity using fvOptions and also with a constant/porosityProperties file in the standard format such as that found in the interFoam/angledDuct tutorial. This format works fine in a simpleFoam case.

I have tried OF 2.4.0, OF 3.0.0, OF 3.0.1, OF 3.0.x and none seem to add in the porosity region. In fact I get this error every time.

--> FOAM Warning :
From function void option::checkApplied() const
in file fvOptions/fvOption.C at line 372
Source porosity1 defined for field U but never used

Do I need to link any extra libs or anything?

If you were successful, please let me know and I can set up a smaller case and upload....

Br

Jason