How to use dynamicRefineFvMesh library using OpenFOAM®
OpenFOAM® has a library (dynamicRefineFvMesh) allowing to realize an automatic mesh refinement during the computation process. This library is used with solvers containing “DyM” in their name (pimpleDyMFoam, interDyMFoam ..). For this article, we will illustrate the possibilities offered by this library with the interDyMFoam solver (version compatible with dynamic mesh libraries).
This post has been written with OpenFOAM 1712+ (www.openfoam.com) before dynamic mesh functionality was merged into standard solvers.
The setting of the automatic mesh refinement is relatively simple. Just fill in the dynamicMeshDict dictionary (or dataset) located in the directory /constant. The dictionary requires several entries including a scalar field for the automatic refinement criterion (volume fraction, concentration of a passive scalar ..), a refinement frequency, a maximum number of cells, levels of refinement etc …
[pastacode lang=”cpp” manual=”dynamicFvMesh%20%20%20dynamicRefineFvMesh%3B%0AdynamicRefineFvMeshCoeffs%0A%7B%0A%20%20%20%20%2F%2F%20How%20often%20to%20refine%0A%20%20%20%20refineInterval%20%201%3B%0A%20%20%20%20%2F%2F%20Field%20to%20be%20refinement%20on%0A%20%20%20%20field%20%20%20%20%20%20%20%20%20%20%20alpha.water%3B%0A%20%20%20%20%2F%2F%20Refine%20field%20inbetween%20lower..upper%0A%20%20%20%20lowerRefineLevel%200.001%3B%0A%20%20%20%20upperRefineLevel%200.999%3B%0A%20%20%20%20%2F%2F%20If%20value%20%3C%20unrefineLevel%20unrefine%0A%20%20%20%20unrefineLevel%20%20%2010%3B%0A%20%20%20%20%2F%2F%20Have%20slower%20than%202%3A1%20refinement%0A%20%20%20%20nBufferLayers%20%20%201%3B%0A%20%20%20%20%2F%2F%20Refine%20cells%20only%20up%20to%20maxRefinement%20levels%0A%20%20%20%20maxRefinement%20%20%202%3B%0A%20%20%20%20%2F%2F%20Stop%20refinement%20if%20maxCells%20reached%0A%20%20%20%20maxCells%20%20%20%20%20%20%20%20400000%3B%0A%20%20%20%20%2F%2F%20Flux%20field%20and%20corresponding%20velocity%20field.%20Fluxes%20on%20changed%0A%20%20%20%20%2F%2F%20faces%20get%20recalculated%20by%20interpolating%20the%20velocity.%20Use%20’none’%0A%20%20%20%20%2F%2F%20on%20surfaceScalarFields%20that%20do%20not%20need%20to%20be%20reinterpolated.%0A%20%20%20%20correctFluxes%0A%20%20%20%20(%0A%20%20%20%20(phi%20none)%0A%20%20%20%20(nHatf%20none)%0A%20%20%20%20(rhoPhi%20none)%0A%20%20%20%20(alphaPhi%20none)%0A%20%20%20%20(alphaPhiUn%20none)%0A%20%20%20%20(alphaPhiCorr0%20none)%0A%20%20%20%20(ghf%20none)%0A%0A%20%20%20%20)%3B%0A%20%20%20%20%2F%2F%20Write%20the%20refinement%20level%20as%20a%20volScalarField%0A%20%20%20%20dumpLevel%20%20%20%20%20%20%20true%3B%0A%7D” message=”Dictionnaire dynamicMeshDict situé dans le répertoire constant” highlight=”” provider=”manual”/]
The nBufferLayers parameter specifies the number of “buffer” cells between two refinement levels. This parameter has an identical effect to nCellsBetweenLevels used by snappyHexMesh. Thus, nBufferLayers 1; indicates that there will be one cell between two mesh levels. The maxRefinement parameter determines the maximum number of times a cell can be cut. Here, with maxRefinement 2; a hexahedral type cell can be cut at most twice, thus giving rise to 2 ^ 3 * 2 ^ 3 (64 cells). The maxCells parameter limits the maximum number of cells in the model. It is also possible to choose to correct some flux after cutting a cell by listing them using the correctFluxes keyword. In this case, the flux will be “updated” by interpolating the velocity field to the face.
It is important to note that automatic remeshing under OpenFOAM® has some limitations :
- it only works for 3D case,
- it only uses the hexRef8 class (hexahedral cells) as cells that can be refined. The cells of this type being cut in 3 directions,
- the only refinements which can be removed are those which have been previously generated. It is therefore not possible to remove the refinements added by snappyHexMesh.
- It is also not possible to refine the levels already added by snappyHexMesh (for example with the use of the distance or inside mode) beyond the maximum level maxRefinement.
Automatic refinement can struggle on cells called protected cells. A message similar to the following appears in the terminal when the solving the numerical model:
[pastacode lang=”bash” manual=”Detected%20200%20cells%20that%20are%20protected%20from%20refinement.%20Writing%20these%20to%20cellSet%20protectedCells.%20″ message=”” highlight=”” provider=”manual”/]
The protected cells appear when you start the calculation and you have used renumberMesh (after decomposePar). A set directory is created in each decomposition directory of the calculation case:
[pastacode lang=”bash” manual=”%2FyourCase%2Fprocessor*%2Fconstant%2FpolyMesh%2Fset” message=”” highlight=”” provider=”manual”/]
The two figures below illustrate the problems encountered in the case where protected cells are present:
Abnormal functioning due to protected cells (gauche) – Normal functioning (droite)
Example case : Water jet on an inclined plate
The simulation of a jet (water) on an inclined plate (at 45 ° for example) is a relevant case to illustrate the dynamic libraryRefineFvMesh possibilities. As mentioned before, the automatic refinement is only possible in 3D (it is possible to use it in 2D but stability problems can occur because the cells will still be cut out of the plane).
The mesh is relatively coarse at the beginning (~ 140k cells) in order to limit the calculation time.
Of course, it’s a pure demonstrative case. The physical validity of the calculation is not sought here.
Numerical methodology
The methodology proposed for this case is as follows:
- Use blockMesh to create the necessary background mesh for snappyHexMesh.
- Meshing with snappyHexMesh by “snapping” the surfaces (.stl format, created under Salome) of the inclined plate and the pipe.
- Define boundary conditions, turbulence model, fluid properties, scheme, linear solvers
- Decompose the case on n processors (here 4) with decomposePar.
- Use the interDyMFoam solver on 4 cores.
- Post-process the results with Paraview®.
The following figures show the surfaces (.stl) of the inclined plate and the pipe boundary as well as the mesh. Three boundary layer elements were added with snappyHexMesh.
Input stl surface for snappyHexMesh (left) – Cartesian mesh generated by snappyHexMesh (right)
The turbulence model used is the k-omega SST model.
The boundary conditions are summarized in the following table:
| Patch | U | k | omega | p | nut | alpha |
|---|---|---|---|---|---|---|
| Inlet | fixedValue (5 m/s) | fixedValue | fixedValue | zeroGradient | calculated | fixedValue |
| Outlet/Top/Bottom/Front/Side | pressureInletOutletVelocity uniform (0.0 0 0); | inletOutlet | inletOutlet | totalPressure; p0 uniform 0; | calculated | inletOutlet |
| Plate | noSlip | kqRWallFunction | omegaWallFunction | fixedFluxPressure; | nutkWallFunction | zeroGradient |
| Symmetry | symmetry | symmetry | symmetry | symmetry | symmetry | symmetry |
Note that for numerical stability arguments it is advisable to combine pressureInletOutletVelocity for velocity with totalPressure for dynamic pressure. The reader may refer to the boundary condition calculator for estimating k and omega (https://cfd-training.com/en/turbulent-boundary-conditions-calculator/) at the inlet patch of the computational domain.
The speed pressure coupling is handled with the PIMPLE loop using only one external corrector. Under these conditions, the PIMPLE loop is equivalent to the PISO loop.
[pastacode lang=”cpp” manual=”PIMPLE%0A%7B%0A%20%20%20%20correctPhi%20%20%20%20%20%20%20%20%20%20yes%3B%0A%20%20%20%20nOuterCorrectors%20%20%20%201%3B%0A%20%20%20%20nCorrectors%20%20%20%20%20%20%20%20%203%3B%0A%20%20%20%20nNonOrthogonalCorrectors%200%3B%0A%7D” message=”PIMPLE settings” highlight=”” provider=”manual”/]
For the pressure resolution three internal correctors are used.
For this case example, 2 buffer cells and 3 levels (max) of refinement are defined in the dynamicMeshDict dictionary (the other settings are identical to the dictionary presented at the beginning of the article):
[pastacode lang=”cpp” manual=”%20%20%20%20%2F%2F%20Have%20slower%20than%202%3A1%20refinement%0A%20%20%20%20nBufferLayers%20%20%202%3B%0A%0A%20%20%20%20%2F%2F%20Refine%20cells%20only%20up%20to%20maxRefinement%20levels%0A%20%20%20%20maxRefinement%20%20%203%3B” message=”” highlight=”” provider=”manual”/]
Results
The GIF below shows the evolution of the jet with the automatic refinement of the free surface. Automatic refinement is achieved by referring to 0 mesh level. Thus, as refinement is possible, level 0 cells of the initial mesh can be refined three times. Level 1 cells, twice, etc. Initially, level 3 cells are present in the vicinity of the plate. These cells will not be refined. It should be noted that the number of cells increases during the calculation until about 870k.
Download the tutorial (free)
OpenFOAM® Tutorial : Automatic mesh refinement
0.00€interFoam solver
But I believe it has problems with reconstruction, decomposition and restarting the case. Ready to discuss that.
You say that there might be a problem if you get bad protected cells. Though in the final case you have been able to remove the protected cells from the simulation.
How did you manage to not produce protected cells?
Dear Magnus,
Protected cells are created if you use renumberMesh. In the tutorial case you will see that renumberMesh is not called.
Best regards,
CFD-training admin
In my case with pimpleFoam they are also becoming protected without doing anything between snappyHexMesh and pimpleFoam. Some of these cells are completely flat cells that have been morphed to snap to the surface, and some are square cells in the middle of the domain.
In my case is it wind flow over an inclined plate, quite similar to your setup.
This is very strange. Please write us at contact@cfd-training.com we will take a look.
Best regards,
CFD-training admin
Hi,
Have you tried this simulate the case “Water jet on an inclined plate” for a parabolic velocity profile at the inlet or any other codedfixedvalue inlet conditions. I have a similar problem, I got some compilation errors while using dynamicRefineFvMesh along with parabolic velocity condition at the inlet. Do you experience similar errors in your case as well?
I haven’t done the test but it should work.
Best regards