CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Natural convection on a vertical cylinder (https://www.cfd-online.com/Forums/openfoam-solving/241734-natural-convection-vertical-cylinder.html)

Roxana_Laspiur March 16, 2022 14:53
Natural convection on a vertical cylinder
 
2 Attachment(s)
Hi!

I am new to OpenFoam. I am trying to simulate internal natural convection on a vertical cylinder, but I have a problem with the result. I am working with the buoyantBoussinesqPimpleFoam solver (transient case).
A time large (approx 30 min) quasi-steady, I get a full stratified fluid full but streamlines aren't expected (toroidal shape). I also compared the vertical velocity with radius, with a paper by Hess and Miller, did not come up to expectations. :confused:

I think my problem is in the Dict of fvSolution and fvSchemes.
I am working in the laminar regime Ra=10¹⁰.
BC:
  • wall top and bottom isolated
  • sidewall T= 296 K
Fluid is water and T_0= 293 K


Attached are the Dict and images of the T-profile and the streamlines.
Thanks for your help and sorry for my English.


blockMeshDict:
Code:
FoamFile
{
    version    2.0;
    format      ascii;
    class      dictionary;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

//lado del cuadrado interno
l          0.30; 
nl        -0.30;
s          0.35; 
sn        -0.35; //longitud de la curvatura del cuadrado interior
// datos del cilindro
R          0.5;  //radio
Rn        -0.5;
H          1;  //altura
//arcos del cilindro
X          0.353553;  //x=R*cos 45º
Xn        -0.353553;
Z          0.353553;  //z
Zn        -0.353553;

scale        1 ;
vertices
(
    //Cuadrado interior
    ($l  0 $l)  //0
    ($nl 0 $l) //1
    ($nl 0 $nl)//2
    ($l  0 $nl) //3
   
    ($l  $H $l)  //4
    ($nl $H $l) //5
    ($nl $H $nl)//6
    ($l  $H $nl) //7
    //cilindro externo
    ($X 0 $Z)  //8
    ($X 0 $Zn)  //9
    ($Xn 0 $Zn) //10
    ($Xn 0 $Z)  //11

    ($X $H $Z)  //12
    ($X $H $Zn)  //13
    ($Xn $H $Zn) //14
    ($Xn $H $Z)  //15
);

blocks
(
    hex (0 3 2 1 4 7 6 5) (20 20 35) simpleGrading (1 1 1)
    //block 1
    hex (8 9 3 0 12 13 7 4) (20 10 35) simpleGrading (1 4 1)
    //block 2
    hex (9 10 2 3 13 14 6 7) (20 10 35) simpleGrading (1 4 1)
    //block 3
    hex (10 11 1 2 14 15 5 6) (20 10 35) simpleGrading (1 4 1)
    //block 4
    hex (11 8 0 1 15 12 4 5) (20 10 35) simpleGrading (1 4 1)
);
edges
(
    arc 0 1 (0 0 $s)
    arc 1 2 ($sn 0 0)
    arc 2 3 (0 0 $sn)
    arc 3 0 ($s 0 0)
   
    arc 4 5 (0 $H $s)
    arc 5 6 ($sn $H 0)
    arc 6 7 (0 $H $sn)
    arc 7 4 ($s $H 0)

    //cilindro externo
    arc 8 9  ($R 0 0)
    arc 9 10  (0 0 $Rn)
    arc 10 11 ($Rn 0 0)
    arc 11 8  (0 0 $R)

    arc 12 13 ($R $H 0)
    arc 13 14 (0 $H $Rn)
    arc 14 15 ($Rn $H 0)
    arc 15 12 (0 $H $R)
);

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


ControlDict:
Code:
FoamFile
{
    version    2.0;
    format      ascii;
    class      dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application    buoyantBoussinesqPimpleFoam;

startFrom      latestTime;

startTime      0;

stopAt          endTime;

endTime        3000;

deltaT          0.03;

writeControl    timeStep;

writeInterval  1000;

purgeWrite      0;

writeFormat    ascii;

writePrecision  6;

writeCompression off;

timeFormat      general;


timePrecision  6;

runTimeModifiable true;

adjustTimeStep  no;

maxCo          0.5;

functions
{
    #includeFunc singleGraph1
    #includeFunc singleGraph2
}


fvSchemes:
Code:
ddtSchemes
{
    default        Euler;
}

gradSchemes
{
    default        Gauss linear;
}

divSchemes
{
    default        none;

    div(phi,U)      Gauss QUICK; //upwind
    div(phi,T)      Gauss QUICK; //upwind
    div(phi,k)      Gauss upwind;
    div(phi,epsilon) Gauss upwind;
    div(phi,R)      Gauss upwind;
    div(R)          Gauss linear;
    div((nuEff*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default        Gauss linear uncorrected;
}

interpolationSchemes
{
    default        linear;
}

snGradSchemes
{
    default        uncorrected;
}



fvSolution:
Code:
solvers
{
    p_rgh
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance      1e-6;
        relTol          0.01;
    }

    p_rghFinal
    {
        $p_rgh;
        relTol          0;
    }

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

    "(U|T|k|epsilon|R)Final"
    {
        $U;
        relTol          0;
    }
}

PIMPLE
{
    momentumPredictor no;
    nOuterCorrectors 1;
    nCorrectors    2;
    nNonOrthogonalCorrectors 0;
    pRefCell        0;
    pRefValue      0;
}
relaxationFactors
{
    equations
    {
        "(U|T|k|epsilon|R)" 1;
        "(U|T|k|epsilon|R)Final" 1;
    }
}

Roman1 March 17, 2022 03:31
Hello. May be the thing is of boundary condition. In the picture one can see thin red lines along the sides of the cylinder.
This means that there is no thermal conduction between solid and fluid. Please double check. In case of the correct regime there must be the smooth
temperature transition.

Roxana_Laspiur March 17, 2022 07:31
Roman
thanks for the help.
Regarding the border conditions:
T: On the sidewall I imposed a fixed T, it seemed the right thing to do.
U: on all walls noSlip
alphat, nut, k and epsilon: I didn't modify them, is that my mistake?
alphat:

Code:

internalField  uniform 0;

boundaryField
{
    wallTop
    {
        type            alphatJayatillekeWallFunction;
        Prt            0.85;
        value          uniform 0;
    }
    wallDown
    {
        type            alphatJayatillekeWallFunction;
        Prt            0.85;
        value          uniform 0;
    }
    sides
    {
        type            alphatJayatillekeWallFunction;
        Prt            0.85;
        value          uniform 0;
    }
    defaultFaces
    {
        type            empty;
    }
}
k:
Code:
dimensions      [0 2 -2 0 0 0 0];


internalField  uniform 0.1;

boundaryField
{
    wallTop
    {
        type            kqRWallFunction;
        value          uniform 0.1;
    }
    wallDown
    {
        type            kqRWallFunction;
        value          uniform 0.1;
    }
    sides
    {
        type            kqRWallFunction;
        value          uniform 0.1;
    }
    defaultFaces
    {
        type            empty;
    }
}

nut:
Code:
dimensions      [0 2 -1 0 0 0 0];


internalField  uniform 0;

boundaryField
{
    wallTop
    {
        type            nutkWallFunction;
        value          uniform 0;
    }
    wallDown
    {
        type            nutkWallFunction;
        value          uniform 0;
    }
    sides
    {
        type            nutkWallFunction;
        value          uniform 0;
    }
    defaultFaces
    {
        type            empty;
    }
}
epsilon:

Code:
internalField  uniform 0.01;

boundaryField
{
    wallTop
    {
        type            epsilonWallFunction;
        value          uniform 0.01;
    }
    wallDown
    {
        type            epsilonWallFunction;
        value          uniform 0.01;
    }
    sides
    {
        type            epsilonWallFunction;
        value          uniform 0.01;
    }
    defaultFaces
    {
        type            empty;
    }
}

Roman1 March 17, 2022 08:04
Mostly I use the solver fireFoam but environment dynamics is similar. Important BC is pressure, try different BC beginning from the simpliest ones (fixedFluxPressure for ex, but it changes from the solver to solver). Try running without wall functions. This permits you to start the process as it is. After that you can make computing more sharp with adding wall functions. Try U: slip, too.

Roman1 March 17, 2022 08:24
Useful BC for T is wallHeatTransfer

Roxana_Laspiur March 17, 2022 09:12
Roman1
Changues the BC of T
Code:
sides
    {
        type            fixedValue;
        value          uniform 296;
    }
for new BC of T

Code:
sides
    {
        type              wallHeatTransfer;
        Tinf                uniform 296;
        alphaWall      uniform 5.53e-07;
    }
When running I get the following error

Code:
--> FOAM FATAL ERROR: (openfoam-2012)

    request for compressibleTurbulenceModel turbulenceProperties from objectRegistry region0 failed
    available objects of type compressibleTurbulenceModel are
0()

    From const Type& Foam::objectRegistry::lookupObject(const Foam::word&, bool) const [with Type = Foam::ThermalDiffusivity<Foam::CompressibleTurbulenceModel<Foam::fluidThermo> >]
    in file /home/usuario/OpenFOAM/OpenFOAM-v2012/src/OpenFOAM/lnInclude/objectRegistryTemplates.C at line 463.

Roman1 March 17, 2022 09:36
{
type wallHeatTransfer;
Tinf uniform 500;
alphaWall uniform 1;
}

I could try running your case on OpenFoam 7 or 9, but it can take 2-3 days. Generally, the fixing is to run with different combinations of the parameters.

Roman1 March 17, 2022 09:50
adjustTimeStep no; Try yes

piu58 March 17, 2022 10:59
I run a simlar case, for me a combination of natural convection and external flow. This arises much more strange "effects". I ended with the schemes:

Code:
ddtSchemes {  default Euler;}

gradSchemes { default Gauss linear; }


divSchemes {
        //// tsable and accurate
        div(U)                Gauss  linear upwind default;
        div(k)                Gauss  linear upwind default;
        div(omega)                Gauss  linear upwind default;
    div(phi,U)      Gauss linear upwind grad(U);
        div(phi,omega) Gauss linear upwind default;
        div(phi,k) Gauss linear upwind default;


    div((nuEff*dev2(T(grad(U))))) Gauss linear;

}



snGradSchemes {
    default        limited 1;  // phi=1: 1=accurate
}


laplacianSchemes {
    default        Gauss linear limited  1; // phi=1: accurate
}
I run RAS k-w simulation:

Code:
simulationType  RAS;

RAS {
    RASModel        kOmega;
    turbulence      on;
    printCoeffs    on;
}
I was not very happy with LES experiments.

Roman1 March 18, 2022 03:15
buoyantBoussinesqPimpleFoam: Transient solver for buoyant, turbulent flow of incompressible fluids

Is it correct using incompressible solver for natural convection? May be use another one?

piu58 March 18, 2022 03:26
Yes, I used buoyantBoussinesqPimpleFoam for free convection.

Roxana_Laspiur March 18, 2022 06:55
Yes, I read that it is the right solver.
I still can't solve my problem.
Thanks for your suggestions

piu58 March 18, 2022 10:52
May be your downstream area is too small. The BC at the exit acts a a large distance.

Roman1 March 18, 2022 14:53
This is the OF case with your geometry (link). The case made on the base of the tutorial of OpenFoam ver. 9. All seems work OK.

https://transfiles.ru/svtdo

piu58 March 18, 2022 15:26
We have the situation here, that the case works as intended by the tutorial. But the result does not fulfil the theory.

I worked with that solver (buoyantBoussinesqPimpleFoam) too and found it has problems with free stream simulation. In my case I tried coupling buoyant with additional forced streaming. I was not succesfu with that coupling and had to separate the effects.

I recommend:
1) analyse the buoyant stream. The thickness of boundary layer should increase with the sqare root of the length. And the total heat transfer should be in the near where the Nußelt number calculates it. If taht is the case, all is ok with your boundary conditions.
2) simulate the stream in the upper part of the cylinder spearately. The last upper part of the cylinder soes no change the boundary stream much. You may identify the velocity in the region of 2/3 to 3/4 of the total height of the cylinder. Take that and set it as boundary condition for another case. Simulate this one without buoyancy, and with an ordninary solver like pimpleFoam.
3) If that wokrs and gives the expected result you may try to combine both effects, buoyancy and free stream in one simulation. As said, I was not lucky with that. You should experiment with solvers especially for U and p. If you change anything look first, if the buoyancy effect is still correct (square root and Nußelt condition). Look what works better and what does not.

I am very intersted to hear if you come closer to a one step simulation.


All times are GMT -4. The time now is 14:38.