[b.simpson]

Hi,

I am trying to run a simulation using the BuoyantBoussinesqSimpleFoam solver in OpenFOAM v-1812. The simulation is of a naturally ventilated building with space on the outside for the external environment.

Currently when I try running BuoyantBoussinesqSimpleFoam the simulation is exiting due to an error. The error is as follows:

-------------------------------------------------------
Primary job  terminated normally, but 1 process returned
a non-zero exit code.. Per user-direction, the job has been aborted.
-------------------------------------------------------
--------------------------------------------------------------------------
mpirun detected that one or more processes exited with non-zero status, thus causing
the job to be terminated. The first process to do so was:

  Process name: [[5190,1],1]
  Exit code:    142
--------------------------------------------------------------------------

Does anyone have an idea of what might be causing the error and/or how I should go about troubleshooting this?

Any help on this issue would be much appreciated.

Kind regards,
Ben


My boundary conditions are as follows:

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p_rgh;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 1e5;

boundaryField
{
    "xmin|xmax|ymin|ymax|zmax"
    {
        type            totalPressure;
        p0              $internalField;
    }
    "zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ|Inlet"
    {
        type            fixedFluxPressure;
        value           $internalField;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      T;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 0 0 1 0 0 0];

internalField   uniform 293;

boundaryField
{
    "xmin|xmax|ymin|ymax|zmin|zmax|Ceiling|Floor|EWall|NWall|SWall|WWall|Inlet"
    {
        type            fixedValue;
        value           uniform 293;
    }
    PC
    {
        type            fixedGradient;
		gradient		uniform 100;
    }
    Occ
    {
		type            fixedGradient;
		gradient		uniform 130;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volVectorField;
    object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{
    "xmin|xmax|ymin|ymax|zmax"
    {
        type            pressureInletOutletVelocity;
        value           uniform (0 0 0);
    }
	"zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ"
    {
        type            noSlip;
    }
	Inlet
    {
        type            fixedValue;
        value           uniform (0 0.00047 0);
		//value           uniform (0 1 0);
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0";
    object      alphat;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -1 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    "xmin|xmax|ymin|ymax|zmax|Inlet"
    {
        type            calculated;
        value           uniform 0;
    }
	"zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ"
    {
        type            alphatJayatillekeWallFunction;
        Prt             0.85;
        value           uniform 0;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0";
    object      CO2;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 0 0 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    "xmin|xmax|ymin|ymax|zmax"
    {
        type            inletOutlet;
        inletValue      $internalField;
        value           $internalField;
    }
    "zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ"
    {
        type            zeroGradient;
    }
	Inlet
    {
        type            fixedValue;
        value           uniform 1.0;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0";
    object      epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

//epsilonInlet  0.03; // Cmu^0.75 * k^1.5 / L ; L =10

dimensions      [0 2 -3 0 0 0 0];

internalField   uniform 0.01;

boundaryField
{
   "xmin|xmax|ymin|ymax|zmax"
    {
        type            inletOutlet;
        inletValue      uniform 0.01;
        value           uniform 0.01;
    }
    "zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ"
    {
        type            epsilonWallFunction;
        value           $internalField;
    }
	Inlet
    {
        type            fixedValue;
        value           $internalField;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0";
    object      k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

kInlet          1.5;   // approx k = 1.5*(I*U)^2 ; I = 0.1

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0.1;

boundaryField
{
   "xmin|xmax|ymin|ymax|zmax"
    {
        type            inletOutlet;
        inletValue      uniform 0.01;
        value           uniform 0.01;
    }
    "zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ"
    {
        type            kqRWallFunction;
        value           uniform 0.1;
    }
	Inlet
    {
        type            fixedValue;
        value           uniform 0.1;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0";
    object      nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -1 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    "xmin|xmax|ymin|ymax|zmax|Inlet"
    {
        type            calculated;
        value           uniform 0;
    }
    "zmin|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ"
    {
        type            nutkWallFunction;
        value           uniform 0;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 1e5;

boundaryField
{
    "xmin|xmax|ymin|ymax|zmin|zmax|Ceiling|Floor|EWall|NWall|SWall|WWall|PC|Occ|Inlet"
    {
        type            fixedFluxPressure;
        rho             rhok;
        value           uniform 0;
    }
	
	#includeEtc "caseDicts/setConstraintTypes"
}

[Lisandro Maders]

Hello,

Can you run other cases in parallel (tutorials?).. Doesn't look like a case problem, but some wrong mpi setup or similar..

Also, did you remember to source OpenFOAM?

Cheers,

[b.simpson]

Thank you Lisandro, for your response.

I have tried running one of my other cases and the windaroundbuildings tutorial case, both run in parallel. Both ran to an endTime of 5000 without any errors.

I tried running my current case not in parallel and the solver stopped at a similar time to the parallel run, but without any error message or convergence message. I have copied the end of the log file below.

Can you explain what you mean by source OpenFOAM?

Kind regards,
Ben

Time = 499

DILUPBiCGStab:  Solving for Ux, Initial residual = 4.2035517e-12, Final residual = 4.2035517e-12, No Iterations 0
DILUPBiCGStab:  Solving for Uy, Initial residual = 1.3809356e-12, Final residual = 1.3809356e-12, No Iterations 0
DILUPBiCGStab:  Solving for Uz, Initial residual = 3.4626842e-11, Final residual = 3.4626842e-11, No Iterations 0
DILUPBiCGStab:  Solving for T, Initial residual = 5.7255095e-07, Final residual = 5.7255095e-07, No Iterations 0
DICPCG:  Solving for p_rgh, Initial residual = 4.8174887e-29, Final residual = 4.8174887e-29, No Iterations 0
time step continuity errors : sum local = 4.7379025e+47, global = 6.7267478e+38, cumulative = 1.3690775e+39
DILUPBiCGStab:  Solving for epsilon, Initial residual = 0.081559159, Final residual = 2.7467555e-07, No Iterations 1
bounding epsilon, min: -5.0589917e+76 max: 6.6333303e+104 average: 6.9901096e+98
DILUPBiCGStab:  Solving for k, Initial residual = 9.7152148e-06, Final residual = 9.7152148e-06, No Iterations 0
ExecutionTime = 12550.98 s  ClockTime = 12575 s

fieldMinMax fieldMinMax write:
    min(mag(U)) = 0 in cell 0 at location (-2.9166667 -2.9166667 0)
    max(mag(U)) = 8.3277972e+51 in cell 58642 at location (-0.28756519 5.8448262 0.38346076)
    min(CO2) = -1.3654693e-08 in cell 2506767 at location (13.100119 16.077201 2.082341)
    max(CO2) = 1 in cell 471004 at location (7.8713678 2.5999999 1.2113401)
    min(p) = -3.6330321e+87 in cell 1506619 at location (-0.14127546 8.297961 3.4)
    max(p) = 1.0130987e+73 in cell 595727 at location (-0.018164488 7.8713315 3.1237291)
    min(T) = 292.99634 in cell 524500 at location (11.752019 15.324934 0.72123229)
    max(T) = 476.00394 in cell 267900 at location (13.908153 16.782254 3.3148686)

scalarTransport execute: CO2
DILUPBiCGStab:  Solving for CO2, Initial residual = 2.1732333e-28, Final residual = 2.1732333e-28, No Iterations 0

Time = 500

DILUPBiCGStab:  Solving for Ux, Initial residual = 1, Final residual = 0.0032024424, No Iterations 1
DILUPBiCGStab:  Solving for Uy, Initial residual = 1, Final residual = 0.0032024424, No Iterations 1
DILUPBiCGStab:  Solving for Uz, Initial residual = 1, Final residual = 0.0032024579, No Iterations 1
DILUPBiCGStab:  Solving for T, Initial residual = 1.1094197e-06, Final residual = 1.1094197e-06, No Iterations 0

[HPE]

Hi

- The maximum magnitude of the velocity field is nonphysically large. Also, that of epsilon field.
- Could you double-check the mesh metrics? `checkMesh -allGeometry -allTopology`?
- Could you run your case with another turbulence closure model, e.g. omega-based model?

[b.simpson]

Thank you for your response.

I have run check mesh and couldnt see any issues in the log file.

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

Create mesh for time = 0

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

Enabling all geometry checks.

Time = 0

Mesh stats
    points:           349461
    faces:            1004544
    internal faces:   961536
    cells:            327680
    faces per cell:   6
    boundary patches: 7
    point zones:      0
    face zones:       0
    cell zones:       0

Overall number of cells of each type:
    hexahedra:     327680
    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.
    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
    xmin                2560     2709     ok (non-closed singly connected)   (-3 -3 0) (18.7 -3 3.4)
    xmax                2560     2709     ok (non-closed singly connected)   (-3 18.7 0) (18.7 18.7 3.4)
    ymin                2560     2709     ok (non-closed singly connected)   (-3 -3 0) (-3 18.7 3.4)
    ymax                2560     2709     ok (non-closed singly connected)   (18.7 -3 0) (18.7 18.7 3.4)
    zmin                16384    16641    ok (non-closed singly connected)   (-3 -3 0) (18.7 18.7 0)
    zmax                7920     8360     ok (non-closed singly connected)   (-3 -3 3.4) (18.7 18.7 3.4)
    zmax_ceiling        8464     8649     ok (non-closed singly connected)   (0 0 3.4) (15.7 15.7 3.4)

Checking faceZone topology for multiply connected surfaces...
    No faceZones found.

Checking basic cellZone addressing...
    No cellZones found.

Checking geometry...
    Overall domain bounding box (-3 -3 0) (18.7 18.7 3.4)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (4.4299191e-16 4.4299191e-16 -7.1450608e-15) OK.
    Max cell openness = 1.2245107e-16 OK.
    Max aspect ratio = 1.0239131 OK.
    Minimum face area = 0.027777776. Maximum face area = 0.029122167.  Face area magnitudes OK.
    Min volume = 0.0047222218. Max volume = 0.0049507683.  Total volume = 1601.026.  Cell volumes OK.
    Mesh non-orthogonality Max: 0 average: 0
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 1.2789791e-13 OK.
    Coupled point location match (average 0) OK.
    Face tets OK.
    Min/max edge length = 0.16666666 0.17065218 OK.
    All angles in faces OK.
    Face flatness (1 = flat, 0 = butterfly) : min = 1  average = 1
    All face flatness OK.
    Cell determinant (wellposedness) : minimum: 0.1249674 average: 0.93518669
    Cell determinant check OK.
    Concave cell check OK.
    Face interpolation weight : minimum: 0.49409237 average: 0.49993708
    Face interpolation weight check OK.
    Face volume ratio : minimum: 0.97664544 average: 0.99975127
    Face volume ratio check OK.

Mesh OK.

End

I have also just tried running the k-omegaSST turbulence model instead of k-epsilon. This ran without error past where the k-epsilon model usually failed. Unfortunately I was unable to run the model until convergence as my computer decided to crash several hours in

Here is the end of the log file before my computer died.

Time = 1613

DILUPBiCGStab:  Solving for Ux, Initial residual = 0.0022768434, Final residual = 3.8302314e-05, No Iterations 1
DILUPBiCGStab:  Solving for Uy, Initial residual = 0.00087106937, Final residual = 1.5127342e-05, No Iterations 1
DILUPBiCGStab:  Solving for Uz, Initial residual = 0.0016276516, Final residual = 2.8328342e-05, No Iterations 1
DILUPBiCGStab:  Solving for T, Initial residual = 0.0012589984, Final residual = 6.7196765e-05, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.0089650297, Final residual = 8.7014216e-05, No Iterations 55
time step continuity errors : sum local = 1.5389332e-05, global = -1.902396e-06, cumulative = -0.0041707033
smoothSolver:  Solving for omega, Initial residual = 0.00039958681, Final residual = 2.8080538e-05, No Iterations 3
DILUPBiCGStab:  Solving for k, Initial residual = 0.00018173884, Final residual = 4.1791532e-06, No Iterations 1
ExecutionTime = 17547.01 s  ClockTime = 17682 s

fieldMinMax fieldMinMax write:
    min(mag(U)) = 0 in cell 0 at location (-2.9166667 -2.9166667 0) on processor 0
    max(mag(U)) = 1.049388 in cell 108330 at location (1.8558424 15.718582 3.3788139) on processor 2
    min(CO2) = 0 in cell 0 at location (-2.9166667 -3 0.085) on processor 0
    max(CO2) = 1 in cell 277850 at location (1.899056 7.5999999 1.2113291) on processor 0
    min(p) = 99965.943 in cell 276497 at location (13.758831 15.718746 3.4) on processor 3
    max(p) = 100000 in cell 135893 at location (12.521604 15.718749 0) on processor 3
    min(T) = 293 in cell 0 at location (-2.9166667 -3 0.085) on processor 0
    max(T) = 1186.3781 in cell 108330 at location (1.8558424 15.718582 3.3788139) on processor 2

scalarTransport execute: CO2
DILUPBiCGStab:  Solving for CO2, Initial residual = 0.00036557033, Final residual = 6.8689503e-06, No Iterations 1

Time = 1614

DILUPBiCGStab:  Solving for Ux, Initial residual = 0.0029013552, Final residual = 4.6281166e-05, No Iterations 1
DILUPBiCGStab:  Solving for Uy, Initial residual = 0.0010731421, Final residual = 1.7715774e-05, No Iterations 1
DILUPBiCGStab:  Solving for Uz, Initial residual = 0.0017262862, Final residual = 2.9371945e-05, No Iterations 1
DILUPBiCGStab:  Solving for T, Initial residual = 0.0012577848, Final residual = 6.7378228e-05, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.0090579121, Final residual = 8.7629514e-05, No Iterations 43
time step continuity errors : sum local = 1.5462643e-05, global = -1.1040759e-06, cumulative = -0.0041718074
smoothSolver:  Solving for omega, Initial residual = 0.00039540712, Final residual = 2.8313504e-05, No Iterations 3
DILUPBiCGStab:  Solving for k, Initial residual = 0.00017238916, Final residual = 4.1214014e-06, No Iterations 1
ExecutionTime = 17560.16 s  ClockTime = 17698 s

fieldMinMax fieldMinMax write:
    min(mag(U)) = 0 in cell 0 at location (-2.9166667 -2.9166667 0) on processor 0
    max(mag(U)) = 1.046601 in cell 108330 at location (1.8558424 15.718582 3.3788139) on processor 2
    min(CO2) = 0 in cell 0 at location (-2.9166667 -3 0.085) on processor 0
    max(CO2) = 1 in cell 277850 at location (1.899056 7.5999999 1.2113291) on processor 0
    min(p) = 99965.943 in cell 271108 at location (1.9411685 15.718749 3.4) on processor 2
    max(p) = 100000 in cell 73613 at location (7.8713315 -0.018749525 0) on processor 1
    min(T) = 293 in cell 0 at location (-2.9166667 -3 0.085) on processor 0
    max(T) = 1180.8241 in cell 108330 at location (1.8558424 15.718582 3.3788139) on processor 2

scalarTransport execute: CO2
DILUPBiCGStab:  Solving for CO2, Initial residual = 0.00036510744, Final residual = 6.8622869e-06, No Iterations 1

Is there a reason why the k-omegaSST model would run without error but the k-epsilon model wont?

Kind regards,
Ben

[HPE]

The culprit is the epsilon. There would be lots of reasons, I'm afraid. Please do search web - there should already be lots of insightful remarks about epsilon's numerical fragility.

Good luck.

PS:Might be wise to check with other epsilon-based models, e.g. kepsilonphitf or rngkepsilon. etc.

[Lisandro Maders]

You could try to under-relax the fields a little bit. The BC's looks fine. This could be a mesh problem, can you share some image of a slice of your mesh?

Also, can you share the fvSchemes and fvSolution files you are using? Did you use the same files for both k-epsilon and kOmega-SST?

[b.simpson]

Thanks for the replies.

I have been looking into better defining my BC's (especially for k and epsilon) previously they were just taken from tutorials.

As for the fvSolutions and fvSchemes, I have copied them below. For the Omega run I used the same solutions and schemes as this but jut added Omega. For that I copied the Omega relevant parts from the motorbike tutorial. I will try a run with under-relaxing the fields.

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default         steadyState;
}

gradSchemes
{
    default         Gauss linear;
}

divSchemes
{
    default         none;
	
    div(phi,U)      bounded Gauss upwind;
    div(phi,T)      bounded Gauss upwind;
    div(phi,k)      bounded Gauss upwind;
    div(phi,epsilon) bounded Gauss upwind;
	div(phi,CO2) 	bounded Gauss upwind;
    div((nuEff*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}


// ************************************************************************* //

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    p_rgh
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0.01;
    }

    "(U|T|h|k|epsilon|CO2)"
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-05;
        relTol          0.1;
    }

    G
    {
        $p_rgh;
        tolerance       1e-05;
        relTol          0.1;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 0;
    pRefCell        0;
    pRefValue       0;

    residualControl
    {
        p_rgh           1e-4;
        U               1e-4;
	CO2             1e-4;
        T               1e-4;
	h               1e-4;
        G               1e-4;

        // possibly check turbulence fields
        "(k|epsilon|omega)" 1e-4;
    }
}

relaxationFactors
{
    fields
	{
	rho             1.0;
        p_rgh           0.9;
    }
    equations
    {
        U               0.3;
	CO2             0.3;
	h               0.9;
        T               0.5;
        "(k|epsilon|R)" 0.7;
	P		0.8;
	G               0.7;
    }
}


// ************************************************************************* //

Thanks again for the support.

Kind regards,

Ben

[b.simpson]

I have been able to get my model running by adjusting the relaxation factors for epsilon as you suggested. A value of 1 allowed my model to run without errors.

I still need to go back a better define my BC's but at least I have a functioning model for the time being.

Thanks for all the support with this issue.

Regards,
Ben