CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (http://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Overestimated temperature values (http://www.cfd-online.com/Forums/openfoam-solving/73744-overestimated-temperature-values.html)

 makaveli_lcf March 16, 2010 04:20

Overestimated temperature values

3 Attachment(s)
Hello foamers!

I came up with a following problem: during simulation of pouring hot melt into the casting mold I got a temperature values bigger than pouring temperature! My BCs are only fixedValue at the inlet and at the walls of mold (where temperature is lower than at the inlet), and zeroGradient at nozzel boundaries, flat fixed top "freesurface" and at the outlet.

Here are the snapshot, where red color indicates overheated cells. From top comes nozzel, flat top interface is with zeroGradient BC, vertical mold is with constant temperature BC:

Attachment 2599
Attachment 2600

So you can see from second picture, that difference in temperature is ~30 K over pouring value, which is not appropriated. The problem arises also from other surfaces with zeroGradient BC. For example from outer wall of submerged part of the nozzle:

Attachment 2601

Where can this problem come from? I don't have any sources of heat, and energy equation contains only energy dissipation parts:

DT/Dt = div(alpha_eff grad(T))

My idea was that error comes from gradient schemes due to non-ortogonality of the mesh. I use following fvSchemes:

Code:

```/*--------------------------------*- C++ -*----------------------------------*\ | =========                |                                                | | \\      /  F ield        | OpenFOAM: The Open Source CFD Toolbox          | |  \\    /  O peration    | Version:  1.6                                  | |  \\  /    A nd          | Web:      www.OpenFOAM.org                      | |    \\/    M anipulation  |                                                | \*---------------------------------------------------------------------------*/ FoamFile {     version    2.0;     format      ascii;     class      dictionary;     location    "system";     object      fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes {     default        Euler; } gradSchemes {     default        Gauss linear; } divSchemes {     default        none;     div(phi,U)      Gauss upwind;     div(phi,T)      Gauss upwind;     div(phi,k)      Gauss upwind;     div(phi,epsilon) Gauss upwind;     div(phi,R)      Gauss upwind;     div(R)          Gauss linear;     div((nuEff*dev(grad(U).T()))) Gauss linear; } laplacianSchemes {     default        none;     laplacian(nuEff,U) Gauss linear corrected;     laplacian((1|A(U)),p) Gauss linear corrected;     laplacian(kappaEff,T) Gauss linear corrected;     laplacian(DkEff,k) Gauss linear corrected;     laplacian(DepsilonEff,epsilon) Gauss linear corrected;     laplacian(DREff,R) Gauss linear corrected; } interpolationSchemes {     default        linear; } snGradSchemes {     default        corrected; } fluxRequired {     default        no;     p              ; } // ************************************************************************* //```
Waiting for your suggestions! It is strange, the OF produces error in such simple case.

 alberto March 17, 2010 03:21

I'd give a try with limiters on the gradients and laplacian.

The mesh is not orthogonal, but it is not skewed either from the detail you showed.

Please let us know what you find.

Best,

 makaveli_lcf March 18, 2010 04:55

1 Attachment(s)
Alberto,
I have followed your advice. First please check if I understood you correctly.
So, I modified fvSchemes

Code:

```gradSchemes {     default        Gauss linear;     grad(T)        limited Gauss linear; } laplacianSchemes {  ....... // !!    laplacian(kappaEff,T) Gauss linear corrected;     laplacian(kappaEff,T) Gauss linear limited 0.5; ....... } snGradSchemes { //  !!    default        corrected;     default        limited 0.5; }```
And I still get production of high temperature values, see animation I attached.

Attachment 2628

It represents same geometry with box slice cutted out to see what is going on
You can see still overheated cell indicated, T error ~10 K.

 makaveli_lcf April 1, 2010 10:11

1 Attachment(s)
I have found that not a diffusion, but convectional part produces wrong values.
So I splitted energy equation following way:

dT1/dt + div(phi, T1) = 0 (1)
dT2/dt = div(alpha_eff * grad(T2)) (2)

Below are the results only transport equation (1) at the left picture, and equations (2), which accounts only for diffusion at the right picture respectively. "dT" value means temperature overestimation value (above pouring temperature). Flow is present in both cases, it's direction is shown with arrows in the picture on the right.

Attachment 2768

All schemes which I used for such test are without limiting.
So what are your suggestions in that case? There is nothing about limiting divergence nonorthogonality corrections in either User or Programming Guides.
Here I use for equation (1) following divergence scheme:

Code:

`    div(phi,T1)                  Gauss linearUpwind Gauss linear;`

 alberto April 1, 2010 10:56

Quote:
 Originally Posted by makaveli_lcf (Post 252747) All schemes which I used for such test are without limiting. So what are your suggestions in that case? There is nothing about limiting divergence nonorthogonality corrections in either User or Programming Guides. Here I use for equation (1) following divergence scheme: Code: `    div(phi,T1)                  Gauss linearUpwind Gauss linear;`
Hello,

in the first fvSchemes you posted, you were using upwind interpolation, which is limited by definition, since you simply use the upwind cell center value on the face.

If you use linearUpwind, you can additionally limit the gradient as shown here ( http://www.cfd-online.com/Forums/ope...earupwind.html ), using

Code:

` div(phi,T1)      Gauss linearUpwind cellLimited  Gauss linear 1;`
Best,

 makaveli_lcf April 1, 2010 11:03

Alberto!

Thank you very much for such immediate response! I will try your suggestions and report here the results.

 alberto April 1, 2010 12:25

Just one additional thought: did you try to re-mesh the part? What mesher did you use to generate the grid that gives you problems?

Best,

 makaveli_lcf April 1, 2010 12:44

It was a Gambit, so actually mesh is the third part product, I should not change it...
Btw, I have already tried
Quote:
 div(phi,T1) Gauss linearUpwind cellLimited Gauss linear 1;
and it produces more errors(((

Alberto, what is your opinion regarding reasons, which can cause such problems with a simple scalar transport equation?

Just for additional information, here is output of checkMesh:

Code:

```Mesh stats     points:          378191     faces:            1625937     internal faces:  1540645     cells:            642651     boundary patches: 10     point zones:      0     face zones:      0     cell zones:      0 Overall number of cells of each type:     hexahedra:    227870     prisms:        138900     wedges:        0     pyramids:      1338     tet wedges:    0     tetrahedra:    274543     polyhedra:    0 Checking topology...     Boundary definition OK.     Point usage OK.     Upper triangular ordering OK.     Face vertices OK.     Number of regions: 1 (OK). Checking patch topology for multiply connected surfaces ...     Patch              Faces    Points  Surface topology     SEN_Outer          5072    2916    ok (non-closed singly connected)     Port_Wall          10588    5414    ok (non-closed singly connected)     SEN_Wall            16692    12887    ok (non-closed singly connected)     WW_c                22256    22896    ok (non-closed singly connected)     WW_b                13272    11878    ok (non-closed singly connected)     WW_a                7734    6059    ok (non-closed singly connected)     NarrowWall          2782    3024    ok (non-closed singly connected)     TopWall            2974    2570    ok (non-closed singly connected)     Outlet              3442    2762    ok (non-closed singly connected)     Inlet              480      519      ok (non-closed singly connected) Checking geometry...     Overall domain bounding box (-0.686 -1.05 -0.0816556) (0.686 0.803 0.0820117)     Mesh (non-empty, non-wedge) directions (1 1 1)     Mesh (non-empty) directions (1 1 1)     Boundary openness (-1.55627e-17 2.36049e-16 1.53059e-17) OK.     Max cell openness = 3.32553e-16 OK.     Max aspect ratio = 15.485 OK.     Minumum face area = 1e-06. Maximum face area = 0.000136551.  Face area magnitudes OK.     Min volume = 8.95047e-10. Max volume = 8.44667e-07.  Total volume = 0.0974925.  Cell volumes OK.     Mesh non-orthogonality Max: 65.4515 average: 13.9421     Non-orthogonality check OK.     Face pyramids OK.     Max skewness = 0.959158 OK. Mesh OK. End```
and here's checkMesh -allGeometry -allTopology
Code:

```Mesh stats     points:          378191     faces:            1625937     internal faces:  1540645     cells:            642651     boundary patches: 10     point zones:      0     face zones:      0     cell zones:      0 Overall number of cells of each type:     hexahedra:    227870     prisms:        138900     wedges:        0     pyramids:      1338     tet wedges:    0     tetrahedra:    274543     polyhedra:    0 Checking topology...     Boundary definition OK.     Point usage OK.     Upper triangular ordering OK.     Topological cell zip-up check OK.     Face vertices OK.     Face-face connectivity OK.   <<Writing 3 cells with with single non-boundary face to set twoInternalFacesCells     Number of regions: 1 (OK). Checking patch topology for multiply connected surfaces ...     Patch              Faces    Points  Surface topology                  Bounding box     SEN_Outer          5072    2916    ok (non-closed singly connected)  (-0.157456 -0.22 -0.0476851) (0.157458 4.17144e-15 0.0476866)     Port_Wall          10588    5414    ok (non-closed singly connected)  (-0.153715 -0.22 -0.021562) (0.153715 1.06581e-17 0.021562)     SEN_Wall            16692    12887    ok (non-closed singly connected)  (-0.0846087 -1.24345e-17 -0.0599263) (0.0846087 0.803 0.0599263)     WW_c                22256    22896    ok (non-closed singly connected)  (-0.686 -1.05 -0.0380647) (0.686 4.09273e-15 0.0378848)     WW_b                13272    11878    ok (non-closed singly connected)  (-0.3 -1.05 -0.0766829) (0.3 8.1554e-15 0.0763984)     WW_a                7734    6059    ok (non-closed singly connected)  (-0.1 -1.05 -0.0816556) (0.1 8.8413e-15 0.0820117)     NarrowWall          2782    3024    ok (non-closed singly connected)  (-0.686 -1.05 -0.027) (0.686 1.36424e-15 0.027)     TopWall            2974    2570    ok (non-closed singly connected)  (-0.686 -1.3667e-14 -0.0816556) (0.686 8.8413e-15 0.0820117)     Outlet              3442    2762    ok (non-closed singly connected)  (-0.686 -1.05 -0.027) (0.686 -1.05 0.027)     Inlet              480      519      ok (non-closed singly connected)  (-0.0599263 0.803 -0.0599263) (0.0599263 0.803 0.0599263) Checking geometry...     Overall domain bounding box (-0.686 -1.05 -0.0816556) (0.686 0.803 0.0820117)     Mesh (non-empty, non-wedge) directions (1 1 1)     Mesh (non-empty) directions (1 1 1)     Boundary openness (-1.55627e-17 2.36049e-16 1.53059e-17) OK.     Max cell openness = 3.32553e-16 OK.     Max aspect ratio = 15.485 OK.     Minumum face area = 1e-06. Maximum face area = 0.000136551.  Face area magnitudes OK.     Min volume = 8.95047e-10. Max volume = 8.44667e-07.  Total volume = 0.0974925.  Cell volumes OK.     Mesh non-orthogonality Max: 65.4515 average: 13.9421     Non-orthogonality check OK.     Face pyramids OK.     Max skewness = 0.959158 OK.     Min/max edge length = 0.001 0.01938 OK.     All angles in faces OK.     Face flatness (1 = flat, 0 = butterfly) : average = 0.999998  min = 0.996924     All face flatness OK.     Cell determinant (wellposedness) : minimum: 0 average: 2.69124  ***Cells with small determinant found, number of cells: 3   <<Writing 3 under-determined cells to set underdeterminedCells Failed 1 mesh checks. End```

 alberto April 1, 2010 13:17

Quote:
 Originally Posted by makaveli_lcf (Post 252770) It was a Gambit, so actually mesh is the third part product, I should not change it... Btw, I have already tried and it produces more errors((( Alberto, what is your opinion regarding reasons, which can cause such problems with a simple scalar transport equation?
Well, just to have a better understanding I would change the mesh, adding boundary layers of hexahedral cells in the problematic zone, so to have a more regular grid there.

Was this mesh working OK in FLUENT if you tried it?

Best,
Alberto

 alberto April 3, 2010 13:40

FYI, similar problems on meshes with tets were described her:

http://www.cfd-online.com/Forums/ope...foam-tets.html

Best,

 makaveli_lcf April 27, 2010 11:38

Alberto, thank you for reference!

I will comment there to join our effort.

 makaveli_lcf May 1, 2010 04:48

A have a following question:

Is there any reason why we solve energy equation in buoyantBoussinesqPisoFoam
after momentum predictor and not after pressure and mass flux correction?

Have a nice weekend!

Best,

Alexander

 alberto May 1, 2010 21:31

Quote:
 Originally Posted by makaveli_lcf (Post 257116) Is there any reason why we solve energy equation in buoyantBoussinesqPisoFoam after momentum predictor and not after pressure and mass flux correction?
I would say it is done to obtain updated thermo data, as a function of h, before solving for the pressure equation.

Best,

 All times are GMT -4. The time now is 10:03.