[raghav.venky]

I'm trying to model a solar chimney based on Boussinesq model. After a few time steps I get this error:

Time = 13

DILUPBiCG:  Solving for Ux, Initial residual = 9.40767e-05, Final residual = 9.32478e-07, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.000374868, Final residual = 7.71871e-06, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.00110305, Final residual = 2.33651e-05, No Iterations 1
DILUPBiCG:  Solving for T, Initial residual = 0.0675763, Final residual = 0.000322582, No Iterations 1
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2   in "/lib/x86_64-linux-gnu/libc.so.6"
#3  Foam::DICPreconditioner::calcReciprocalD(Foam::Field<double>&, Foam::lduMatrix const&) at ??:?
#4  Foam::DICPreconditioner::DICPreconditioner(Foam::lduMatrix::solver const&, Foam::dictionary const&) at ??:?
#5  Foam::lduMatrix::preconditioner::addsymMatrixConstructorToTable<Foam::DICPreconditioner>::New(Foam::lduMatrix::solver const&, Foam::dictionary const&) at ??:?
#6  Foam::lduMatrix::preconditioner::New(Foam::lduMatrix::solver const&, Foam::dictionary const&) at ??:?
#7  Foam::PCG::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#8  Foam::fvMatrix<double>::solveSegregated(Foam::dictionary const&) at ??:?
#9  Foam::fvMatrix<double>::solve(Foam::dictionary const&) at ??:?
#10  
 at ??:?
#11  
 at ??:?
#12  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#13  
 at ??:?
Floating point exception (core dumped)

This is my checkMesh:

Build  : 2.3.0-f5222ca19ce6
Exec   : checkMesh
Date   : Oct 26 2014
Time   : 18:57:25
Host   : "raghav"
PID    : 5903
Case   : /home/raghav/OpenFOAM/raghav-2.3.0/run/solarChimney
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:           1449
    internal points:  0
    faces:            9305
    internal faces:   6411
    cells:            3929
    faces per cell:   4
    boundary patches: 5
    point zones:      0
    face zones:       0
    cell zones:       0

Overall number of cells of each type:
    hexahedra:     0
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    3929
    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: 1 (OK).

Checking patch topology for multiply connected surfaces...
    Patch               Faces    Points   Surface topology                  
    blackWall           1250     692      ok (non-closed singly connected)  
    glassWall           1328     731      ok (non-closed singly connected)  
    outlet              67       64       ok (non-closed singly connected)  
    inlet               75       69       ok (non-closed singly connected)  
    wall                174      173      ok (non-closed singly connected)  

Checking geometry...
    Overall domain bounding box (0 0 0) (4.2 4 6)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (1.67772e-16 -9.04982e-17 -9.41181e-19) OK.
    Max cell openness = 2.08637e-16 OK.
    Max aspect ratio = 5.98181 OK.
    Minimum face area = 0.00145422. Maximum face area = 0.129964.  Face area magnitudes OK.
    Min volume = 3.4493e-05. Max volume = 0.00620352.  Total volume = 4.8.  Cell volumes OK.
    Mesh non-orthogonality Max: 62.114 average: 21.9416
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.583302 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

I think the mesh is okay. When I see the solution, the velocities are very high, indicating the solution has diverged. Does anyone know what the error means and where the problem is? Case attached.

Thanks

[alexeym]

Hi,

Could you describe the case in more details? From your case files:

1. Boundaries walls, glassWall, and blackWall overlap. Not quite sure it's good decision to have them this way.
2. Thickness of the mesh is 1 cell, are you quite sure this is what you want? And this cell is in direction of temperature gradient.
3. Mesh is tetrahedral, use leastSquares for gradient schemes.
4. You'd like to use k-epsilon turbulence model, are you sure your case is turbulent?

[raghav.venky]

@alexeym

Thanks for your input.

1. I don't understand what you mean by that. The three surfaces are distinct and are either parallel or perpendicular to each other. Can you please what you mean by overlap?

2. Yes I thought that was not good. I used Salome to generate the mesh. I don't know how to refine the mesh in one direction and keep it coarse in the other directions. I'll make the entire mesh fine and then see what happens.

3. I did not know this. Is there a place where I can read up on this?

4. I thought buoyantBoussinesqSimpleFoam has to be turbulent. Is there a way to turn it off? I'll simulate non turbulence and then go on to turbulence.


I'll post the results after I make the changes you've suggested in a few days.

Thanks!

[alexeym]

Hi,

1. It seems, this was my mistake, paraview messed up wall as patch name and wall as a group of patches, so it displayed wall patch as group of patches (wall + backWall + glassWall), that's why I thought they overlap.

2. The geometry of the mesh is very simple, I don't know why you'd like to use Salome instead of blockMesh.

3. http://www.cfd-online.com/Forums/ope...mesh-pipe.html, https://www.google.fr/search?q=openf...s+grad+schemes

4. There's turbulence model called "laminar", if you take a look at the code

tmp<volScalarField> laminar::nut() const
{
    return tmp<volScalarField>
    (
        new volScalarField
        (
            IOobject
            (
                "nut",
                runTime_.timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::NO_WRITE
            ),
            mesh_,
            dimensionedScalar("nut", nu()().dimensions(), 0.0)
        )
    );
}

so nuEff == nu + nut is just nu, i.e. molecular viscosity.

[raghav.venky]

Used blockMesh to generate mesh. Hence not using leastSquares for gradient schemes. Looks like there is no option (ie turbulenceProperties dictionary) to make it laminar for buoyantBoussinesqSimpleFoam. I tried turning off turbulence in the RASProperties (in the constant directory). The solution still diverges. I've attached the case files.

Could the boundary conditions be a problem?

[raghav.venky]

I played around. Set the velocity boundary condition at the inlet and outlet as 0. Now no air can come in or go out. I get a nice circulation of the air inside ie the air near the hot boundary rises and the air near cold boundary falls.

Now, to simulate air entering and and leaving ie driven due to buoyancy, I had given a zeroGradient boundary condition. But this led divergent solutions. Any suggestions?

[alexeym]

You can try starting with transient solver instead of steady state. Also you can increase under-relaxation, though it led to inconsistency between inlet and outlet (i.e. amount of mater added to the system wasn't equal to amount of removed matter) with zeroGradient BCs.

[hanness]

Hi Raghav,

you could try with a fixedValue pressure on the inlet and a zeroGradient on the outlet or vice versa which would be a common setup for flows through a domain. Usually you would provide the velocity accordingly (fixedValue for velocity where pressure is zeroGradient and vice versa)
For open boundaries you could also try a totalPressure boundary condition for p_rgh and inletOutlet or pressureInletOutletVelocity for U.
For buoyancy driven flows exposed to the ambient I usually use a modified totalPressure BC for p_rgh on both inlet and outlet (modified to account for the barometric pressure differences) and pressureInletOutletVelocity for U. That works fine for large scale simulations.

Hopefully one of those combinations works for you.
Hannes

[raghav.venky]

@alexeym and @hanness

Thank you very much for your suggestions. Sorry for the delayed response, I've been busy with other stuff. The closed system (0 velocity inlet and outlet) is very interesting and I'll be working on that for a while. Is there way to get the heat transfer rate or the heat transfer coefficient?

I believe it starts with getting the gradient of the temperature on the surface where I want to find the heat transfer. Is there a way to do that?

Also, I will have to get the heat flux or the heat transfer coefficient. I have no idea how to do that.

[hanness]

Hi Raghav,

have a look at the wallHeatFlux utility coming with OpenFOAM.

Hannes

[raghav.venky]

Thanks hanness, that was helpful

[mathew1105]

hi sorry , how do you calculate the Totalpressure at the outlet and inlet for a buoyancy flow?

I am stock with a similar problem and am facing convergene problems


thank you for your help