CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Running, Solving & CFD

Energy balance in solid using chtMRSF

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree2Likes
  • 1 Post By Bloerb
  • 1 Post By Bloerb

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   January 15, 2016, 08:42
Default Energy balance in solid using chtMRSF
  #1
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
Hi all,

I simulate a heat exchanger using chtMultiRegionSimpleFoam and OF 2.4.0.
The results are somewhat realistic but I think there might be a problem with thermal convergence in the SOLID because the heat balance doesn't sum up.

heat balance solid:

The heat transferred at water side (blue) fluctuates and tops the heat transferred from exhaust gas side (red) which makes no sense for me



Other parameters as global temperature increase, pressure drop and the residuals point to a converging solution (become constant) although the maximum temperature in the fluid rises continously.

I measure the transferred heat with a functionObject

Code:
    energyAbsorbedSOLID
    {
        type patchExpression;
        outputControlMode timeStep;
        outputInterval 1;
        accumulations (
            sum
        );
        patches (".*");
        region SOLID;
        variables ("kvar=140;"
                   );
    expression "kvar*(snGrad(T))*area()";
        verbose true;
    }
hcl734 is offline   Reply With Quote

Old   January 15, 2016, 08:47
Default
  #2
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
Another thing are the high time step continuity errors, which decrease but are still quite high



Could somebody point me to a good definition of those values so I can understand them better?
hcl734 is offline   Reply With Quote

Old   January 26, 2016, 08:47
Default
  #3
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
It seems to be a matter of convergence.
If I continue the run it becomes constant with a small difference.
Increasing nNonOrthCorrectors in the Solid was also helpful, since more iterations are done for the solid then, which let the case converge faster.
hcl734 is offline   Reply With Quote

Old   January 26, 2016, 11:44
Default
  #4
Senior Member
 
Join Date: Sep 2013
Posts: 228
Rep Power: 12
Bloerb is on a distinguished road
Quote:
Increasing nNonOrthCorrectors in the Solid was also helpful, since more iterations are done for the solid then, which let the case converge faster.
hmm that is true but there is a better way...Please post your fvSchemes and fvSolution for your solid region. Heat transfer calculations need a very low residual to converge. A few tips about this in general:

Check your residuals (v2.4.0 has a neat function object for this). If there is no change check the log file if the h equation does iterate. You might want to use minIter 1 to be safe. The residuals alone however are often misleading. Always use the wallHeatFlux utility to check if the heatflux from one region to the next is the same. Then check again in paraview. If the heat flux isn't smooth your mesh might just be to coarse (always check the y+ value. I'd recommend a value below 1). Problems arise if the mesh isn't identical on bothsides which can be significant. Using prism layers helps alot. Temperature on a few patches or min/max values should also be checked. For all of those there are function objects to use.

In your case it is probably the relaxation factor.
The relaxation factor for h (or e depending on your thermoPhysicalProperties) has a huge influence on convergence speed. I only use them to stop the simulation from crashing in the first few timesteps. I'd fully remove them if possible (in both regions). Another reason might be your mesh or discretisation.
hcl734 likes this.
Bloerb is offline   Reply With Quote

Old   January 27, 2016, 03:58
Default
  #5
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
Hi Stephan,

I am using OF 3.0 now, just to avoid confusion.


fvSchemes for WATER

Code:
ddtSchemes
{
    default steadyState;
}

gradSchemes
{
    default         Gauss linear;
    grad(U)         cellLimited Gauss linear 1;
}

divSchemes
{
    default         none;
    div(phi,U)      bounded Gauss upwind;
    div(phi,h)      bounded Gauss upwind;
    div(phi,k)      bounded Gauss upwind;
    div(phi,K)      bounded Gauss upwind;
    div(phi,epsilon) bounded Gauss upwind;
    div(phi,omega) bounded Gauss upwind;
    div(phi,R)      bounded Gauss upwind;
    div(R)          Gauss linear;
    div((muEff*dev2(T(grad(U))))) Gauss linear;
    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear limited 0.333;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         limited 0.333;
}

fluxRequired
{
    default         no;
    p_rgh;
}

wallDist
{
    method meshWave;
}
fvSolution for WATER

Code:
solvers
{
    rho
    {
        solver          PCG
        preconditioner  DIC;
        tolerance       1e-7;
        relTol          0;
    }

    p_rgh
    {
        solver           GAMG;
        tolerance        1e-7;
        relTol           0.01;

        smoother         DIC;

        cacheAgglomeration true;
        nCellsInCoarsestLevel 3126;
        agglomerator     faceAreaPair;
        mergeLevels      1;

        maxIter          100;
    }

    "(U|h|k|epsilon|omega)"
    {
        solver           PBiCG;
        preconditioner   DILU;
        tolerance        1e-7;
        relTol           0.1;
    }
}

SIMPLE
{
    momentumPredictor on;
    nNonOrthogonalCorrectors 2;
    pRefCell        0;
    pRefValue       100000;
    rhoMin          rhoMin [1 -3 0 0 0] 1000;
    rhoMax          rhoMax [1 -3 0 0 0] 1100;
}

relaxationFactors
{
    fields
    {
        rho             1;
        p_rgh           0.5;
    }
    equations
    {
        U               0.3;
        h               0.7;
        nuTilda         0.7;
        k               0.7;
        epsilon         0.7;
        omega           0.7;
        "ILambda.*"     0.7;
    }
}


relaxationFactors
{
    fields
    {
        rho             1;
        p_rgh           0.4;
    }
    equations
    {
        U               0.5;
        h               0.9;
        nuTilda         0.9;
        k               0.9;
        epsilon         0.9;
        omega           0.9;
        "ILambda.*"     0.9;
    }
}

fvSchemes for SOLID

Code:
ddtSchemes
{
    default     steadyState;
}

gradSchemes
{
    default         Gauss linear;
}

divSchemes
{
    default         none;
}

laplacianSchemes
{
    default         none;
    laplacian(alpha,h)  Gauss linear corrected;
    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         uncorrected;
}

fluxRequired
{
    default         no;
}

wallDist
{
    method meshWave;
}
fvSolution for SOLID
Code:
solvers
{
    h
    {
        solver           PCG;
        preconditioner   DIC;
        tolerance        1e-06;
        relTol           0.1;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 10;
}

relaxationFactors
{
    fields
    {
    }
    equations
    {
        h               0.7;
    }
}
hcl734 is offline   Reply With Quote

Old   January 27, 2016, 04:11
Default
  #6
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
My residuals go down and become constant. For h at 1e-6 but it takes sometime.
I also monitored pressure drop, temperature increase and maximum temperature in fluid and solid which become constant as well.
Heat flux I checked with a functionObject for the SOLID region

Code:
    energyAbsorbedWATER
    {
        type patchExpression;
        outputControlMode timeStep;
        outputInterval 1;
        accumulations (
            sum
        );
        patches ("SOLID_TO_WATER");
        region WATER;
        variables ("kvar=0.5;"
                   );
    expression "(snGrad(T))*area()";
        verbose true;
    }
Meshes on SOLID and WATER side don't match and I am using mappedWall and AMI to couple them.
The average yPlus is below one for all walls.
hcl734 is offline   Reply With Quote

Old   January 27, 2016, 04:15
Default
  #7
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
I also wanted to check the heat balance in the fluid.
And to obtain the enthalpy flow difference between IN- and OUTLET I used

Code:
   energyOutlet
   {
        type patchExpression;
        outputControlMode timeStep;
        outputInterval 1;
        accumulations (
            sum
        );
        patches ("OUTLET");
        region WATER;
    variables (
                   );
    expression "(A+B*T-C*sqr(T))*T*phi";
        verbose true;
   }
witch A, B, C being the coefficients of my cp-polynomial.
But there seems to be some difference between what is going into the water and what I get from subtracting enthalpy flow for inlet and outlet.
hcl734 is offline   Reply With Quote

Old   January 27, 2016, 12:31
Default
  #8
Senior Member
 
Join Date: Sep 2013
Posts: 228
Rep Power: 12
Bloerb is on a distinguished road
If your residuals become constant (in one region) make sure that this is not due to not solving anymore.

Code:
Solving for h, Initial residual = .., Final residual = .., No Iterations 0
Lowering your tolerance and maybe adding a minimum number of iteration avoids this. The solver stops if the resdiual 1e-06 is reached but this is not always already fully converged. The relaxation factor has a huge impact on performance and might just block changes. Therefore I'd set them to 1.0 in all Regions. This will also speed up your simulation quite a bit. (maybe 1000 iterations instead of 10000) If it crashes you might want to increase them to 1.0 a few iterations into solving. This should be done in all regions.
Code:
solvers
{
    h
    {
        solver           PCG;
        preconditioner   DIC;
        tolerance        1e-08;
        relTol           0.1;
        minIter     1;
    }
}

SIMPLE
{
    nNonOrthogonalCorrectors 3;
}

relaxationFactors
{
    fields
    {
    }
    equations
    {
        h               1.0;
    }
}
The fvSchemes and fvSolution files however do look decent. Your y+ value is also absolutely fine. If these changes did not help I suspect some problem with the mesh. Since your mesh differs between regions interpolation losses have to be expected. One way I check if this is the case is looking if the heatlfux calculated with the wallHeatFlux utility looks decent. Often you can see dots or stripes if there is a big difference. This should look smooth like a temperature field.

Do the same problems occur whitout using polynomial transport?
Bloerb is offline   Reply With Quote

Old   January 28, 2016, 04:10
Default
  #9
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
1) Number of iterations

Quote:
Originally Posted by Bloerb View Post
If your residuals become constant (in one region) make sure that this is not due to not solving anymore.

Code:
Solving for h, Initial residual = .., Final residual = .., No Iterations 0
Lowering your tolerance and maybe adding a minimum number of iteration avoids this. The solver stops if the resdiual 1e-06 is reached but this is not always already fully converged.

You are right, thats what happens after some time.

Code:
Solving for solid region SOLID
DICPCG:  Solving for h, Initial residual = 9.999899e-07, Final residual = 9.999899e-07, No Iterations 0
DICPCG:  Solving for h, Initial residual = 9.999774e-07, Final residual = 9.999774e-07, No Iterations 0
DICPCG:  Solving for h, Initial residual = 9.999774e-07, Final residual = 9.999774e-07, No Iterations 0
DICPCG:  Solving for h, Initial residual = 9.999774e-07, Final residual = 9.999774e-07, No Iterations 0
DICPCG:  Solving for h, Initial residual = 9.999774e-07, Final residual = 9.999774e-07, No Iterations 0
DICPCG:  Solving for h, Initial residual = 9.999774e-07, Final residual = 9.999774e-07, No Iterations 0
DICPCG:  Solving for h, Initial residual = 9.999774e-07, Final residual = 9.999774e-07, No Iterations 0
[...]
Do you know how to reduce the tolerance of residual because something like

Code:
    residualControl
    {
        p_rgh          1e-6;
        U              1e-8; 
        T              1e-6; 

        // possibly check turbulence fields
        "(k|epsilon|omega)"    1e-7; 
    }
Doesn't work according to http://www.cfd-online.com/Forums/ope...implefoam.html

Or how I can add a minimum number of iterations?
hcl734 is offline   Reply With Quote

Old   January 28, 2016, 13:24
Default
  #10
Senior Member
 
Join Date: Sep 2013
Posts: 228
Rep Power: 12
Bloerb is on a distinguished road
This is done as mentioned above. You just need to lower your tolerance or add minIter (again highlighted below) in your fvSolution files. This is one of the reasons heat transfer simulations should not be judged by residual alone. Much more import are the temperature values or the heat fluxes. Residual Control only aborts a calculation once the inital residual reaches a certain value. The solver on the other side stops solving once the tolerance or the relative tolerance from one time step to the next is reached . It is however correct that residual control is not implemented in this solver, maybe in part for this reason.
Code:
solvers
{
    h
    {
        ....
        tolerance        1e-08;
        relTol           0.1;
         minIter     1;
    }
}

relaxationFactors
{
    fields
    {
    }
    equations
    {
        h               1.0;
    }
}
hcl734 likes this.
Bloerb is offline   Reply With Quote

Old   February 8, 2016, 12:17
Default
  #11
Member
 
Join Date: May 2015
Posts: 68
Rep Power: 6
hcl734 is on a distinguished road
Hi,

your suggestions have been pure gold as they worked out for my case pretty good.
But now I am having another question:
It seems like I am using a significant amount of energy (17 W) when coupling solid and fluid region with compressible::turbulentTemperatureCoupledBaffleMix ed.
Some loss makes sense since I have Arbitrary Mesh Interfaces (AMI) which I couple using nearestPatchFaceAMI. But 17 W is too much and not confirmed by CFX calculations.
So I wonder whether it is a matter of thermophysicalProperties
I am using temperature dependent values for heat capacity (hPolynomial) and coupling is specified as
Code:
    {
        type            compressible::turbulentTemperatureCoupledBaffleMixed;
        value           uniform 348.15;
        Tnbr            T;
        kappa           solidThermo; //fluidThermo
        kappaName       none;
    }
Is this correct?
Is solidThermo temperature dependent or should I use the lookup entry?


Description of coupling condition

Code:
Description                                    Mixed boundary condition for temperature, to be used for heat-transfer                                    on back-to-back baffles.                                
                                   Specifies gradient and temperature such that the equations are the same                                    on both sides:                                    - refGradient = zero gradient                                    - refValue = neighbour value                                    - mixFraction = nbrKDelta / (nbrKDelta + myKDelta())                                
                                   where KDelta is heat-transfer coefficient K * deltaCoeffs                                
                                   Example usage:                                        myInterfacePatchName                                        {                                            type        compressible::turbulentTemperatureCoupledBaffleMixed;                                            neighbourFieldName  T;                                            kappa               lookup;                                            kappaName           kappa;                                            value               uniform 300;                                        }                                
                                   Needs to be on underlying mapped(Wall)FvPatch.                                
                                   Note: kappa : heat conduction at patch. Gets supplied how to lookup                                        calculate kappa:                                    - 'lookup' : lookup volScalarField (or volSymmTensorField) with name                                    - 'fluidThermo' : use fluidThermo and compressible::RASmodel to calculate                                        kappa                                    - 'solidThermo' : use solidThermo kappa()                                    - 'directionalSolidThermo' directionalKappa()                                
                                   Note: runs in parallel with arbitrary decomposition. Uses mapped                                    functionality to calculate exchange.
hcl734 is offline   Reply With Quote

Old   February 8, 2016, 15:03
Default
  #12
Senior Member
 
Join Date: Sep 2013
Posts: 228
Rep Power: 12
Bloerb is on a distinguished road
solidThermo is correct for solid regions and fluidThermo correct for fluid regions. This means that kappa is determined via the thermophysical model used. (In your case the polynomial ones) Since the "turbulentTemperatureCoupledBaffleMixed" boundary condition can be used with other solvers and solvers written by yourself it has the lookup functionality. lookup sets kappa to the value specified in the kappa file. (You need to create one and probably modify the chtmultiregion solvers to use it). The name of the file needs to be set by the kappaName entry. Which is why it right now is set to none. If no file is present lookup should however be the same as solidThermo.

One way to check where your looses occur is to use the mapFields utility to map one boundary field to the neighbour region and afterwards subtract them with foamCalc. With that you should see problematic areas. The syntax for this is a bit complicated but you could give it a try.

Are the results using the nearestPatchFaceAMI better than with the other options available or why did you switch? And can you post a closeup of the boundary between the regions so I can get a feeling for the differences between the mesh points of the two regions?

I do not expect your polynomial to influence your losses.
Bloerb is offline   Reply With Quote

Reply

Tags
chtmrsf, chtmultiregionsimplefoam, convergence, multi regions

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
HELP needed: Energy balance in gaseous combustion James Willie FLUENT 11 August 10, 2005 05:58
UDS for solid enegy and solid mass balance Arvind Phate FLUENT 2 July 20, 2005 23:53
UDF Scalar Code: HT 1 Greg Perkins FLUENT 8 October 20, 2000 13:40
UDFs for Scalar Eqn - Fluid/Solid HT Greg Perkins FLUENT 0 October 14, 2000 00:03
Why FVM for high-Re flows? Zhong Lei Main CFD Forum 23 May 14, 1999 14:22


All times are GMT -4. The time now is 01:23.