CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Bugs

Velocity/pressure/Temperature blow-up while dealing with unstructured meshes

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

Reply
 
LinkBack Thread Tools Display Modes
Old   October 6, 2011, 16:10
Default Velocity/pressure/Temperature blow-up while dealing with unstructured meshes
  #1
Senior Member
 
Pawel Sosnowski
Join Date: Mar 2009
Location: Trieste, Italy
Posts: 105
Rep Power: 9
psosnows is on a distinguished road
Good day!

We would like to present here a very disturbing issue regarding solving cases with unstructured meshes with OpenFOAM.
We are holding from calling it a "major bug" just because we need someone to confirm the presented results.
The issue arises with usage of non-regular meshes for solving buoyant problems.

At certain moment during a simulation, the pressure, temperature and velocity field blows up. This always starts in a single, random, unstructured element within the fluid.

THE SOLVER
For presentation of the bug, we use an incompressible PISO-based, transient, buoyant solver (attached to this thread in the .zip).
The Velocity equation:
Code:
fvVectorMatrix UEqn
(
    fvm::ddt(U)
  + fvm::div(phi, U)
  - fvm::laplacian(mu/rho, U)
  + g * betaT * (T - T0)
);
The Temperature Equation:
Code:
volScalarField k = (K/(rho*Cp));
k.correctBoundaryConditions();
fvScalarMatrix TEqn
(
    fvm::ddt(T)
  + fvm::div(phi,T)
 ==
    fvm::laplacian(k,T)
);
where U is a volVectorField, phi is standard-calculated, mu, rho, Cp, T, K are volScalarFields, g is a dimmensionedVector and T0 is a dimmensionedScalar.

THE TEST CASE
Geometry:
- the test case consists of a 3D fluid box of dimension 0.5 x 0.5 x 0.5;
- front and back are cyclic.
Velocity:
- top, bottom, left and right patches have fixedValue (0 0 0) no-slip condition;
Temperature:
- top and bottom are zeroGradient;
- left has fixed value 1;
- right has fixed value 0;
- for each new mesh, the steady, linear distribution of temperature field was obtained using laplacianFoam for initial condition;
Other fields' amd properties set-up:
- top, bottom, left and right are zeroGradient;
- density rho [1 -3 0 0 0 0 0] 1000;
- specific heat Cp [0 2 -2 -1 0 0 0] 1000;
- dynamic visc. mu [1 -1 -1 0 0 0 0] 0.1;
- thermal cond. K [1 1 -3 -1 0 0 0] 1;
- therm. expansion coeff betaT [0 0 0 -1 0 0 0] 0.01;
- reference temp. T0 [0 0 0 1 0 0 0] 0.5;
- gravity acc. g [0 1 -2 0 0 0 0] ( 0 -9.81 0 );
- max Courant number maxCoNum 0.3;
Schemes and precisions can be found in the attached cases.

All test cases are attached in the .zip

Steps to reproduce:
- unpack
- blockMesh;
- for "Pyramid" case, run modifyMesh;
- change fvSolution to allow running laplacianFoam;
- run laplacianFoam to acquire linear temperature distribution;
- change back fvSolution;
- use obtained linear T field as initial condition for attached controlBuoyantFoam solver;
- run the cases;

REFERENCE RUN
The reference case was meshed using blockMesh, with regular mesh with number of cells 40 x 40 x 40. The mesh is not fine enough near the walls and will produce errors near the walls, but considering problems presented later, those errors can be omitted (but a question arises- why do they occur anyway?).

CROSS 2x2x2 REFINEMENT
The mesh was refined in single cell lines in 3D, to acquire splitting of one cell into 2x2x2 hexes (see hex2x2x2Mesh.jpg). The simulation was performed giving same results as Reference run.

PIRAMID SPLITTING
One cell of Reference mesh was was split into 6 identical, regular pyramids using modifyMesh tool. The simulation was performed, and after working for a while, the simulation broke. Cells near the refinement developed non-physical large temperature and velocity (see TempThreshold.jpg).
After initial test, the mesh was treated with renumberMesh tool, but the simulation gave same result.

CONCLUSIONS AND REMARKS
The presented cases were designed to visualise the error. We encountered this strange behaviour many times, always while using unstructured mesh. The simulations blew up in random mesh positions for no particular reasons. We have run them with increased number of PISO and non-orthogonal corrections (even up to 16 PISO and 8 non-ort). And we always encountered problems.

Industrial applications in most cases require unstructured grids. It is true, that maybe one can play with the mesh and pray that this time the simulation will not blow up... but if a simple solver does not handle simple mesh, that means something is wrong.

Right now we are investigating the "deep" OF code checking for the reason of the problem. We are looking forward to start some collaboration to eliminate it.

Anxious for Your kind replays,
Pawel S and Colleagues
Attached Images
File Type: jpg TempThreshold.jpg (82.2 KB, 42 views)
File Type: jpg hex2x2x2Mesh.jpg (84.0 KB, 29 views)
Attached Files
File Type: zip BUG_UNSTR.zip (36.1 KB, 4 views)
psosnows is offline   Reply With Quote

Old   October 18, 2011, 09:14
Default
  #2
Senior Member
 
kmooney's Avatar
 
Kyle Mooney
Join Date: Jul 2009
Location: Amherst, MA USA - San Diego, CA USA
Posts: 279
Rep Power: 9
kmooney is on a distinguished road
Unstructured grids tend to be sensitive to your choice of discretization schemes.

Post #5 by Dr. Patterson in the thread below contains a few recommended schemes for skewed grids. You might want to give them a try if you haven't yet.

http://www.cfd-online.com/Forums/ope...red-grids.html
kmooney is offline   Reply With Quote

Old   October 18, 2011, 09:21
Default
  #3
Senior Member
 
Pawel Sosnowski
Join Date: Mar 2009
Location: Trieste, Italy
Posts: 105
Rep Power: 9
psosnows is on a distinguished road
Thank you for the answer!
I am investigating this subject right now.
psosnows is offline   Reply With Quote

Old   October 19, 2011, 14:48
Default
  #4
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,894
Rep Power: 26
alberto will become famous soon enoughalberto will become famous soon enough
Since it might interest others, there is a bug-report http://www.openfoam.com/mantisbt/view.php?id=310
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as live DVD/USB, hard drive image and virtual image.
OpenQBMM - An open-source implementation of quadrature-based moment methods
alberto is offline   Reply With Quote

Old   October 19, 2011, 15:48
Default
  #5
Senior Member
 
Alberto Passalacqua
Join Date: Mar 2009
Location: Ames, Iowa, United States
Posts: 1,894
Rep Power: 26
alberto will become famous soon enoughalberto will become famous soon enough
To put in practice what Henry suggested, use a limited scheme for the divergence. If your Pe is high, central differencing will give you problems.

Additionally, on unstructured grids, it is more appropriate to use leastSquares for gradients, since it preserves second-order accuracy.

You find a working set of schemes which fix the problem in the attachment.
__________________
Alberto Passalacqua

GeekoCFD - A free distribution based on openSUSE 64 bit with CFD tools, including OpenFOAM. Available as live DVD/USB, hard drive image and virtual image.
OpenQBMM - An open-source implementation of quadrature-based moment methods
alberto is offline   Reply With Quote

Reply

Thread Tools
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 On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Merge 2 unstructured meshes vinc ANSYS Meshing & Geometry 0 May 19, 2011 16:14
structured vs unstructured hexahedral meshes Usman FLUENT 0 January 19, 2008 05:20
Structure vs unstructured meshes Des Aubery Main CFD Forum 2 July 25, 2003 23:11
artifical dissipation on unstructured meshes Frederic Main CFD Forum 0 October 31, 2001 18:42
Higher Order FV Schemes for unstructured meshes Apurva Shukla Main CFD Forum 4 December 15, 2000 10:17


All times are GMT -4. The time now is 13:49.