need help with setting up buoyantSimpleFoam case

boffin5 · May 12, 2024, 13:18

Hi all,

For a long time I've been working on a case with a radiator in a duct, mounted on an airplane. An image is attached. Of late, I have been trying to use buoyantSimpleFoam to simulate it. As usual, I use tutorials for guidance, but I am finding that boundary conditions for the various tutorials are all over the map.

To construct the geometry, I used ennova, and the output was 2 meshes, one for the airplane/duct and the other for the heat exchanger. They were created simultaneously, therefore the meshes interface perfectly. In OpenFoam, I used mergeMesh to combine them, yielding 2 regions, and then stitchMesh to relate the faces at the radiator inlet and outlet faces, resulting in one region. Then I used topoSet to collect the radiator cells into a porous zone. Here is the output from topoSet:

Code:

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

Create polyMesh for time = 0

Reading topoSetDict

Time = 0
    mesh not changed.
Created cellSet porousCells
    Applying source boxToCell
    Adding cells with center within boxes 1((2.02542 -0.22352 -0.64202) (2.08257 0 -0.47057))
    cellSet porousCells now size 4197
Created cellZoneSet porousZone
    Applying source setToCellZone
    Adding all cells from cellSet porousCells ...
    cellZoneSet porousZone now size 4197
End

In paraView, I click on 'Read zones' and 'Copy data to cell zones' to visualize the porous zone. But another way to do it is by this command: 'foamToVTK -cellSet porousZone', whereupon the output is:

Code:

Create time

Create mesh for time = 0

Keeping old VTK files in "/home/boffin5/cfdaero/meredith-buoyant-stitch/VTK"

--> FOAM FATAL ERROR: 
Problem

    From function void Foam::fvMeshSubset::setLargeCellSubset(const labelList&, Foam::label, Foam::label, bool)
    in file fvMeshSubset/fvMeshSubset.C at line 1172.

FOAM aborting

#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::error::abort() at ??:?
#2  Foam::fvMeshSubset::setLargeCellSubset(Foam::List<int> const&, int, int, bool) at ??:?

So I'm not sure what is going on. The output from checkMesh shows 1 highly skewed face; however this is not necessarily fatal, and I have edited the fvSchemes file to include 'skewCorrected' to account for the problem face.

Code:

Exec   : checkMesh

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

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           811893
    faces:            975893
    internal faces:   925868
    cells:            148777
    faces per cell:   12.7826
    boundary patches: 11
    point zones:      1
    face zones:       3
    cell zones:       1

Overall number of cells of each type:
    hexahedra:     1924
    prisms:        8
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    polyhedra:     146845
    Breakdown of polyhedra by number of faces:
        faces   number of cells
            5   4
            6   449
            7   3188
            8   5529
            9   10506
           10   14141
           11   16019
           12   15209
           13   18605
           14   20894
           15   15960
           16   11265
           17   7344
           18   3913
           19   2047
           20   904
           21   399
           22   191
           23   93
           24   63
           25   39
           26   22
           27   11
           28   10
           29   8
           30   3
           31   2
           32   2
           33   1
           34   2
           35   2
           36   3
           37   3
           39   1
           40   1
           41   1
           42   1
           43   2
           52   2
           54   1
           56   1
           66   1
           67   1
           83   1
           86   1

Checking topology...
    Boundary definition OK.
    Cell to face addressing 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
                    flap     5531    10013  ok (non-closed singly connected)
                fuselage     6860    14868  ok (non-closed singly connected)
                    lips     4865     7908  ok (non-closed singly connected)
                 fairing    20535    32634  ok (non-closed singly connected)
                 spinner      359      945  ok (non-closed singly connected)
            trailingedge      635      826  ok (non-closed singly connected)
                   edges     1640     2251  ok (non-closed singly connected)
                   inlet      135      272  ok (non-closed singly connected)
                  outlet      135      272  ok (non-closed singly connected)
                frontier     2317     4506  ok (non-closed singly connected)
                 clplane     6207    12369  ok (non-closed singly connected)

Checking geometry...
    Overall domain bounding box (-15 -7 -7) (50 3.98651e-09 7)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (-1.55353e-18 5.92364e-16 2.52647e-17) OK.
    Max cell openness = 3.46311e-16 OK.
    Max aspect ratio = 34.163 OK.
    Minimum face area = 5.32353e-10. Maximum face area = 1.18133.  Face area magnitudes OK.
    Min volume = 9.63107e-13. Max volume = 1.45954.  Total volume = 6368.32.  Cell volumes OK.
    Mesh non-orthogonality Max: 84.2436 average: 10.4766
   *Number of severely non-orthogonal (> 70 degrees) faces: 145.
    Non-orthogonality check OK.
  <<Writing 145 non-orthogonal faces to set nonOrthoFaces
    Face pyramids OK.
 ***Max skewness = 4.63047, 1 highly skew faces detected which may impair the quality of the results
  <<Writing 1 skew faces to set skewFaces
    Coupled point location match (average 0) OK.

Failed 1 mesh checks.

End

Optimistically, I choose to believe that the porous zone cellZoneSet is okay, per the topoSet output. Now for the boundary conditions:

Code:

object      alphat;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

//internalField   uniform 0;
internalField   uniform 1e-3;

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
    
    inlet
    {
        type            calculated;
        value           $internalField;
    }
    
    outlet
    {
        type            calculated;
        value           $internalField;
    }
    
    fuselage
    {
        type            compressible::alphatWallFunction;
        value           $internalField;
    }
    
    spinner
    {
        type            compressible::alphatWallFunction;
        value           $internalField;
    }
    
    lips
    {
        type            compressible::alphatWallFunction;
        value           $internalField;
    }
    
    fairing
    {
        type            compressible::alphatWallFunction;
        value           $internalField;
    }
    
    flap
    {
        type            compressible::alphatWallFunction;
        value           $internalField;
    }
   
    edges {
        type            compressible::alphatWallFunction;
        value           $internalField;
        }
    trailingedge {
        type            compressible::alphatWallFunction;
        value           $internalField;
        }
}

After 'fuselage' are the various pieces of the geometry that make up the model. For brevity I'll only show 'fuselage' for the remaining BCs.

Code:

object      epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

//internalField   uniform 0.54;
internalField   uniform 200;

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
    
    inlet
    {
        type            turbulentMixingLengthDissipationRateInlet;
        mixingLength    0.005;
        value           uniform 200;    
     }
    
    outlet
    {
        type            inletOutlet;
        inletValue      $internalField;
        value           $internalField;
    }
    
    fuselage
    {
        type            epsilonWallFunction;
        value           $internalField;
    }

Code:

object      k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 1;    // 0.375

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
    
    inlet
    {
        type            turbulentIntensityKineticEnergyInlet;
        intensity       0.05;
        value           uniform 1;
    }
    
    outlet
    {
        type            inletOutlet;
        inletValue      uniform 1;
        value           uniform 1;
    }
    
    fuselage
    {
        type            kqRWallFunction;
        value           uniform 1;
    }

Code:

object      nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 0;

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
    
    inlet
    {
        type            calculated;
        value           uniform 0;
    }
    
    outlet
    {
        type            calculated;
        value           uniform 0;
    }
    
    fuselage
    {
        type            nutkWallFunction;
        value           uniform 0;
    }

Code:

 object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 1e5;  //90812

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
       
    inlet
    {
        type            zeroGradient;
    }
    
    outlet
    {
        type            fixedValue;
        value           $internalField;
    }
    
    fuselage
    {
        type             zeroGradient;
    }

Code:

object      p_rgh;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 1e5;  //84559

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
    
    inlet
    {
        type            fixedFluxPressure;
        value           uniform 0;
     }
    
    outlet
    {
        type            prghPressure;
        p               uniform 1e5;
        value           uniform 1e5;
    }
    
    fuselage
    {
        type            fixedFluxPressure;
        value           uniform 0;
    }

Code:

 object      T;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 0 0 1 0 0 0];

internalField   uniform 282.214;

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
        
    inlet
    {
        type            fixedValue;
        value           $internalField;
    }
    
    outlet
    {
        type            zeroGradient;
    }
    
    fuselage
    {
        type            zeroGradient;
    }

Code:

object      U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform (59.161 0 0);

boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
    clplane
    {
        type            slip;
    }
    
    frontier
    {
        type            slip;
    }
    
    inlet
    {
        type            fixedValue;
        value           uniform (59.161 0 0);
    }
    
    outlet
    {
        type            zeroGradient;
     }
    
    fuselage
    {
        type            noSlip;
    }

Okay, now when I try to run it, it fails at the second time step:

Code:

Create time

Create mesh for time = 0


SIMPLE: Convergence criteria found
        p_rgh: tolerance 0.0001
        U: tolerance 1e-05
        e: tolerance 1e-05
        "(k|epsilon|omega)": tolerance 0.001

Reading thermophysical properties

Selecting thermodynamics package 
{
    type            heRhoThermo;
    mixture         pureMixture;
    transport       const;
    thermo          hConst;
    equationOfState perfectGas;
    specie          specie;
    energy          sensibleEnthalpy;
}

Reading field U

Reading/calculating face flux field phi

Creating turbulence model

Selecting turbulence model type RAS
Selecting RAS turbulence model kEpsilon
RAS
{
    model           kEpsilon;
    turbulence      on;
    printCoeffs     on;
    Cmu             0.09;
    C1              1.44;
    C2              1.92;
    C3              0;
    sigmak          1;
    sigmaEps        1.3;
}

Creating thermophysical transport model

Selecting thermophysical transport type RAS
Selecting default RAS thermophysical transport model eddyDiffusivity

Reading g

Reading hRef
Calculating field g.h


Reading pRef
Reading field p_rgh

No MRF models present

Radiation model not active: radiationProperties not found
Selecting radiationModel none
Creating finite volume options from "constant/fvOptions"

Selecting finite volume options model type explicitPorositySource
    Source: porosity1
    - selecting cells using cellZone porousZone
    - selected 4197 cell(s) with volume 0.00219013
Porosity region porosity1:
    selecting model: DarcyForchheimer
    creating porous zone: porousZone

Starting time loop

Time = 0.25

smoothSolver:  Solving for Ux, Initial residual = 1, Final residual = 0.00181065, No Iterations 1
smoothSolver:  Solving for Uy, Initial residual = 1, Final residual = 0.0963853, No Iterations 1
smoothSolver:  Solving for Uz, Initial residual = 1, Final residual = 0.00167042, No Iterations 1
smoothSolver:  Solving for h, Initial residual = 1, Final residual = 0.0197347, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.999649, Final residual = 0.00553735, No Iterations 95
time step continuity errors : sum local = 0.0254991, global = -1.95086e-05, cumulative = -1.95086e-05
smoothSolver:  Solving for epsilon, Initial residual = 0.586636, Final residual = 0.0554639, No Iterations 1
bounding epsilon, min: -704323 max: 612853 average: 458.189
smoothSolver:  Solving for k, Initial residual = 1, Final residual = 0.0828106, No Iterations 2
bounding k, min: -450.635 max: 413.426 average: 2.56551
ExecutionTime = 2.54 s  ClockTime = 2 s

Time = 0.5

smoothSolver:  Solving for Ux, Initial residual = 0.995391, Final residual = 6.33002e+237, No Iterations 1000
smoothSolver:  Solving for Uy, Initial residual = 0.999018, Final residual = 6.78445e+237, No Iterations 1000
smoothSolver:  Solving for Uz, Initial residual = 0.996324, Final residual = 7.79098e+237, No Iterations 1000
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2  ? in "/lib64/libc.so.6"
#3  void Foam::magSqr<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh>(Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&, Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh> const&) at ??:?
#4  Foam::tmp<Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh> > Foam::magSqr<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh>(Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh> const&) at ??:?
#5  ? at ??:?
#6  __libc_start_main in "/lib64/libc.so.6"
#7  ? at /home/abuild/rpmbuild/BUILD/glibc-2.31/csu/../sysdeps/x86_64/start.S:122

My question is simply: I know there are problems with the BCs, but I need help to resolve them. Any and all guidance would be much appreciated!

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Setting up a transient case with laminar incompressible flow	VentinS	OpenFOAM Running, Solving & CFD	10	May 19, 2019 12:26
Is Playstation 3 cluster suitable for CFD work	hsieh	OpenFOAM	9	August 16, 2015 14:53
Setting up Case for different turbulence models with buoyancy-driven flow	MisterX	OpenFOAM Pre-Processing	1	November 19, 2012 05:10
Setting up a turbine nozzle case	cfdengineer	CFX	2	August 2, 2010 05:49
Free surface boudary conditions with SOLA-VOF	Fan	Main CFD Forum	10	September 9, 2006 12:24