[arvhult]

Hi guys! I am totally new to this forums, and am looking for some input on a matter of simulating purely buoyancy-driven flow.

I have a model of a vertical pipe with inversely stratified water, the interface being in the centre. The top water is 370 K and the bottom is 470 K. I have used NIST to produce my own polynomials for water properties. The pressure is 7 MPa, and I am absolutely certain the polynomials are valid outside of my temperature interval. The solver is buoyantPimpleFoam, so that heat transfer is also modeled.

So when relaxed, the water should mix by buoyancy forces. My 250 cm long pipe with 20 cm diameter doesn't diverge, but multiplying the length by 8 (and thus also the number of cells) does. The pressure field reaches 1001 iterations in every time step, and after about 0.9 seconds the residual is too great. Any other parameter is solved for in 1 iteration. It is apparent to me that most probably the buoyancy forces become too great, and the vast acceleration of the fluid parcels simply make the physics problematic. I myself understand the physics and fluid mechanics, and recieved help from a teacher in setting the numerics, e.g. correctors and relaxation, etc. My only set boundary conditions except zeroGradient in temperature, and no-slip in U, is 7 MPa pressure on top, and fixedFluxPressure on the walls.

I have set nOuterCorrectors to 3, and nonOrthogonalCorrectors 0 (changing the latter didn't make a difference, I found). I have a relaxation factor of 0.5 for all parameters. in controlDict i juse adjustable time steps, with a maximum Co of 0.5. With 8 CPUs, I simply decompose the pipe vertically and solve in parallel.

Other than the lengthening, the model is completely identical to the one that did not diverge. QUESTION: Is there any way to make this model work, other than refining the mesh so much that it becomes unpractical to finish the simulation within reasonable times? I am running the finer mesh now (3 million cells), but so far there are only slightly lower residuals in p_rgh.

(I am validating the model against an experiment, so I need at least 150 seconds of the simulations. With 8 CPU's at hand, that time is precious to my project.)

Thank you for any input

[arvhult]

Also, I have k-epsilon turbulence, RASmodel and blockMesh. I have temporarily doubled the maxIter for p_rgh, to see how that goes.

[jherb]

I think it would help if you could post you boundary conditions here (the files in the "0" directory). And also the files in the system directory.

[nsf]

Hi Carl-Magnus,

Welcome to the forums!

I think this thread might help you
http://www.cfd-online.com/Forums/ope...watertank.html

Although the title says buoytantSimpleFoam they end up talking about the pimple version of the solver. One thing you could try is to force the dpdt term to zero as suggested here.

Personally I have some doubts about using icoPolynomial with a compressible solver. In your case suspect one would have to develop a compressible version of icoPolynomial. If you attempt that then maybe this thread is of use
http://www.cfd-online.com/Forums/ope...ion-state.html

Best Regards
Nicolas

[arvhult]

Originally Posted by jherb

I think it would help if you could post you boundary conditions here (the files in the "0" directory). And also the files in the system directory.

Hi. My simulation crashed around the same time step, even with increased maxIter. Unfortunately i have no written output to show you, since it was set to begin after 50 seconds.

Courant Number mean: 1.20916 max: 1.74202
deltaT = 0.000133197
Time = 0.948441

diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
PIMPLE: iteration 1
DILUPBiCG:  Solving for Ux, Initial residual = 0.139108, Final residual = 0.00152522, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.34399, Final residual = 0.00131887, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.143764, Final residual = 0.00206552, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.0166886, Final residual = 4.20692e-05, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.00840404, Final residual = 8.39006e-05, No Iterations 1465
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 2.38231e-05, global = 1.97348e-07, cumulative = -3.12036e-06
DICPCG:  Solving for p_rgh, Initial residual = 0.00149478, Final residual = 1.48895e-05, No Iterations 1186
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.43263e-05, global = 4.67573e-07, cumulative = -2.65278e-06
DICPCG:  Solving for p_rgh, Initial residual = 0.00219061, Final residual = 1.88537e-05, No Iterations 1488
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.17406e-05, global = 1.43209e-09, cumulative = -2.65135e-06
PIMPLE: iteration 2
DILUPBiCG:  Solving for Ux, Initial residual = 0.130218, Final residual = 0.00288864, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.270003, Final residual = 0.00209802, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.135979, Final residual = 0.00270099, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.0580598, Final residual = 0.000321248, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.0168534, Final residual = 0.00015614, No Iterations 1473
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 9.48016e-05, global = 5.11128e-08, cumulative = -2.60024e-06
DICPCG:  Solving for p_rgh, Initial residual = 0.000767346, Final residual = 6.73006e-06, No Iterations 1492
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 5.39433e-05, global = 1.50556e-09, cumulative = -2.59873e-06
DICPCG:  Solving for p_rgh, Initial residual = 0.00157791, Final residual = 1.36561e-05, No Iterations 1484
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 5.88858e-05, global = -6.39112e-09, cumulative = -2.60513e-06
PIMPLE: iteration 3
DILUPBiCG:  Solving for Ux, Initial residual = 0.168643, Final residual = 1.9668e-08, No Iterations 4
DILUPBiCG:  Solving for Uy, Initial residual = 0.30858, Final residual = 1.14325e-07, No Iterations 3
DILUPBiCG:  Solving for Uz, Initial residual = 0.160013, Final residual = 9.34459e-07, No Iterations 3
DILUPBiCG:  Solving for h, Initial residual = 0.811487, Final residual = 1.55269e-07, No Iterations 6
DICPCG:  Solving for p_rgh, Initial residual = 0.315944, Final residual = 0.00313795, No Iterations 1464
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.00948333, global = -3.08738e-05, cumulative = -3.34789e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.001084, Final residual = 9.0015e-06, No Iterations 1488
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.0160053, global = -2.87319e-07, cumulative = -3.37662e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.334231, Final residual = 8.98454e-09, No Iterations 1536
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 3.23113e-08, global = 3.7752e-13, cumulative = -3.37662e-05
DILUPBiCG:  Solving for epsilon, Initial residual = 0.985638, Final residual = 9.23328e-07, No Iterations 10
DILUPBiCG:  Solving for k, Initial residual = 0.999651, Final residual = 1.0427e-07, No Iterations 9
bounding k, min: -16.2075 max: 48108.4 average: 6.53818
ExecutionTime = 28571 s  ClockTime = 28812 s

fieldMinMax fieldMinMax1 output:
    min(mag(U)) = 0 at position (0.000140286 0 -0.000392805) on processor 0
    max(mag(U)) = 18873.4 at position (0.00972534 0.448438 -0.00851755) on processor 1
    min(p) = 6.98115e+06 at position (0.0102747 2 0.00851755) on processor 7
    max(p) = 1.11691e+12 at position (0.0096 0 0.0008) on processor 0
    min(T) = 375.435 at position (0.0102747 1.99844 0.00851755) on processor 7
    max(T) = 686.871 at position (0.01 0.0015625 -0.0004) on processor 0
    min(k) = 1e-15 at position (0.000140286 0.0015625 -0.000392805) on processor 0
    max(k) = 48108.4 at position (0.0181196 1.01719 0.00447374) on processor 4
    min(epsilon) = 5.14027e-05 at position (0.0096309 0.0015625 -0.00959039) on processor 0
    max(epsilon) = 55069.9 at position (0.0184835 1.02031 0.00382936) on processor 4
    min(rho) = 60.3145 at position (0.01 0.0015625 -0.0004) on processor 0
    max(rho) = 960.914 at position (0.0102747 1.99844 0.00851755) on processor 7

Courant Number mean: 582.627 max: 1054.19
deltaT = 6.31753e-08
Time = 0.948442

There should be no problems with the boundary conditions, but here is the pressure.

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

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

internalField   uniform 7e+06;

boundaryField
{

    // when setting everything in p_rgh, these settings should all be "calculated"
    // this way conflicts are avoided
    bottom
    {
        type            calculated;
        value           $internalField; // internalField is just a variable from line 20
        // value can always be used and is a useful when having 
        // types that are not supported by paraview, like calculated, 
        // because then there is at least a value at the beginning
    }
    top
    {
        type            calculated;
        value           $internalField;
    }
    fixedWalls
    {
        type            calculated;
        value           $internalField;
    }
}


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

and p_rgh

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

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

internalField   uniform 7e6;

boundaryField
{
    bottom
    {
        type            fixedFluxPressure;
        value           uniform 7e6;
    }

    // one surface should be set to include the pressure if there is no outlet
    top
    {
        type            fixedValue;
        value           uniform 7e6;
    }

    fixedWalls
    {
        type            fixedFluxPressure;
        value           uniform 7e6;
    }
}

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

The temperature is simply mapped with setFields, and boundary conditions to zeroGradient (so no problem there; no heating from walls). U is uniform (0 0 0), and the others (alphat, epsilon, k, mut) are default as in the hotRoom tutorial.

Should I perhaps invest time in implementing the water property data libraries into my OpenFOAM distribution? (Don't recall the name of it, but XSteam utilizes it.)

[arvhult]

Originally Posted by nsf

Hi Carl-Magnus,

Welcome to the forums!

I think this thread might help you
http://www.cfd-online.com/Forums/ope...watertank.html

Although the title says buoytantSimpleFoam they end up talking about the pimple version of the solver. One thing you could try is to force the dpdt term to zero as suggested here.

Personally I have some doubts about using icoPolynomial with a compressible solver. In your case suspect one would have to develop a compressible version of icoPolynomial. If you attempt that then maybe this thread is of use
http://www.cfd-online.com/Forums/ope...ion-state.html

Best Regards
Nicolas

Hi Nicolas,
Thank you for all advice, I will definitely look into it after lunch.
BTW, look at the courant number (code in my post above) - can it not just be that my case has a too great buoyant force? Such a courant number can result from a sudden acceleration, right? I do not know how good the self-adjustable time stepping works in OpenFOAM, but I presume there are limitations.

[nsf]

Hi,

My guess is that the high courant numbers are a symptom, not a cause. Although I have seen cases where the variable time stepping caused problems I don't think it's a problem for your case. However it's so easy to rule out. You could try with a fixed, small time step. Restart the solution from some time where the solution seems stable (maybe t=0.8) and see what happens.

Should I perhaps invest time in implementing the water property data libraries into my OpenFOAM distribution? (Don't recall the name of it, but XSteam utilizes it.)

Although I haven't looked at it my self I, it looks like that's just what Roman did in the last link i posted.

/Nicolas

[arvhult]

Originally Posted by nsf

Hi,

My guess is that the high courant numbers are a symptom, not a cause. Although I have seen cases where the variable time stepping caused problems I don't think it's a problem for your case. However it's so easy to rule out. You could try with a fixed, small time step. Restart the solution from some time where the solution seems stable (maybe t=0.8) and see what happens.



Although I haven't looked at it my self I, it looks like that's just what Roman did in the last link i posted.

/Nicolas

As a nice coincidence, Roman is the one who actually made my model converge in the first place (He was a teacher of mine at KTH), so I have already asked about the property libraries.

Yes, I also understand the courant number is a symptom. Thank you for confirming that the adjustable time step should be good enough!
Since I have no written data due to the early divergence, I cannot restart it close to the point of divergence. But I will definitely re-run it for about 1.5 seconds, and write out every 0.1 seconds, to see if the physics are messed up.

Thank you all for the kind feedback! I will hit you back with any updates.

[arvhult]

Perhaps I should also supply you with my blockMeshDict, although I would not hope it is the problem since I have lengthened a functioning mesh, without ruining the aspect ratio.

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

convertToMeters 1;

vertices
(
    //Bottom circle
    (0 0 0)
    (0.01 0 0.01)
    (0.02 0 0)
    (0.01 0 -0.01)
    //Bottom square
    (0.002 0 0)
    (0.01 0 0.008)
    (0.018 0 0)
    (0.01 0 -0.008)

    //Top circle
    (0 2 0)
    (0.01 2 0.01)
    (0.02 2 0)
    (0.01 2 -0.01)
    //Top square
    (0.002 2 0)
    (0.01 2 0.008)
    (0.018 2 0)
    (0.01 2 -0.008)
);

blocks
(
    //Blocks with grading
    hex (0 4 7 3 8 12 15 11) (6 20 640) simplegrading (2 1 1)
    hex (1 5 4 0 9 13 12 8) (6 20 640) simplegrading (2 1 1)
    hex (2 6 5 1 10 14 13 9) (6 20 640) simplegrading (2 1 1)
    hex (3 7 6  2 11 15 14 10) (6 20 640) simplegrading (2 1 1)
    hex (4 5 6 7 12 13 14 15) (20 20 640) simplegrading (1 1 1)
);

edges
(
    arc 3 0 (0.00292893219 0 -0.00707106781)
    arc 11 8 (0.00292893219 2 -0.00707106781)
    arc 0 1 (0.00292893219 0 0.00707106781)
    arc 8 9 (0.00292893219 2 0.00707106781)
    arc 1 2 (0.01707106781 0 0.00707106781)
    arc 9 10 (0.01707106781 2 0.00707106781)
    arc 2 3 (0.01707106781 0 -0.00707106781)
    arc 10 11 (0.01707106781 2 -0.00707106781)
);

boundary
(
    bottom
    {
        type wall;
        faces
        (
            (0 3 7 4)
            (1 0 4 5)
            (2 1 5 6)
            (3 2 6 7)
            (4 7 6 5)
        );
    }
    top
    {
        type wall;
        faces
        (
            (11 8 12 15)
            (8 9 13 12)
            (9 10 14 13)
            (10 11 15 14)
            (12 13 14 15)
        );
    }
    fixedWalls
    {
        type wall;
        faces
        (
            (0 8 11 3)
            (1 9 8 0)
            (2 10 9 1)
            (3 11 10 2)
        );
    }
);

mergePatchPairs
(
);

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

[arvhult]

Hi all,

I have now some written data on the divergence. This time it stopped at 0.8, again with a high Courant number.

diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
PIMPLE: iteration 1
DILUPBiCG:  Solving for Ux, Initial residual = 0.2464, Final residual = 0.00335823, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.593198, Final residual = 0.00190324, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.247312, Final residual = 0.00349747, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.00869132, Final residual = 2.18017e-05, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.00398637, Final residual = 3.61218e-05, No Iterations 1489
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 2.03653e-05, global = -1.40403e-11, cumulative = -1.02108e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00357569, Final residual = 3.54487e-05, No Iterations 182
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 9.35465e-06, global = 2.01011e-07, cumulative = -1.00097e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00233777, Final residual = 2.24788e-05, No Iterations 1489
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 9.29375e-06, global = -1.3286e-09, cumulative = -1.00111e-05
PIMPLE: iteration 2
DILUPBiCG:  Solving for Ux, Initial residual = 0.0864367, Final residual = 0.00149477, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.224922, Final residual = 0.00107517, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.0844334, Final residual = 0.00146596, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.00691049, Final residual = 3.44991e-05, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.00641169, Final residual = 5.38555e-05, No Iterations 1477
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.93211e-05, global = -2.24201e-08, cumulative = -1.00335e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.000772992, Final residual = 6.98352e-06, No Iterations 1492
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 3.72301e-06, global = -1.50077e-10, cumulative = -1.00336e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00844094, Final residual = 7.77021e-05, No Iterations 1485
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 6.92087e-06, global = -5.63935e-10, cumulative = -1.00342e-05
PIMPLE: iteration 3
DILUPBiCG:  Solving for Ux, Initial residual = 0.0579123, Final residual = 8.39988e-09, No Iterations 3
DILUPBiCG:  Solving for Uy, Initial residual = 0.508397, Final residual = 5.73844e-09, No Iterations 3
DILUPBiCG:  Solving for Uz, Initial residual = 0.056135, Final residual = 4.83073e-09, No Iterations 3
DILUPBiCG:  Solving for h, Initial residual = 0.0361884, Final residual = 2.59179e-08, No Iterations 3
DICPCG:  Solving for p_rgh, Initial residual = 0.0186181, Final residual = 0.00015513, No Iterations 1469
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 3.53444e-05, global = -2.98665e-09, cumulative = -1.00372e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00162671, Final residual = 1.35096e-05, No Iterations 1490
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 6.90039e-05, global = 1.57893e-09, cumulative = -1.00356e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.033054, Final residual = 7.48061e-09, No Iterations 1527
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.73547e-09, global = -5.41923e-14, cumulative = -1.00356e-05
DILUPBiCG:  Solving for epsilon, Initial residual = 0.0695951, Final residual = 1.982e-07, No Iterations 4
DILUPBiCG:  Solving for k, Initial residual = 0.558717, Final residual = 2.36136e-08, No Iterations 5
ExecutionTime = 25219 s  ClockTime = 25547 s

fieldMinMax fieldMinMax1 output:
    min(mag(U)) = 0 at position (0.000140286 0 -0.000392805) on processor 0
    max(mag(U)) = 15.4614 at position (0.016613 1.97656 -0.00715851) on processor 7
    min(p) = -1.2251e+08 at position (0.0096 0.0015625 -8.41452e-21) on processor 0
    max(p) = 6.98107e+06 at position (0.01 2 0.0004) on processor 7
    min(T) = 367.115 at position (0.016613 1.02344 0.00715851) on processor 4
    max(T) = 466.916 at position (0.0198597 0.976562 0.000392805) on processor 3
    min(k) = 2.59436e-05 at position (0.0194437 0.0015625 0.00264287) on processor 0
    max(k) = 0.0688622 at position (0.0103691 0.995312 -0.00959039) on processor 3
    min(epsilon) = 7.90773e-05 at position (0.0096309 0.0015625 -0.00959039) on processor 0
    max(epsilon) = 4.34332 at position (0.0103691 0.995312 -0.00959039) on processor 3
    min(rho) = 876.878 at position (0.0198597 0.976562 0.000392805) on processor 3
    max(rho) = 966.759 at position (0.016613 1.02344 0.00715851) on processor 4

Courant Number mean: 1.61501 max: 2.36562
deltaT = 8.85805e-05
Time = 0.806548

diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
PIMPLE: iteration 1
DILUPBiCG:  Solving for Ux, Initial residual = 0.116203, Final residual = 0.00157145, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.283288, Final residual = 0.00107667, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.116585, Final residual = 0.00160151, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.0216613, Final residual = 5.94309e-05, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.0176146, Final residual = 0.000147573, No Iterations 1491
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 2.20453e-05, global = 1.00842e-09, cumulative = -1.00346e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00166656, Final residual = 1.51823e-05, No Iterations 1493
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.30494e-05, global = 2.15517e-12, cumulative = -1.00346e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00252207, Final residual = 2.37896e-05, No Iterations 1491
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.19989e-05, global = -1.98796e-10, cumulative = -1.00348e-05
PIMPLE: iteration 2
DILUPBiCG:  Solving for Ux, Initial residual = 0.137127, Final residual = 1.64292e-05, No Iterations 2
DILUPBiCG:  Solving for Uy, Initial residual = 0.282967, Final residual = 0.00240918, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.136695, Final residual = 1.75124e-05, No Iterations 2
DILUPBiCG:  Solving for h, Initial residual = 0.243831, Final residual = 0.00152883, No Iterations 1
DICPCG:  Solving for p_rgh, Initial residual = 0.0631348, Final residual = 0.000625358, No Iterations 1484
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.000246716, global = 1.16196e-08, cumulative = -1.00232e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00104445, Final residual = 1.00774e-05, No Iterations 1495
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.000126416, global = -1.36375e-09, cumulative = -1.00245e-05
DICPCG:  Solving for p_rgh, Initial residual = 0.00192377, Final residual = 1.89664e-05, No Iterations 1491
diagonal:  Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.000122684, global = 6.97997e-10, cumulative = -1.00239e-05
PIMPLE: iteration 3
DILUPBiCG:  Solving for Ux, Initial residual = 0.237038, Final residual = 3.90856e-08, No Iterations 3
DILUPBiCG:  Solving for Uy, Initial residual = 0.208437, Final residual = 5.63899e-08, No Iterations 3
DILUPBiCG:  Solving for Uz, Initial residual = 0.232736, Final residual = 4.06552e-08, No Iterations 3
DILUPBiCG:  Solving for h, Initial residual = 0.979527, Final residual = 6.20112e-07, No Iterations 4

The picture shows a couple of slices at the interface between the hot and cold water. Keep in mind a diameter of 2 cm, and pipe height of 2 m. There are maxima of the mixing velocity as seen in sharp red. The opaque volume colour code is hot/cool temperature, but not that important.


greenshot screen capture

Here is a picture of my mesh. I chose the weird design just to get a better aspect ratio of the cells.


how to do a screenshot on a pc

I am preparing a fine mesh, and attempting to solve the coarse one from 0.8 s with a manually set, much smaller time step.

[arvhult]

Hi,
As i suspected, there was nothing wrong with the adjustable time step; even with the smallest manual stepping (1e-7 s) the simulation would not converge, starting from 0.8 s. Still > 1400 iterations for p_rgh.

I am queuing running a simulation with easier properties; atmospheric pressure, thus more incompressible.
Please give feedback if you see something problematic with my mesh; meshing is not my experienced side, but this one (above) seems nice enough for OpenFOAM to handle.

Otherwise I will keep you posted on the progress with the fine mesh run.

[nsf]

Hi,

Did you try setting dpdt off yet? I'm feeling quite confident that your case will run if you do so.

I did try to implement a compressiblePolynomial. It's based on icoPolynomial but additionally a polynomial for the compressibility has to be entered.

I tried it on the hot room tutorial, but even with compressibility added, it crashes. Both icoPolynomial and my compressiblePolynomial runs just fine when dpdt is turned off. The simulation will also diverge as soon as dpdt is turn on if initialized with dpdt off.

I suspect that dpdt will diverge for incomressible liquids as a direct consequence of the definition of incompressible. I mentioned this here.

Now how, one motivates that dpdt can be safely neglected is another thing. I suspect that the same thing is done in fluent. At least from what I understand of the theory manual.

Good luck.
Nicolas

[arvhult]

Originally Posted by nsf

Hi,

Did you try setting dpdt off yet? I'm feeling quite confident that your case will run if you do so.

I did try to implement a compressiblePolynomial. It's based on icoPolynomial but additionally a polynomial for the compressibility has to be entered.

I tried it on the hot room tutorial, but even with compressibility added, it crashes. Both icoPolynomial and my compressiblePolynomial runs just fine when dpdt is turned off. The simulation will also diverge as soon as dpdt is turn on if initialized with dpdt off.

I suspect that dpdt will diverge for incomressible liquids as a direct consequence of the definition of incompressible. I mentioned this here.

Now how, one motivates that dpdt can be safely neglected is another thing. I suspect that the same thing is done in fluent. At least from what I understand of the theory manual.

Good luck.
Nicolas

Dear Nicolas,

Thanks for the heads up on that. When I now came to work, the "easier" simulation above mentioned had converged for 2 s, i.e. past the before time where it diverged. I will restart it now from 2 s; if it diverges, I will turn dpdt off.

Next will also be to implement the water property libraries; I just recieved them from Roman.

I also started two simulations of the actual case my work is focused on, which has a nicer aspect ratio (diameter 2 dm, height 2-3 m). They didn't diverge up to 150 seconds. This is good news, but I still need to validate my model against the nasty one.

I will keep you posted on how the continued run goes - and thank you again for the great help.
Best regards, Carl-Magnus

[arvhult]

Hi everyone;
Posting an update now as promised. My simulation is still running, coming up on its 11th second (in a week...).

So, ehrm, this might take a while. Do you know if a maximum courant number larger than 0.5 would work? (That is kind of a default value I have used from tutorials and other models). Changing that would be the only way to speed up the simulation I believe.

Other than that, thank you for all the help; hopefully the problem is now solved. I will still look into employing the property libraries, to make the physics nicer to the solver.

Best regards,
Carl-Magnus

[nsf]

Hi,

I happy to hear it worked for you. Did you have to turn off dpdt or did it work with it switch on and varying density?

You could get a simulation to "work" with high Co. By increasing the outer iterations and introducing under relaxation. (will this improve the speed of the simulation?). I've seen post in this forums of people running with max Co as high as 2-4.

It's however difficult to say if it's a good idea. The best way to find out is to try.
Why not try to run with Co > 0.5 and compare it to your present case? Just restart with a higher Co. Let it run for a day or two and compare. Was it faster? Did the solution change?

Since you have a simple geometry you probably have high Co in large parts of the domain and any errors introduced by the high Co will affect a large part of the domain. And using high Cos is probably not the best approach. If that's not the case and you have high Co in just a few cells consider remeshing instead of running with high Co.

Best Regards
Nicolas

[arvhult]

Originally Posted by nsf

Hi,

I happy to hear it worked for you. Did you have to turn off dpdt or did it work with it switch on and varying density?

You could get a simulation to "work" with high Co. By increasing the outer iterations and introducing under relaxation. (will this improve the speed of the simulation?). I've seen post in this forums of people running with max Co as high as 2-4.

It's however difficult to say if it's a good idea. The best way to find out is to try.
Why not try to run with Co > 0.5 and compare it to your present case? Just restart with a higher Co. Let it run for a day or two and compare. Was it faster? Did the solution change?

Since you have a simple geometry you probably have high Co in large parts of the domain and any errors introduced by the high Co will affect a large part of the domain. And using high Cos is probably not the best approach. If that's not the case and you have high Co in just a few cells consider remeshing instead of running with high Co.

Best Regards
Nicolas

Dear Nicolas,

I changed the model to atmospheric pressure, and the polynomials as well. Then it became more incompressible, and icoPolynomial probably works better.

The simulation is still in the most transient phase (10-20 s) where most accelerations occur. I may wait until monday before taking further measures; I would expect the water to slow down by then.

Can I turn dpdt off while restarting, or will that mess things up?

Otherwise, I will look into the Co and the property libraries to make things more realistic.

Best regards, Carl

[nsf]

Hi,

I don't think you can turn it off while running. Well maybe you can but it's not a good idea. As far as I understand the code, it will stop updating dpdt if you turn it off while the simulation is running. It won't reset dpdt=0. On the other hand, I don't think you need to restart in order to turn of dpdt. Just stop the simulation.

/Nicolas