[fusij]

Hello,

I have a free-surface flow that is similar to the waterchannel tutorial. I am using Local Time Stepping (LTS) but I have convergence problems. I will below go through my case step by step and provide the case files if someone is interested to take a look and give me some pointers.

A schematic of the geometry can be seen below



I have used cfMesh to create the mesh, here is the output of the checkMesh -allGeometry -allTopology:

Create time

Create polyMesh for time = 0

Enabling all (cell, face, edge, point) topology checks.

Enabling all geometry checks.

Time = 0

Mesh stats
    points:           385360
    internal points:  288677
    edges:            1074519
    internal edges:   882306
    internal edges using one boundary point:   93367
    internal edges using two boundary points:  0
    faces:            994222
    internal faces:   898692
    cells:            305063
    faces per cell:   6.20499
    boundary patches: 5
    point zones:      0
    face zones:       0
    cell zones:       0

Overall number of cells of each type:
    hexahedra:     268418
    prisms:        2112
    wedges:        0
    pyramids:      5269
    tet wedges:    0
    tetrahedra:    2336
    polyhedra:     26928
    Breakdown of polyhedra by number of faces:
        faces   number of cells
            6   6170
            7   3188
            8   599
            9   11939
           10   1
           12   3697
           13   1
           15   1273
           16   1
           18   54
           21   4
           24   1

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
    Topological cell zip-up check OK.
    Face-face connectivity OK.
    Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces...
                   Patch    Faces   Points                  Surface topology Bounding box
                    atmo    29308    30048  ok (non-closed singly connected) (34.1243 215.677 108.957) (71.5834 236.001 114)
                   inlet     3800     4181  ok (non-closed singly connected) (51.0844 217.546 107.302) (71.5835 226.81 110.267)
                  outlet     3283     3598  ok (non-closed singly connected) (21.1025 219.023 104.886) (25.9396 238.029 108.5)
             sula-fleyti     9304     9450  ok (non-closed singly connected) (31.7173 224.968 105.022) (35.0795 226.499 108.936)
                   walli    49835    51290  ok (non-closed singly connected) (21.1025 215.677 104.886) (71.5835 238.029 110.267)

Checking geometry...
    Overall domain bounding box (21.1025 215.677 104.886) (71.5835 238.029 114)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (1.89714e-16 1.40128e-15 1.1662e-15) OK.
    Max cell openness = 3.27551e-16 OK.
    Max aspect ratio = 10.2315 OK.
    Minimum face area = 0.000174025. Maximum face area = 3.46417.  Face area magnitudes OK.
    Min volume = 1.50399e-06. Max volume = 8.48143.  Total volume = 3955.37.  Cell volumes OK.
    Mesh non-orthogonality Max: 64.2831 average: 9.91544
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 2.10038 OK.
    Coupled point location match (average 0) OK.
 ***Error in face tets: 16 faces with low quality or negative volume decomposition tets.
  <<Writing 16 faces with low quality or negative volume decomposition tets to set lowQualityTetFaces
    Min/max edge length = 0.00357058 1.46541 OK.
   *There are 2062 faces with concave angles between consecutive edges. Max concave angle = 63.3745 degrees.
  <<Writing 2062 faces with concave angles to set concaveFaces
    Face flatness (1 = flat, 0 = butterfly) : min = 0.80978  average = 0.999301
    All face flatness OK.
    Cell determinant (wellposedness) : minimum: 0.0415301 average: 12.492
    Cell determinant check OK.
 ***Concave cells (using face planes) found, number of cells: 2444
  <<Writing 2444 concave cells to set concaveCells
    Face interpolation weight : minimum: 0.0647744 average: 0.45461
    Face interpolation weight check OK.
    Face volume ratio : minimum: 0.0177947 average: 0.801404
    Face volume ratio check OK.

Failed 2 mesh checks.

End

My interpretation of the failed mesh checks is that I should be able to solve a case on this mesh, low number of lowQualityTetFaces should not have a detrimental effect on the convergence. I have tried several configurations for cfMesh and it always seems to create a certain amount of lowQualityTetFaces.

Here is an image showing some of the lowQualityTetFaces in the domain



This is a list of how I set the boundary conditions, provided by patchSummary

Valid fields:
    volScalarField    nut
    volScalarField    alpha.water
    volVectorField    U
    volScalarField    k
    volScalarField    alpha.water.org
    volScalarField    p_rgh
    volScalarField    omega

patch    : atmo
    scalar        nut        calculated
    scalar        alpha.water        inletOutlet
    scalar        k        inletOutlet
    scalar        alpha.water.org        inletOutlet
    scalar        p_rgh        totalPressure
    scalar        omega        inletOutlet
    vector        U        pressureInletOutletVelocity

patch    : inlet
    scalar        nut        calculated
    scalar        alpha.water        variableHeightFlowRate
    scalar        k        fixedValue
    scalar        alpha.water.org        variableHeightFlowRate
    scalar        p_rgh        fixedFluxPressure
    scalar        omega        fixedValue
    vector        U        variableHeightFlowRateInletVelocity

patch    : outlet
    scalar        nut        calculated
    scalar        alpha.water        zeroGradient
    scalar        k        inletOutlet
    scalar        alpha.water.org        zeroGradient
    scalar        p_rgh        fixedFluxPressure
    scalar        omega        inletOutlet
    vector        U        zeroGradient

group    : wall
    scalar        nut        generic
    scalar        alpha.water        zeroGradient
    scalar        k        generic
    scalar        alpha.water.org        zeroGradient
    scalar        p_rgh        fixedFluxPressure
    scalar        omega        generic
    vector        U        fixedValue

End

I use setFields to have a initial water level in my domain as seen below

I have basically adopted the fvSchemes and fvSolution from the DTCHull tutorial that uses localEuler (LTS) as its ddt scheme. I squeezed down the tolerances on the solvers to lower the continuity error in each timestep and also played with the damping/smoothing coefficients. You can find the controlDict, fvSchemes and fvSolution files attached to this post.

I am confused now if the mesh, B.C. or schemes/solvers are causing the divergence. If you can give me comments on this case that would be appreciated.

You can find the case files here: case-files

[KateBradbrook]

What type of divergence appears first?

If bounding message for turbulence parameters - in my experience, this implies grid problems - look at grid in area of high turbulence parameters.

If you quickly get strange velocity vectors or alpha behavious at/near boundaries, then it is probably a BC issue. At outlet you have fixedFluxPressure for P and zeroGradient for U. Whilst this can be consistent for simple cases, there is also the potential for small inconsistencies (particularly -ve velocities) to feedback into strange behaviour. Especially when starting from a zero velocity condition with LTS. Try a fixed value for one of them (e.g. calculate U based on inflow/area at outlet for your initial water level and keep P as fixedFluxPressure; or use a free-outfall (P=0) or hydrostatic pressure condition and keep U as zeroGradient) - at least to get a stable solution and then you can adjust later....

Hope these thoughts help. If things blow up quickly it IS normally grid, BC, or initial conditions - changing solution schemes or parameters in my experience only delays divergence if there is an underlying issue with problem specification, so best to go back to less restrained tolerances and know more quickly if it is going to work or not!

[fusij]

Hey Kate, thank you for your answer.

After fixing the pressure at the outlet, everything started to run smoothly. Now I converge in 10000-12000 timesteps using LTS with Co=10.

So having few lowQualityTetFaces in your domain is not prohibiting convergence.

My issues came from the BC's. Letting the pressure float at the outlet and only fixing pressure for atmosphere BC was not enough.

[franjo_j]

Originally Posted by fusij

Hey Kate, thank you for your answer.

After fixing the pressure at the outlet, everything started to run smoothly. Now I converge in 10000-12000 timesteps using LTS with Co=10.

So having few lowQualityTetFaces in your domain is not prohibiting convergence.

My issues came from the BC's. Letting the pressure float at the outlet and only fixing pressure for atmosphere BC was not enough.

My experiences show that lowQualityTetFaces cause problems to lagrangian simulations. Furthermore, the mesh is useless if you have any negative pyramids. In other cases, convergence depends on your settings and the problem at hand.

[fusij]

Originally Posted by franjo_j

My experiences show that lowQualityTetFaces cause problems to lagrangian simulations. Furthermore, the mesh is useless if you have any negative pyramids. In other cases, convergence depends on your settings and the problem at hand.

Thanks for the input. I have noticed that cfMesh creates these lowQualTetFaces quite often, at least in my case. They are spread at regular intervals and at a certain distance from the surface. This distance changes when I thicken/decrease my refinement area but does not remove these faces. Do you know why this could be? Has it something to do with cfMesh or is it problem dependent?

I really like using cfMesh, especially the ease of use compared to other OS meshers. It gives you good meshes in a much shorter time, both because of easy setup and quick run-times.

[w.chaabouni]

Hello everyone,

I am running a simulation based on the tutorial "/tutorials/multiphase/interFoam/laminar/waveExampleStokesI". When I run simulation, I have a similar ERROR


/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v1812 |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : v1812 OPENFOAM=1812
Arch : "LSB;label=32;scalar=64"
Exec : interFoam -dry-run-write
Date : May 07 2019
Time : 14:03:22
Host : default
PID : 1007
I/O : uncollated
Case : /home/ofuser/workingDir/OpenFoam/tests/case5
nProcs : 1
trapFpe: Floating point exception trapping enabled (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 10)
allowSystemOperations : Allowing user-supplied system call operations

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

Create mesh for time = 0

Operating in 'dry-run' mode: case will run for 1 time step. All checks assumed OK on a clean exit
Selecting simplified mesh model staticFvMesh
Creating simplified mesh using "/home/ofuser/workingDir/OpenFoam/tests/case5/constant/polyMesh"
Mesh bounds: (0 0 0) (10 0.04 30)

PIMPLE: Operating solver in PISO mode

Reading field p_rgh

Reading field U

Reading/calculating face flux field phi

Reading transportProperties

Selecting incompressible transport model Newtonian
Selecting incompressible transport model Newtonian
Selecting turbulence model type laminar
Selecting laminar stress model Stokes

Reading g

Reading hRef
Calculating field g.h

No MRF models present

No finite volume options present
DICPCG: Solving for pcorr, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0, global = 0, cumulative = 0
Courant Number mean: 0 max: 0

Starting time loop

Courant Number mean: 0 max: 0
Interface Courant Number mean: 0 max: 0
Time = 0.01

PIMPLE: iteration 1
Selecting waveModel StokesI
#0 Foam::error:rintStack(Foam::Ostream&) at ??:?
#1 Foam::sigFpe::sigHandler(int) at ??:?
#2 ? in /lib64/libc.so.6
#3 Foam::waveModel::initialiseGeometry() at ??:?
#4 Foam::waveModel::readDict(Foam::dictionary const&) at ??:?
#5 Foam::waveModels::waveGenerationModel::readDict(Fo am::dictionary const&) at ??:?
#6 Foam::waveModels::irregularWaveModel::readDict(Foa m::dictionary const&) at ??:?
#7 Foam::waveModels::regularWaveModel::readDict(Foam: :dictionary const&) at ??:?
#8 Foam::waveModels::StokesI::readDict(Foam::dictiona ry const&) at ??:?
#9 Foam::waveModels::StokesI::StokesI(Foam::dictionar y const&, Foam::fvMesh const&, Foam:olyPatch const&, bool) at ??:?
#10 Foam::waveModel::addpatchConstructorToTable<Foam:: waveModels::StokesI>::New(Foam::dictionary const&, Foam::fvMesh const&, Foam:olyPatch const&) at ??:?
#11 Foam::waveModel::New(Foam::word const&, Foam::fvMesh const&, Foam:olyPatch const&) at ??:?
#12 Foam::waveModel::lookupOrCreate(Foam:olyPatch const&, Foam::fvMesh const&, Foam::word const&) at ??:?
#13 Foam::waveAlphaFvPatchScalarField::updateCoeffs() at ??:?
#14 Foam::fvPatchField<double>::evaluate(Foam::UPstrea m::commsTypes) at ??:?
#15 Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>::Boundary::evaluate() at ??:?
#16 ? at interFoam.C:?
#17 ? at ??:?
#18 __libc_start_main in /lib64/libc.so.6
#19 ? at ??:?
Floating point exception




Can anyone help me with it please?