[Hawxliquid]

Hi,


I am currently to do a simple case, the tutorial example "simpleCarSolarPanel" from the heattransfer tutorials.


In the tuorial, it seems the temperature is constant everywhere in the domain. But what I would like to simulate is actually the temperature in the domain, as well as its effect on the natural convection it generates.


I thought of using the "qr" as a an input for heat in my domain, but I can not find how to use any boundary conditions that would let me do that..


Would anyone have experience on that? I do not need a very realistic simulation, just a way of having the boundary cells with a high qr to input heat at a higher rate in the domain, that the ones with a lower qr.


Thank you!

[Yann]

Hello,

You can replace the zeroGradient temperature BC by something like externalWallHeatFluxTemperature so your walls get heated up with solar radiation.

Have a look at the documentation: https://doc.openfoam.com/2306/tools/...uxTemperature/

You can also check the tutorial linked in the documentation.

Hope this helps,
Yann

[Hawxliquid]

Hi Yann,


Great. Thanks, it seems it sent me on the good direction


I now have results, but have some issue still.


Namely, I get zone with unrealistic high T°. I tried playing around with the radiation propeties, boundaryRadiationProperties, 0/T, etc.
I manage to lower a bit the high T°, but it is still unrealistically high..


I tried to use fvOptions with limitT, but somehow it does not do anything ? I get a

limitTemperature limitTAir Lower limited 0 (0%) of cells, with min limit 290
limitTemperature limitTAir Upper limited 6911 (19.88%) of cells, with max limit 350
limitTemperature limitTAir Unlimited Tmin 299.642
limitTemperature limitTAir Unlimited Tmax 350

but it does not seem to limit it.


If anyone had any idea why I get such high temperature, or knows a way to apply the limitT (not the best option I know ), It will be of great help!


I am sharing the case here:

[Yann]

It looks like your limitT is working:

limitTemperature limitTAir Upper limited 6911 (19.88%) of cells, with max limit 350

There are 6911 cells were max temperature is limited, and it corresponds to 19.88% of the total number of cells in the mesh.

About the temperature values you get, I guess it must be related to the radiation received, the radiative properties of your walls and the other boundary conditions used in the tutorial.

[Hawxliquid]

Ah, that is peculiar, because when I post-process the case using paraview, I see temperatures way above the limit, as shown in the attached picture (https://imgur.com/a/MjpQIvT). Any idea why?




Additionally, I looked around but could not find good documentation about the boundaryRadiationProperties. Especially about the different types (opaqueDiffusive vs opaqueReflective) and why the absorptivity and emissivity are defined by two floats. (In the example (0.3 0.7)).


Any idea where I can get information about this?




EDIT :


Maybe regarding the limitT, I found a weird issue.
At the beginning of my log, I get the sentence:
"
Selecting finite volume options type limitTemperature
Source: limitTAir
- selecting all cells
- selected 34761 cell(s) with volume 39.8291
"


But it seems that, looking at paraview, my mesh is made of 89843 cells and 51980 points ?!

[Yann]

Originally Posted by Hawxliquid

Ah, that is peculiar, because when I post-process the case using paraview, I see temperatures way above the limit, as shown in the attached picture (https://imgur.com/a/MjpQIvT). Any idea why?

I think this is because limitT limits temperature in the cells. You are visualizing temperature on boundary faces which are probably not limited. The max temperature in your internalMesh should correspond to the max temperature defined in limitT.

Originally Posted by Hawxliquid

Additionally, I looked around but could not find good documentation about the boundaryRadiationProperties. Especially about the different types (opaqueDiffusive vs opaqueReflective) and why the absorptivity and emissivity are defined by two floats. (In the example (0.3 0.7)).

Any idea where I can get information about this?

There are 2 floats because the wall absorption/emission model uses a multiband model. In radiationProperties you also have a spectralDistribution for each band.

If information is not available in the documentation, you have to search into the source code:

grep -r "opaqueDiffusive" $FOAM_SRC

This should get your to something like: $FOAM_SRC/thermophysicalModels/radiation/submodels/boundaryRadiationProperties/opaqueDiffusive/opaqueDiffusive.H
You can check the file header to hopefully gather information about it, and if this is not enough you'll have to dig into the code itself, which might be pretty tedious.

Note that this is only and example, you can also use find to search for specific files and directories, while grep will search a pattern into files. It can be interesting to use both, depending on the situation.

In addition, you might find interesting things in the release notes. For instance the initial release of solarLoad: https://www.openfoam.com/news/main-n...ical-modelling
There must have been some additions in the next releases.

I find radiation models in OpenFOAM pretty messy, good luck getting your way through the jungle!
Yann

[Yann]

Originally Posted by Hawxliquid

EDIT :


Maybe regarding the limitT, I found a weird issue.
At the beginning of my log, I get the sentence:
"
Selecting finite volume options type limitTemperature
Source: limitTAir
- selecting all cells
- selected 34761 cell(s) with volume 39.8291
"


But it seems that, looking at paraview, my mesh is made of 89843 cells and 51980 points ?!

Your mesh is made of 34761 cells, as stated by OpenFOAM. By default, ParaView has a "Decompose polyhedra" option. You have to untick it when loading your case to see the actual mesh stats.

Cheers,
Yann

[Hawxliquid]

Thank you so much for all the information Yann, they are very useful!


I manage to make the simple case work, but not the case I want to apply it to. I have been trying for a while to find where my error could arise from, but to no avail.


When I launch the case (even with radiation turned off), my velocity diverges, as shown here :

https://pastebin.com/23G5mQ6w


I also see at the first iteration that the enthalpy equation seems to diverge, right at the start? Why would that happen? (DILUPBiCGStab: Solving for h, Initial residual = 1, Final residual = 78.6652, No Iterations 1000)










For information, my fvSolution :

 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
 |  \\    /   O peration     | Version:  v2206                                 |
 |   \\  /    A nd           | Website:  www.openfoam.com                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
     version     2.0;
     format      ascii;
     class       dictionary;
     object      fvSolution;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 
 solvers
 {
     p_rgh
     {
         solver          PCG;
         preconditioner  DIC;
         tolerance       1e-08;
         relTol          0.01;
     }
 
 
     "(U|h|k|epsilon)"
     {
         solver           PBiCGStab;
         preconditioner   DILU;
         tolerance        1e-7;
         relTol           0.05;
     }
 }
 
 
 SIMPLE
 {
     nNonOrthogonalCorrectors 10;
 
 
 
 
     rhoMin          0.5;
     rhoMax          2.0;
     residualControl
     {
         p_rgh           1e-2;
         U               1e-4;
         T               1e-2;
 
 
         // possibly check turbulence fields
         "(k|epsilon|omega)" 1e-3;
     }
 }
 
 
 
 
 relaxationFactors
 {
     fields
     {
         p_rgh          0.3;
         p               0.3;
         U               0.7;
         T             0.7;
     }
     equations
     {
         U               0.3;
         T               0.3;
         "(k|epsilon|R)" 0.2;
         rho            0.3;
 
 
     }
 }
 
 
 
 
 // ************************************************************************* //

fvSchemes :

 
 
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
 |  \\    /   O peration     | Version:  v2206                                 |
 |   \\  /    A nd           | Website:  www.openfoam.com                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
     version     2.0;
     format      ascii;
     class       dictionary;
     object      fvSchemes;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 
 ddtSchemes
 {
     default         steadyState;
 }
 
 
 gradSchemes
 {
     default         Gauss linear;
 }
 
 
 divSchemes
 {
     default         none;
 
 
     div(phi,U)      bounded Gauss upwind;
 
 
     energy          bounded Gauss upwind;
     div(phi,h)      $energy;
     div(phi,K)      $energy;
 
 
     turbulence      bounded Gauss upwind;
     div(phi,k)      $turbulence;
     div(phi,epsilon) $turbulence;
 
 
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
 }
 
 
 laplacianSchemes
 {
     default         Gauss linear limited corrected 0.333;
 }
 
 
 interpolationSchemes
 {
     default         linear;
 }
 
 
 snGradSchemes
 {
     default         limited corrected 0.333;
 }
 
 
 fluxRequired
 {
     default         no;
     p_rgh;
 }
 
 
 
 
 // ************************************************************************* //

If anyone had any idea of why this happen, that would be very useful!




Thanks

EDIT:


It seems like that my issue comes from my pressure and p_rgh BC. But I can't seem to find how to fix it. I have an open domain (dome_big), and my BC are:


p

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2006                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 101325;
boundaryField
{
    #includeEtc "caseDicts/setConstraintTypes"
    
   ".*"
   {
       type     calculated;
       value    $internalField;
   }

}

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

prgh

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2206                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p_rgh;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

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

internalField   uniform 101325;

boundaryField
{
    dome_big
    {
        type            prghTotalPressure;
        p0              $internalField;



    }

    ground
    {
        type            fixedFluxPressure;
        value           $internalField;


    }

    buildings
    {
        type            fixedFluxPressure;
        value           $internalField;



    }

    water
    {
        type            fixedFluxPressure;
        value           $internalField;



    }
}


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

[Yann]

Hello,

Originally Posted by Hawxliquid

I also see at the first iteration that the enthalpy equation seems to diverge, right at the start? Why would that happen? (DILUPBiCGStab: Solving for h, Initial residual = 1, Final residual = 78.6652, No Iterations 1000)

This usually happens due to initialization, and it can happen on any variable but pressure and energy are probably the most sensitive. It's not really a big deal as long as it gets back to normal within few iterations, which is what is happening for you, since you get this on the second iteration:

DILUPBiCGStab:  Solving for h, Initial residual = 0.580014, Final residual = 0.0151477, No Iterations 4

So basically the solver struggles to solve h on first iteration, probably due to the initialization of the simulation, but gets back on track and starts converging on the next iterations. So far so good, I wouldn't mind too much about that. Initial conditions are important, but it's not always possible to have a good first guess before starting the simulation.
(By the way, you can monitor residual by plotting it, it makes it easier to follow what's going on in the simulation. I would also recommend to plot variables of interest such as min/max temperature, max velocity, etc...)

Back to your diverging velocity: it can be VERY tricky to find proper boundary conditions for buoyant flows, especially for open domains. When dealing with natural convection, the flow conditions at the boundaries are usually not known and the flow is driven by buoyant forces inside the domain. Numerically speaking, this is not something easy to deal with.

I personally still don't have reliable set of boundary conditions for such cases. You can have a look at this thread which basically discuss this matter: Common setup for p_rgh in complete open domains

Another recommendation: try avoiding 1st order schemes when it comes to natural convection. It can lead to weird stuff like reversed convection flow.

Aside from numerical schemes, I think your main challenge will be to find a proper set of boundary conditions for your case.

Good luck!
Yann

[Hawxliquid]

Originally Posted by Yann

Hello,



This usually happens due to initialization, and it can happen on any variable but pressure and energy are probably the most sensitive. It's not really a big deal as long as it gets back to normal within few iterations, which is what is happening for you, since you get this on the second iteration:

DILUPBiCGStab:  Solving for h, Initial residual = 0.580014, Final residual = 0.0151477, No Iterations 4

So basically the solver struggles to solve h on first iteration, probably due to the initialization of the simulation, but gets back on track and starts converging on the next iterations. So far so good, I wouldn't mind too much about that. Initial conditions are important, but it's not always possible to have a good first guess before starting the simulation.
(By the way, you can monitor residual by plotting it, it makes it easier to follow what's going on in the simulation. I would also recommend to plot variables of interest such as min/max temperature, max velocity, etc...)

Back to your diverging velocity: it can be VERY tricky to find proper boundary conditions for buoyant flows, especially for open domains. When dealing with natural convection, the flow conditions at the boundaries are usually not known and the flow is driven by buoyant forces inside the domain. Numerically speaking, this is not something easy to deal with.

I personally still don't have reliable set of boundary conditions for such cases. You can have a look at this thread which basically discuss this matter: Common setup for p_rgh in complete open domains

Another recommendation: try avoiding 1st order schemes when it comes to natural convection. It can lead to weird stuff like reversed convection flow.

Aside from numerical schemes, I think your main challenge will be to find a proper set of boundary conditions for your case.

Good luck!
Yann

Again, merci for your time and such great information!


I think I am starting to get a okay convergence. If anyone is interested, previously my " far field" model was made of the side and the sky as the same patch, hence making the convergence problematic regarding BC.
I separated the two, and used a prghTotalPressure for the "sky", and fixedFluxPressure on the rest of the BCs. That seem to work so far.


Now, just have to play with the radiation model to get "realistic" temperature on the BCs! Currently they reach thousands of K, for a qr equivalent to a sunny day