next up previous contents
Next: 4 Going further Up: 3 A more complex Previous: 3.3 Visualisation   Contents

3.4 Using dynamic adaptive mesh refinement

For the moment our simulation is not very well resolved. We could always change the GfsRefine 6 line to something bigger but it would not make really good use of the quadtree approach used in Gerris. A code using a simple regular Cartesian grid approach would be faster and would produce the same results. Instead we are going to use dynamic adaptive mesh refinement, where the quadtree structure of the discretisation is used to adaptively follow the small structures of the flow, thus concentrating the computational effort on the area where it is most needed.

This is done using yet another object class: GfsAdapt, also derived from GfsEvent. Various criteria can be used to determine where refinement is needed. In practice, each criterium will be defined through a different object derived from GfsAdapt. If several GfsAdapt objects are specified in the same simulation file, refinement will occur whenever at least one of the criteria is verified.

For this first example, we will use a simple criterium based on the local value of the vorticity. A cell will be refined whenever

$\displaystyle {\vert\nabla\times{\bf v}\vert\Delta x\over\max\vert{\bf v}\vert}$ > $\displaystyle \delta$,

where $ \Delta$x is the size of the cell and $ \delta$ is a user-defined threshold which can be interpreted as the maximum angular deviation (caused by the local vorticity) of a particle traveling at speed max|$ \bf v$| across the cell. This criterium is implemented by the GfsAdaptVorticity object.

The general syntax for an GfsAdapt object is:

[GfsEvent] { mincells = 1 maxcells = 100000 minlevel = 1 maxlevel = 10 cmax = 1e-2 }
where mincells specifies the minimum number of cells in the domain, maxcells the maximum number of cells, minlevel the level below which it is not possible to coarsen a cell, maxlevel the level above which it is not possible to refine a cell and cmax the maximum cell cost. The default values are 0 for minlevel and mincells and infinite for maxlevel and maxcells. An important point is that, for the moment, it is not possible to dynamically refine solid boundaries. A simple solution to this restriction is to always refine the solid boundary with the maximum resolution at the start of the simulation and to restrict the refinement using the maxlevel identifier in GfsAdapt.

What happens if the maximum number of cells is reached? The refinement algorithm will keep the number of cells fixed but will minimize the maximum cost over all the cells. This can be used for example to run a constant-size simulation where the cells are optimally distributed across the simulation domain. This would be done by setting maxcells to the desired number and cmax to zero.

Following this we can modify our simulation file:

4 3 GfsSimulation GfsBox GfsGEdge {} {
  GfsTime { end = 9 }
  GfsRefine 7
  GtsSurfaceFile half-cylinder.gts
  GfsInit {} { U = 1 }
#  GfsOutputBoundaries {} boundaries
  GfsAdaptVorticity { istep = 1 } { maxlevel = 7 cmax = 1e-2 }
  GfsOutputTime { step = 0.02 } stdout
  GfsOutputBalance { step = 0.02 } stdout
  GfsOutputProjectionStats { step = 0.02 } stdout
  GfsOutputPPM { step = 0.02 } vorticity.ppm {
    min = -100 max = 100 v = Vorticity
  GfsOutputSimulation { step = 0.1 } half-cylinder-%3.1f.gfs {
    variables = U,V,P
  GfsOutputTiming { start = end } stdout
GfsBox { left = GfsBoundaryInflowConstant 1 }
GfsBox {}
GfsBox {}
GfsBox { right = GfsBoundaryOutflow }
1 2 right
2 3 right
3 4 right
We have added two lines and commented out (using #) the line outputting the boundaries (we don't need that anymore, we have the simulation files).

The first line we added says that we want to refine dynamically the mesh through the GfsAdaptVorticity object applied every timestep (istep = 1). The $ \delta$ parameter (cmax) is set to 10-2 .

The second line we added is a new GfsOutput object which displays the ``balance'' of the domain sizes across the different processes (when Gerris is ran in parallel). We will use this to monitor how the number of cells evolves with time as the simulation refines or coarsens the mesh according to our vorticity criterium.

We can now run this new simulation. If the previous simulation did not complete yet, don't be afraid to abort it (Ctrl-C), this one is going to be better (and faster).

% gerris2D half-cylinder.gfs
If we now look at the balance summary written by GfsOutputBalance, we see that initially (step: 0) the total number of cells (on all levels) is 86966, which corresponds to a constant resolution of 4 x 27 x 27 = 512 x 128 . At step 10 the number of cells is down to 990 with a corresponding increase in computational speed. If we now look at the first simulation file we saved, using:
% gfs2oogl2D < half-cylinder-0.1.gfs > boundaries
% gfs2oogl2D -S -z 0 -c Vorticity < half-cylinder-0.1.gfs > squares.oogl
we obtain figure 7 showing not only the domain boundaries as usual, but also the boundaries (thin black lines) between different levels of refinement.
Figure 7: Dynamic adaptive mesh refinement t = 0.1
We see that the mesh is very refined around the solid and around the two vortices developing at the trailing edge and very coarse (one cell per box only) on the downstream part of the domain. If you are not sure what these thin black lines represent, just switch on the edge representation in Geomview (using the Inspect $ \rightarrow$ Appearance menu). You will get a picture looking like figure 8, showing all the cells used for the discretisation.
Figure 8: Dynamic adaptive mesh refinement t = 0.1 . Detail of the cells.
As the simulation goes on, you can see the number of cells in the domain increase as the trailing vortices develop. With a bit of patience you will get to figure 9 showing the fully developed Von Karman vortex street with patches of increased resolution following each vortex. Even when the flow is fully developed using adaptive mesh refinement still saves a factor of ~6 in time and memory use. The advantage of adaptive mesh refinement is even more obvious in situations where it is necessary to use very large domains to avoid any contamination of the solution by the boundary conditions.
Figure 9: Dynamic adaptive mesh refinement t = 9 .
You should also try to animate vorticity.ppm which by now should give you a nice animation of the developing trailing vortices becoming unstable and generating the Von Karman street. If ImageMagick is properly installed on your system you can also try:
% convert vorticity.ppm vorticity.mpg
which will produce a much smaller MPEG video file, suitable for distribution through the network.

next up previous contents
Next: 4 Going further Up: 3 A more complex Previous: 3.3 Visualisation   Contents