[fst-student]

Hello,

I'm using buoyantSimpleFoam for a simulation of trough a channel. In simplefoam the simulation works so I transfered it to buoyantSimpleFoam because the temperature has to be calculated. The simulations runs but after a view interations the maximum value of the density increases very fast and results an error.
I only can run these interations by using a relaxtionfactor of 0.1 for the density.
The BC for the walls, the inlet and outlet I copied from the circuitBoardCooling test case.
The last 2 time steps are shown below:

Time = 135.6

DILUPBiCG:  Solving for Ux, Initial residual = 0.011047, Final residual = 6.65335e-05, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.0125556, Final residual = 0.000104728, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.0146263, Final residual = 0.00126821, No Iterations 1
GAMG:  Solving for p_rgh, Initial residual = 0.0283624, Final residual = 0.000223, No Iterations 6
time step continuity errors : sum local = 0.00459438, global = -0.000824157, cumulative = 0.00564771
rho max/min : 12.461 0.832619
DILUPBiCG:  Solving for nuTilda, Initial residual = 0.00074201, Final residual = 2.89444e-05, No Iterations 1
ExecutionTime = 255.28 s  ClockTime = 259 s

Time = 135.7

DILUPBiCG:  Solving for Ux, Initial residual = 0.00981575, Final residual = 5.17544e-05, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.0084177, Final residual = 9.80715e-05, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.0165464, Final residual = 0.00143128, No Iterations 1


--> FOAM FATAL ERROR: 
Maximum number of iterations exceeded

    From function thermo<Thermo, Type>::T(scalar f, scalar T0, scalar (thermo<Thermo, Type>::*F)(const scalar) const, scalar (thermo<Thermo, Type>::*dFdT)(const scalar) const, scalar (thermo<Thermo, Type>::*limit)(const scalar) const) const
    in file /home/openfoam/OpenFOAM/OpenFOAM-2.3.1/src/thermophysicalModels/specie/lnInclude/thermoI.H at line 76.

FOAM aborting

#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::error::abort() at ??:?
#2  Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy>::THs(double, double, double) const at ??:?
#3  Foam::heRhoThermo<Foam::rhoThermo, Foam::pureMixture<Foam::constTransport<Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > >::calculate() at ??:?
#4  Foam::heRhoThermo<Foam::rhoThermo, Foam::pureMixture<Foam::constTransport<Foam::species::thermo<Foam::hConstThermo<Foam::perfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > >::correct() at ??:?
#5  
 at ??:?
#6  __libc_start_main in "/lib/i386-linux-gnu/libc.so.6"
#7  
 at ??:?
Aborted (core dumped)

Do anyone have an idea to solve that problem?

[Hip2BL7]

Hi everyone,

I've also come across this issue a couple of times, and managed to fix it by underrelaxing rho to quite low levels (<0.05). This seems like the incorrect solution though, and I'm wondering if I've missed anything in terms of the calculation which may explain why my density explodes? It seems to interplay with the calculation for hydrostatic pressure p_rgh, just from my messing about.

My process to solving these is now turning off any update in density by setting the relaxation to 0, then slowly opening the floodgates once the solution stabilises. I have not been able to run this case buoyantPimpleFoam.

For reference, my case is simple 2D incompressible air flowing between two parallel plates under forced convection, where the bottom wall is "cold" and the top wall is adiabatic.

If anyone out there knowledgeable in the inner workings of OpenFoam could help us out, that would be kindly appreciated!

e.g.

Starting time loop

Time = 0.01

DILUPBiCGStab:  Solving for h, Initial residual = 1, Final residual = 3.91784e-09, No Iterations 1
GAMG:  Solving for p_rgh, Initial residual = 1, Final residual = 9.91787e-05, No Iterations 13
GAMG:  Solving for p_rgh, Initial residual = 1.15545e-05, Final residual = 9.92366e-09, No Iterations 12
GAMG:  Solving for p_rgh, Initial residual = 1.00347e-08, Final residual = 5.36769e-09, No Iterations 1
time step continuity errors : sum local = 7.43954e-11, global = -7.51005e-18, cumulative = -7.51005e-18
rho min/max : 0.607709 1.27515
DILUPBiCGStab:  Solving for epsilon, Initial residual = 0.999998, Final residual = 3.12551e-09, No Iterations 1
DILUPBiCGStab:  Solving for k, Initial residual = 1, Final residual = 3.65876e-10, No Iterations 1
ExecutionTime = 0.99 s  ClockTime = 2 s

Time = 0.02

DILUPBiCGStab:  Solving for h, Initial residual = 0.482545, Final residual = 3.07279e-09, No Iterations 3
GAMG:  Solving for p_rgh, Initial residual = 0.999935, Final residual = 9.74004e-05, No Iterations 62
GAMG:  Solving for p_rgh, Initial residual = 4.85622e-05, Final residual = 1.66795e-06, No Iterations 100
GAMG:  Solving for p_rgh, Initial residual = 1.66941e-06, Final residual = 5.74829e-07, No Iterations 100
time step continuity errors : sum local = 1.18707e-07, global = -4.26698e-17, cumulative = -5.01798e-17
rho min/max : 0.607626 1.27515
DILUPBiCGStab:  Solving for epsilon, Initial residual = 0.108103, Final residual = 2.0814e-13, No Iterations 1
DILUPBiCGStab:  Solving for k, Initial residual = 0.287233, Final residual = 5.50281e-09, No Iterations 2
ExecutionTime = 2.98 s  ClockTime = 4 s

Time = 0.03

DILUPBiCGStab:  Solving for h, Initial residual = 0.0637923, Final residual = 5.56083e-10, No Iterations 2
GAMG:  Solving for p_rgh, Initial residual = 0.814782, Final residual = 5.55235e-05, No Iterations 13
GAMG:  Solving for p_rgh, Initial residual = 0.187916, Final residual = 1.57453e-05, No Iterations 19
GAMG:  Solving for p_rgh, Initial residual = 0.000180221, Final residual = 1.15654e-08, No Iterations 19
time step continuity errors : sum local = 2.71784e-09, global = -5.1816e-18, cumulative = -5.53614e-17
rho min/max : 0.60756 311.75
DILUPBiCGStab:  Solving for epsilon, Initial residual = 1, Final residual = 2.53892e-14, No Iterations 1
DILUPBiCGStab:  Solving for k, Initial residual = 0.410838, Final residual = 1.1645e-12, No Iterations 1
ExecutionTime = 3.41 s  ClockTime = 5 s

Time = 0.04

DILUPBiCGStab:  Solving for h, Initial residual = 0.999803, Final residual = 1.13088e-09, No Iterations 4
GAMG:  Solving for p_rgh, Initial residual = 0.999998, Final residual = 0.000159061, No Iterations 100
GAMG:  Solving for p_rgh, Initial residual = 7.80899e-05, Final residual = 1.35153e-06, No Iterations 100
GAMG:  Solving for p_rgh, Initial residual = 1.35242e-06, Final residual = 5.66887e-08, No Iterations 100
time step continuity errors : sum local = 0.0398109, global = 3.56945e-13, cumulative = 3.5689e-13
rho min/max : -48.388 560.131
DILUPBiCGStab:  Solving for epsilon, Initial residual = 0.166943, Final residual = 4.04283e-12, No Iterations 1
DILUPBiCGStab:  Solving for k, Initial residual = 0.00298645, Final residual = 1.04669e-11, No Iterations 1
ExecutionTime = 5.65 s  ClockTime = 7 s

Time = 0.05

DILUPBiCGStab:  Solving for h, Initial residual = 0.975413, Final residual = 2.18205e-09, No Iterations 4


--> FOAM FATAL ERROR: 
Negative initial temperature T0: -20348.4

[Tobermory]

Interesting - I haven't played with the buoyant[S|P]impleFoam solvers, but from the brief snippet of the run log that you included it looks like the energy solver gets out of control, resulting in poor temperature predictions and therefore poor density predictions.

Strongly buoyant flows can be a real pig to solve in steady state - in other codes, I have had to resort to a pseudo-transient approach (ie LTS in OpenFOAM) ... but this should be stable since it is the lower wall that is cold, so I don't think that that is the problem.

Are you underrelaxing enthalpy sufficiently? Have you tried playing with the divSchemes on h? Do the temperatures look sensible? If you stop your example simulation after t=0.03 and save, where is the solution going nuts? At the cold wall, I presume, since the densities are already too high there - but is it adjacent to the domain inlet/outlet (ie is it caused by your BCs)?

This might give you a few things to check out, but is probably not much help I am afraid. You could always upload a simple test case that demonstrates the problem and someone on the forum might have time to try debug it for you.

[Hip2BL7]

Originally Posted by Tobermory

Interesting - I haven't played with the buoyant[S|P]impleFoam solvers, but from the brief snippet of the run log that you included it looks like the energy solver gets out of control, resulting in poor temperature predictions and therefore poor density predictions.

Strongly buoyant flows can be a real pig to solve in steady state - in other codes, I have had to resort to a pseudo-transient approach (ie LTS in OpenFOAM) ... but this should be stable since it is the lower wall that is cold, so I don't think that that is the problem.

Are you underrelaxing enthalpy sufficiently? Have you tried playing with the divSchemes on h? Do the temperatures look sensible? If you stop your example simulation after t=0.03 and save, where is the solution going nuts? At the cold wall, I presume, since the densities are already too high there - but is it adjacent to the domain inlet/outlet (ie is it caused by your BCs)?

This might give you a few things to check out, but is probably not much help I am afraid. You could always upload a simple test case that demonstrates the problem and someone on the forum might have time to try debug it for you.

I have lowered enthalpy and T independently, but no further than ~0.5, I'll give it a go and report back I have in fact limited the maximum temperature through fvOptions but to no avail, it ran for a couple more simulations but still broke.

What divSchemes do you suggest I give a try? I had them all set to upwind on account of the errors.

Its definitely something to do with the inlet. The pressure (p, p_rgh) seems to go mental in the first couple of iterations which then seems to affect other parameters downstream. I've tried dozens of different BC's but none seem to have any effect on the stability directly (the issue seems to source on the corners of the inlet).

EDIT: Changing enthalpy relaxation all the way down to 0 had little effect on convergence.

I'm in the process of writing a fvOptions correction to limit the density after each iteration, maybe this'll work...

[Hip2BL7]

So just thought I'd come back to close this thread off with my findings.

My pressure fields seemed quite wild and varied significantly between timesteps, and as the system initialised in the first couple of time steps, fast oscillations in the pressure value would eventually trip U, and the turbulence fields eventually causing some floating point errors in the energy equation, leading to the solver to break. My solution includes:

• PBiCGStab solving for p_rgh
• maxIter=1 on p_rgh field
• upwind divschemes (though I can make them work barely with linear schemes)
• choice initial conditions
• lowered the relaxation on all fields

Even with these changes, I still have quite large oscilations in the pipe - but they settle down over time. I can see reflections occuring at the end of the pipe, I was hoping there would be some transmissible boundary condition which could help get rid of these oscillations but couldn't find anything on it.