Tips and tricks

From Gerris

(Difference between revisions)
Jump to: navigation, search
Revision as of 04:34, 30 January 2008
Rohtav (Talk | contribs)
(Tecplot Snapshots)
← Previous diff
Current revision
Ray (Talk | contribs)
(How to move in 3D in a 2D space ?)
Line 1: Line 1:
== Emacs mode for Gerris files == == Emacs mode for Gerris files ==
-Well, not really but something approaching. Add the following to your <code>.emacs</code>+The default Gerris installation comes with an emacs [http://www.gnu.org/software/emacs/manual/html_node/emacs/Major-Modes.html major mode] to edit Gerris simulation files. To enable it, you need to find where Gerris is installed. For example do:
- (setq auto-mode-alist (cons '("\\.gfs\\'" . shell-script-mode) auto-mode-alist))+ % which gerris2D
 + /usr/bin/gerris2D
-== Generating several movies on-the-fly ==+Then edit your .emacs file and add the following lines:
-While it is fairly simple to use the scripting mode of gfsview and unix pipes to [http://gfs.sourceforge.net/examples/examples/boussinesq.html generate a movie] on the fly from a running simulation, how does one generate several movies simultaneously?+ ;; gerris mode
 + (add-to-list 'load-path "/usr/share/gerris")
 + (require 'gfs-mode)
-Using named unix fifos and the <code>tee</code> utility it is fairly easy too. For example if one has three gfsview files called <code>wide.gfv</code>, <code>closeup.gfv</code> and <code>overview.gfv</code> and want to generate the three corresponding movies <code>wide.mpg</code>, <code>closeup.mpg</code> and <code>overview.mpg</code> in one go, one could use the following script:+This will setup emacs so that opening files with the ".gfs" extension automatically uses the "gfs-mode" major mode. Alternatively you can manually enable this mode on any file using (within emacs):
- #!/bin/sh+ M-x gfs-mode
- +
- movies="wide closeup overview"+
- rm -f $movies+
- mkfifo $movies+
- +
- gerris3D mysimulation.gfs | tee $movies > /dev/null &+
- for movie in $movies; do+
- gfsview-batch3D $movie.gfv < $movie | ppm2mpeg > $movie.mpg &+
- done+
- sleep 10+
- rm -f $movies+
-of course the simulation file <code>mysimulation.gfs</code> should contain lines looking like:+where M- stands for the "emacs Meta key" (usually mapped to "esc" on your keyboard).
- OutputSimulation { istep = 10 } stdout+If you now restart emacs and open a .gfs file, you will see that emacs highlights Gerris keywords. These keywords are also "clickable". Clicking on a keyword will open the corresponding documentation in your web browser. If you want, you can customize the way links are opened by configuring the 'browse-url' emacs function which is used by gfs-mode.
- EventScript { istep = 10 } { echo "Save stdout { width = 1024 height = 768 }" }+
-== Compressing simulation files ==+===Indentation===
-When it is useful to save simulation results at regular intervals, the size of the files can be reduced by using on-the-fly compression. This can be done like this:+The gfs-mode also knows how to properly indent Gerris simulation files. A line can easily be indented by pressing the 'Tab' key with the cursor positioned anywhere on the line.
- OutputSimulation { istep = 100 } sim-%ld.gfs+To automatically indent a block of text, use Ctrl-space (or the mouse) to select the block and type:
- EventScript { istep = 100 } { gzip -f -q sim-*.gfs }+
-GfsView can read compressed GFS files directly.+ M-x indent-region
 +(note that you can use the Tab key to autocomplete command names).
 +===Keyword auto-completion===
-== Writing generic or customized gerris output ==+gfs-mode configures [http://www.gnu.org/software/emacs/manual/html_node/emacs/Dynamic-Abbrevs.html dynamic abbrevs] to provide auto-completion for Gerris keywords. Dynamic abbrevs are usually enabled by default in emacs. They work with any buffer. To use auto-completion, type the beginning of the Gerris keyword (e.g. "Out") and type
 + M-/
-Since most of visualization package on the market (or as open-source) do not support quadtree/octree data format, to use benefits of third party visualization packages rather than GFSView we should convert Gerris results to a general unstructured data. +repeatedly. This will cycle through all the Gerris keywords starting with "Out".
-The function "gfs_write_generic_output" is written for this purpose. +===Comments===
-This function could write output data in [http://www.tecplot.com/ Tecplot] or [http://www.cacr.caltech.edu/~slombey/asci/vtk/vtk_formats.simple.html VTK format]. We mean generic output because Tecplot format is almost a general unstructured data type and you could easily modify its (code or resulted file) for your own need. Alternatively you could contact me and I could do its for you (just describe your format). All data are written in ASCII format.+You can comment out (resp. uncomment) selected blocks of text using
-VTK format could be visualized by powerful free package: [http://www.paraview.org/HTML/Index.html Paraview].+ M-x comment-region (resp. M-x uncomment-region)
-Tecplot format format could be visualized by commersial package+== Generating several movies on-the-fly ==
-[http://www.tecplot.com/ Tecplot]+
-Tecplot format is briefly as follows:+''Note that using [[GfsOutputView]] is simpler and more efficient than the technique described in this section.''
-1) Header (include list of field variables, number of vertexes, number of elements and elements type)+While it is fairly simple to use the scripting mode of gfsview and unix pipes to [http://gfs.sourceforge.net/examples/examples/boussinesq.html generate a movie] on the fly from a running simulation, how does one generate several movies simultaneously?
-2) Spatial coordinates plus filed variables related to listed variables in header+Using named unix fifos it is fairly easy too. For example if one has three gfsview files called <code>wide.gfv</code>, <code>closeup.gfv</code> and <code>overview.gfv</code> and want to generate the three corresponding movies <code>wide.mpg</code>, <code>closeup.mpg</code> and <code>overview.mpg</code> in one go, one could use the following command:
-3) Element connectivity+ % gerris3D mysimulation.gfs | gfsview-batch3D
-For more details regarding to VTK format refer to VTK file format: [http://www.vtk.org/pdf/file-formats.pdf PDF] or [http://www.cacr.caltech.edu/~slombey/asci/vtk/vtk_formats.simple.html HTML]+with a simulation file <code>mysimulation.gfs</code> containing the lines:
 + EventScript { start = 0 } {
 + movies="wide closeup overview"
 + rm -f $movies
 + mkfifo $movies
 + for movie in $movies; do
 + ppm2mpeg < $movie > $movie.mpg &
 + done
 + }
 + OutputSimulation { step = 0.01 } stdout
 + EventScript { step = 0.01 } {
 + movies="wide closeup overview"
 + for movie in $movies; do
 + echo "Clear"
 + cat $movie.gfv
 + echo "Append $movie { width = 1024 height = 768 }"
 + done
 + }
-=== Brief algorithm ===+== Compressing simulation files ==
-1) Construct mesh connectivity 4/8 vertex per cell.+When it is useful to save simulation results at regular intervals, the size of the files can be reduced by using on-the-fly compression. This can be done like this:
-2) Optimize connectivity by removing repeated vertexes, this considerably reduce data size and improve considerably speed of rendering (very important for 3d or large scale simulation).+ OutputSimulation { istep = 100 } sim-%ld.gfs
 + EventScript { istep = 100 } { gzip -f -q sim-*.gfs }
-3) Interpolation of cell-wise data to vertexes.+GfsView can read compressed GFS files directly.
-4) Writing results on a text file.+== Adding objects after a simulation has completed ==
 +If you want to add and execute Gerris objects after a simulation has completed, you can use something like:
-=== Function Prototype ===+ % gerris2D -e "[[GfsOutputScalarStats|OutputScalarStats]] stderr { v = P }" end.gfs > /dev/null
-Prototype of this function is as follows (6 input arguments):+Another example would be computing a distance function (i.e. a levelset function) from a VOF interface description for visualisation using isosurfaces in GfsView. This is easily done using something like:
-<source lang="c">+ % gerris3D -e "[[GfsVariableDistance|VariableDistance]] D T" end.gfs | gfsview3D -s iso.gfv
- void gfs_write_generic_output (gchar *type, +
- gint n_filed, +
- gchar **filed_name, +
- gint plot_depth, +
- GtsBBox * box, +
- GfsDomain * domain);+
-</source>+
 +== Masking out parts of the domain in GfsView ==
-1. type: type of desired output, "VTK" or "TECPLOT" at the moment are supported.+If you want to display only a subset of the mesh in GfsView you can use the <code>NODATA</code> value in a [[GfsFunction]]. For example, setting the scalar field of "Squares" in GfsView to:
-2. n_filed: number of filed data to be visualized, 0 mean just visualization of mesh.+ (P > 0 ? P : NODATA)
-3. filed_name: string array, each item is correspond to each filed, you should include your desired field name, i.e., "U", "V", "W", "P" and desired tracer, e.g., "T", ..., also you could include "Vorticity" to compute and write vorticity field. If for any field, function can not find data pointer from name of variable (except for vorticity), this field is skiped.+will only display the pressure field where it is positive. This trick also works for [[GfsOutputPPM]] and [[GfsOutputGRD]].
-4. plot_depth: depth of cell which are desired to visualize, -1 mean visiting all levels.+== VIM features and Gerris ==
-5. box: if box=NULL, all domain are considered for visualization, but user could define a box to visualize only portion of domain located in the box, to define box you should specify its two corner (by six coordinate, for 2d simulation define two z-value arbitrary).+There are plenty of vim features that can be customized to make your life easier when you work with gerris.
-6. domain: simulation domain.+Using tab, you can autocomplete the keywords of gerris. To do that, add the following lines to your .vimrc:
- + function! Tab_Or_Complete()
-At the moment this function could be used to write data on-the-fly but it is easy to write a short code to convert "*.sim" data (result of Gerris) to described format. Maybe in future somebody write a tool to convert "*.sim" files to desired output, you could write its yourself, for this purpose follow instruction same as "gfs2oogl.c", to load simulation object from "*.sim" file).+ if col('.')>1 && strpart( getline('.'), col('.')-2, 3 ) =~ '^\w'
 + return "\<C-N>"
 + else
 + return "\<Tab>"
 + endif
 + endfunction
 + :inoremap <Tab> <C-R>=Tab_Or_Complete()<CR>
 + if has("autocmd")
 + "set complete+=k/etc/dictionaries-common/words isk+=.,(
 + set complete+=k/usr/share/gerris/gerris.dic isk+=.,(
 + endif " has("autocmd"
-=== How to Use ===+where <code>/usr/share/gerris/gerris.dic</code> may vary according to how gerris was installed on your system.
- +
- +
-To call this function on-the-fly it is sufficient to include this +
-file in "simulation.c" file, i.e. add statement +
- +
-<source lang="c">+
- #include "generic_output.c"+
-</source>+
- +
-after include statements in the "simulation.c" file of your gerris and add its desired function call in time step loop of "simulation_run" function e.g. at a place like this (determined by "######")+
- +
-<source lang="c">+
- +
- static void simulation_run (GfsSimulation * sim)+
- {+
- ...+
- ...+
- ...+
- +
- while (sim->time.t < sim->time.end &&+
- sim->time.i < sim->time.iend) {+
- +
- ...+
- ...+
- ...+
- +
- sim->time.t = sim->tnext;+
- sim->time.i++;+
- +
- gts_range_add_value (&domain->timestep, +
- gfs_clock_elapsed(domain->timer) - tstart);+
- gts_range_update (&domain->timestep);+
- gts_range_add_value (&domain->size, +
- gfs_domain_size (domain, FTT_TRAVERSE_LEAFS, -1));+
- gts_range_update (&domain->size);+
- +
- "######"+
- }+
- +
-</source>+
- +
-assume we like to run "heated.gfs" example +
- +
-use the following statements at the specified place (writing output file after each 100 time steps):+
- +
-<source lang="c">+
- +
- {+
- gint n_field = 5, plot_depth = -1;+
- gchar *filed_name[5] = {"U", "V", "P", "T", "Vorticity"};+
- GtsBBox *box= NULL; +
- +
- if(sim->time.i % 100 == 0) +
- gfs_write_generic_output ("TECPLOT", n_field,+
- filed_name,plot_depth, +
- box, domain);+
- +
- // or +
- +
- if(sim->time.i % 100 == 0) +
- gfs_write_generic_output ("VTK", n_field, +
- filed_name, plot_depth, +
- box, domain);+
- +
- }+
- +
-</source>+
- +
-or if you like to limit visulization within a box use this:+
- +
-<source lang="c">+
- +
- {+
- gint n_field = 5, plot_depth = -1;+
- gchar *filed_name[5] = {"U", "V", "P", "T", "Vorticity"};+
- +
- GtsBBox *box= gts_bbox_new (gts_bbox_class (), NULL, +
- -.5, -.5, -.5,+
- .5, .5, .5);+
- +
- if(sim->time.i % 100 == 0) +
- gfs_write_generic_output ("TECPLOT", n_field, +
- filed_name, plot_depth, +
- box, domain);+
- +
- // or +
- +
- if(sim->time.i % 100 == 0) +
- gfs_write_generic_output ("VTK", n_field, +
- filed_name, plot_depth, +
- box, domain);+
- +
- }+
- +
-</source>+
- +
-note that resulted files are names as:+
- +
- tecplot%i.plt or paraview%i.vtk+
- +
-where %i determine file number, it is automatically incremented after each output+
- +
-Note that you could write various type or a type with various proberties (e.g. writing both hole data and limited data by a box) simultaneously. For this for this purpose you just include suitable call of function in simulation loop.+
- +
-After modification of code, you should recompile your gerris code to activate this added feature. Of course it is better to include such data in input "*.gfs", i.e. treat its like an event, but prefarably this job should be performed by Gerris mantainer, if he detect that this feature is useful for a portion of Gerris users. +
- +
-(in the case of interest, you could announce your interest to this feature in GSF user mailing list).+
- +
- +
-=== Tecplot Snapshots ===+
- +
- +
-[[Image:tecplot1.jpg|middle|thumb|300px]]+
- +
- +
-[[Image:tecplot2.jpg|middle|thumb|300px]]+
- +
-=== Paraview Snapshots ===+
- +
- +
-[[Image:paraview1.jpg|thumb|300px]]+
- +
- +
-[[Image:paraview2.jpg|thumb|300px]]+
- +
- +
-[[Image:paraview3.jpg|thumb|600px]]+
- +
-=== Source Code Listing ===+
- +
-<source lang="c">+
- +
-#include <stdlib.h>+
-#include <math.h>+
-#include <gmodule.h>+
-#include "config.h"+
-#include "gfsconfig.h"+
-#include "simulation.h"+
-#include "output.h"+
-#include "refine.h"+
-#include "solid.h"+
-#include "version.h"+
- +
- +
- +
-#define SQR(a) ((a)*(a))+
- +
-#define my_alloc(pt, n, type, message) {\+
- if(((pt) = (type *) calloc(n,sizeof(type)))==NULL) {\+
- printf("\n\n Assertion: There is not sufficient memory to allocate: %s \n\n",message);\+
- exit(1);\+
- }\+
- }+
- +
- +
-typedef struct {+
- gint node[(FTT_DIMENSION - 1)*4]; +
-} Conectivity;+
- +
- +
-static gint tecplot_counter = 0;+
-static gint vtk_counter = 0;+
- +
-void gfs_write_generic_output (gchar *type, gint n_filed, gchar **filed_name, gint plot_depth, GtsBBox * box, GfsDomain * domain);+
- +
-/********************************************************************+
- +
- count number of leaf cells and compute some spatial data+
- +
-********************************************************************/+
- +
-void count_life_cells (FttCell * cell, gpointer * data)+
-{+
- gint *ncells;+
- gdouble *scene_data, h, x_corner, y_corner, z_corner;+
- FttVector pos;+
- +
- g_return_if_fail (cell != NULL);+
- g_return_if_fail (data != NULL);+
- +
- ncells = data[0];+
- scene_data = data[1];+
- +
- h = ftt_cell_size (cell);+
- +
- scene_data[0] = MIN(h, scene_data[0]);+
- scene_data[1] = MAX(h, scene_data[1]);+
- +
- ftt_cell_pos (cell, &pos);+
- +
- h *= 0.5;+
- +
- x_corner = pos.x - h;+
- scene_data[3] = MIN(x_corner, scene_data[3]);+
- scene_data[6] = MAX(x_corner, scene_data[6]);+
- +
- y_corner = pos.y - h;+
- scene_data[4] = MIN(y_corner, scene_data[4]);+
- scene_data[7] = MAX(y_corner, scene_data[7]);+
- +
-#if FTT_3D+
- +
- z_corner = pos.z - h;+
- scene_data[5] = MIN(z_corner, scene_data[5]);+
- scene_data[8] = MAX(z_corner, scene_data[8]);+
- +
-#endif+
- +
- ++ *ncells;+
-}+
- +
-/********************************************************************+
- +
- fill mesh coordinates+
- +
-********************************************************************/+
- +
-void fill_mesh_coordinate (FttCell * cell, gpointer * data)+
-{+
- gchar *type;+
- gint i, id, *cell_id;+
- FttVector pos, *coordinate;+
- gdouble h, *scene_data, *distance;+
- gdouble scale_y, scale_z, offset_x, offset_y, offset_z;+
- +
- g_return_if_fail (cell != NULL);+
- g_return_if_fail (data != NULL);+
- +
- cell_id = data[0];+
- type = data[1];+
- scene_data = data[2];+
- coordinate = data[3];+
- distance = data[4]; +
- +
- /* scaling and offseting is performed to have only +
- same value of distance for exactly same points */+
- offset_x = scene_data[3];+
- offset_y = scene_data[4];+
- offset_z = scene_data[5];+
- +
- scale_y = scene_data[9];+
- scale_z = scene_data[9] * scene_data[10];+
- +
- id = *cell_id * (FTT_DIMENSION - 1)*4;+
- +
- h = 0.5 * ftt_cell_size (cell);+
- +
- ftt_cell_pos (cell, &pos);+
- +
-#if FTT_2D+
- +
- if(strcmp(type,"TECPLOT") == 0) {+
- coordinate[id+0].x = pos.x + h; coordinate[id+0].y = pos.y + h; +
- coordinate[id+1].x = pos.x - h; coordinate[id+1].y = pos.y + h;+
- coordinate[id+2].x = pos.x - h; coordinate[id+2].y = pos.y - h;+
- coordinate[id+3].x = pos.x + h; coordinate[id+3].y = pos.y - h;+
- }+
- +
- if(strcmp(type,"VTK") == 0) {+
- coordinate[id+0].x = pos.x - h; coordinate[id+0].y = pos.y - h; +
- coordinate[id+1].x = pos.x + h; coordinate[id+1].y = pos.y - h;+
- coordinate[id+2].x = pos.x - h; coordinate[id+2].y = pos.y + h;+
- coordinate[id+3].x = pos.x + h; coordinate[id+3].y = pos.y + h;+
- }+
- +
- for (i = 0; i < (FTT_DIMENSION - 1)*4; i++)+
- distance[id+i] = (coordinate[id+i].x - offset_x) + +
- (coordinate[id+i].y - offset_y) * ( scale_y + 1 ) ;+
- +
-#else /* 3D */+
- +
- if(strcmp(type,"TECPLOT") == 0) {+
- coordinate[id+0].x = pos.x + h; coordinate[id+0].y = pos.y + h; coordinate[id+0].z = pos.z + h; +
- coordinate[id+1].x = pos.x - h; coordinate[id+1].y = pos.y + h; coordinate[id+1].z = pos.z + h; +
- coordinate[id+2].x = pos.x - h; coordinate[id+2].y = pos.y - h; coordinate[id+2].z = pos.z + h; +
- coordinate[id+3].x = pos.x + h; coordinate[id+3].y = pos.y - h; coordinate[id+3].z = pos.z + h; +
- coordinate[id+4].x = pos.x + h; coordinate[id+4].y = pos.y + h; coordinate[id+4].z = pos.z - h; +
- coordinate[id+5].x = pos.x - h; coordinate[id+5].y = pos.y + h; coordinate[id+5].z = pos.z - h; +
- coordinate[id+6].x = pos.x - h; coordinate[id+6].y = pos.y - h; coordinate[id+6].z = pos.z - h; +
- coordinate[id+7].x = pos.x + h; coordinate[id+7].y = pos.y - h; coordinate[id+7].z = pos.z - h;+
- }+
- +
- if(strcmp(type,"VTK") == 0) {+
- coordinate[id+0].x = pos.x - h; coordinate[id+0].y = pos.y - h; coordinate[id+0].z = pos.z - h; +
- coordinate[id+1].x = pos.x + h; coordinate[id+1].y = pos.y - h; coordinate[id+1].z = pos.z - h; +
- coordinate[id+2].x = pos.x - h; coordinate[id+2].y = pos.y + h; coordinate[id+2].z = pos.z - h; +
- coordinate[id+3].x = pos.x + h; coordinate[id+3].y = pos.y + h; coordinate[id+3].z = pos.z - h;+
- coordinate[id+4].x = pos.x - h; coordinate[id+4].y = pos.y - h; coordinate[id+4].z = pos.z + h; +
- coordinate[id+5].x = pos.x + h; coordinate[id+5].y = pos.y - h; coordinate[id+5].z = pos.z + h; +
- coordinate[id+6].x = pos.x - h; coordinate[id+6].y = pos.y + h; coordinate[id+6].z = pos.z + h; +
- coordinate[id+7].x = pos.x + h; coordinate[id+7].y = pos.y + h; coordinate[id+7].z = pos.z + h; +
- }+
- +
- for (i = 0; i < (FTT_DIMENSION - 1)*4; i++)+
- distance[id+i] = (coordinate[id+i].x - offset_x) + +
- (coordinate[id+i].y - offset_y) * ( scale_y + 1 ) + +
- (coordinate[id+i].z - offset_z) * ( scale_z + 1 ) ;+
- +
-#endif /* 3D */+
- +
- ++ *cell_id;+
-}+
- +
-/********************************************************************+
- +
- compute cell vorticity+
- +
-********************************************************************/+
- +
-void compute_cell_vorticity (FttCell * cell, gpointer * data) {+
- +
- GfsVariable **velocity, *vorticity;+
- +
- g_return_if_fail (cell != NULL);+
- g_return_if_fail (data != NULL);+
- +
- velocity = data[0];+
- vorticity = data[1];+
- +
- GFS_VARIABLE (cell, vorticity->i) = gfs_vorticity (cell, velocity);+
- +
-}+
- +
-/********************************************************************+
- +
- interpolate cell data to vertexes and compute interpolation weight+
- +
-********************************************************************/+
- +
-void cell_to_vertex (FttCell * cell, gpointer * data) {+
- +
- static FttDirection d[(FTT_DIMENSION - 1)*4][FTT_DIMENSION] = {+
-#if FTT_2D+
- {FTT_RIGHT, FTT_TOP}, {FTT_LEFT, FTT_TOP}, {FTT_LEFT, FTT_BOTTOM}, {FTT_RIGHT, FTT_BOTTOM}+
-#else /* 3D */+
- {FTT_RIGHT, FTT_TOP, FTT_FRONT}, {FTT_LEFT, FTT_TOP, FTT_FRONT}, +
- {FTT_LEFT, FTT_BOTTOM, FTT_FRONT}, {FTT_RIGHT, FTT_BOTTOM, FTT_FRONT},+
- {FTT_RIGHT, FTT_TOP, FTT_BACK}, {FTT_LEFT, FTT_TOP, FTT_BACK}, +
- {FTT_LEFT, FTT_BOTTOM, FTT_BACK}, {FTT_RIGHT, FTT_BOTTOM, FTT_BACK},+
-#endif /* 3D */+
- };+
- +
- gint i, j, k, n, *cell_id, *n_real_filed;+
- gdouble **field_data;+
- gint id_vertex_cell[(FTT_DIMENSION - 1)*4];+
- +
- Conectivity *conectivity;+
- +
- GfsVariable *variable;+
- +
- g_return_if_fail (cell != NULL);+
- g_return_if_fail (data != NULL);+
- +
- cell_id = data[0];+
- n_real_filed = data[1];+
- conectivity = data[2]; +
- field_data = data[4];+
- k = 5;+
- +
- n = *n_real_filed;+
- +
- for (i = 0; i < (FTT_DIMENSION - 1)*4; i++) +
- id_vertex_cell[i] = conectivity[*cell_id].node[i];+
- +
- for (j = 0; j < n; j++) {+
- variable = data[j+k];+
- for (i = 0; i < (FTT_DIMENSION - 1)*4; i++)+
- field_data[j][id_vertex_cell[i]] = gfs_cell_corner_value (cell, d[i], variable, -1);+
- }+
- +
- ++ *cell_id;+
-}+
- +
-/****************************************************+
- +
- quick sort of double array "a" with length "n" +
- and return integer permutation array "t" +
- +
-*****************************************************/+
- +
-void my_fqsort (gint n, gdouble *a, gint *t) {+
- +
- static gint i, j, k, l, r, s, stackl[100], stackr[100], ww, flg1, flg2;+
- static gdouble w, x;+
- +
- --a; --t;+
- s = 1; stackl[1] = 1; stackr[1] = n;+
- +
-L10: l = stackl[s]; r = stackr[s--];+
- +
-L20: i = l; j = r; k = ( l + r ) / 2; x = a[k];+
- +
- /* repeat until i > j */+
- flg1=1;+
- while (flg1) {+
- flg2=1;+
- while(flg2)+
- if (a[i]<x) ++i; /* search from lower end */ +
- else flg2=0;+
- flg2=1;+
- while (flg2) +
- if (x<a[j]) --j; /* search from upper end */ +
- else flg2=0;+
- +
- if (i<=j) { // swap positions i, j+
- w = a[i]; ww = t[i]; a[i] = a[j]; t[i++] = t[j]; a[j] = w; t[j--] = ww;+
- if (i>j) flg1=0;+
- } else flg1=0;+
- };+
- +
- if ( j - l >= r - i ) {+
- if (l < j ) {+
- stackl[++s] = l; stackr[s] = j;+
- }+
- l = i;+
- } else {+
- if ( i < r ) {+
- stackl[++s]= i; stackr[s] = r;+
- }+
- r = j;+
- }+
- if (l < r) goto L20; +
- if ( s != 0 ) goto L10;+
- return;+
-}+
- +
-/****************************************************+
- +
- quick sort of integer array "a" with length "n" +
- and return integer permutation array "t" +
- +
-****************************************************/+
- +
-void my_iqsort (gint n, gint *a, gint *t) {+
- +
- static gint i, j, k, l, r, s, stackl[100], stackr[100], ww, flg1, flg2;+
- static gint w, x;+
- +
- --a; --t;+
- s = 1; stackl[1] = 1; stackr[1] = n;+
- +
-L10: l = stackl[s]; r = stackr[s--];+
- +
-L20: i = l; j = r; k = ( l + r ) / 2; x = a[k];+
- +
- /* repeat until i > j */+
- flg1=1;+
- while (flg1) {+
- flg2=1;+
- while(flg2)+
- if (a[i]<x) ++i; /* search from lower end */ +
- else flg2=0;+
- flg2=1;+
- while (flg2) +
- if (x<a[j]) --j; /* search from upper end */ +
- else flg2=0;+
- +
- if (i<=j) { // swap positions i, j+
- w = a[i]; ww = t[i]; a[i] = a[j]; t[i++] = t[j]; a[j] = w; t[j--] = ww;+
- if (i>j) flg1=0;+
- } else flg1=0;+
- };+
- +
- if ( j - l >= r - i ) {+
- if (l < j ) {+
- stackl[++s] = l; stackr[s] = j;+
- }+
- l = i;+
- } else {+
- if ( i < r ) {+
- stackl[++s]= i; stackr[s] = r;+
- }+
- r = j;+
- }+
- if (l < r) goto L20; +
- if ( s != 0 ) goto L10;+
- return;+
-}+
- +
-/********************************************************************+
- +
- remove repeated vertex from conectivity and re-new conectivity+
- +
-********************************************************************/+
- +
-void remove_repeated_vertex ( gint n_cells, +
- gint *n_vertex,+
- gdouble *scene_data,+
- gdouble *distance,+
- FttVector *coordinate,+
- Conectivity *conectivity+
- ) {+
- +
- gint i, j, n, n_vertexes_new, nv, *perm, *iperm, *label, *id_vertex;+
- gdouble small = 0.25 * scene_data[0];+
- gchar *flg;+
- FttVector *new_coordinate;+
- +
- nv = 4*(FTT_DIMENSION - 1);+
- n = nv*n_cells;+
- +
- /* allocate array for permutation array and label of repeated vertexes */+
- +
- my_alloc(perm, n, gint, "temporary permutation array");+
- my_alloc(iperm, n, gint, "temporary ipermutation array");+
- my_alloc(label, n, gint, "temporary label array");+
- my_alloc(flg, n, gchar, "temporary flag array");+
- +
- for (i=0; i<n; i++) {+
- perm [i] = label [i]= i;+
- flg [i] = 1;+
- }+
- +
- /* +
- quick sort algorithm, note that qsort of ANCI C is not suitable+
- for our purpose, since we need permutation array rather than sorted list +
- */+
- +
- my_fqsort (n, distance, perm);+
- +
- n_vertexes_new = 1;+
- +
- for (i=0; i<(n-1); i++) {+
- +
- if( fabs( distance[i] - distance[i+1] ) < small ) {+
- +
- gint iv1 = perm[i];+
- gint iv2 = perm[i+1];+
- +
- if ( +
- fabs( coordinate[iv1].x - +
- coordinate[iv2].x ) < small +
- &&+
- fabs( coordinate[iv1].y - +
- coordinate[iv2].y ) < small+
-#if FTT_3D+
- &&+
- fabs( coordinate[iv1].z - +
- coordinate[iv2].z ) < small+
-#endif+
- ) {+
- label[iv2] = label[iv1];+
- flg[iv2] = 0;+
- } else ++n_vertexes_new;+
- } else ++n_vertexes_new;+
- }+
- +
- my_alloc(id_vertex, n_vertexes_new, gint, "temporary id vertex array");+
- +
- n_vertexes_new = 0;+
- +
- for (i=0; i<n; i++) if( flg [i] > 0 ) {+
- perm [n_vertexes_new] = i;+
- id_vertex [n_vertexes_new++] = label[i];+
- }+
- +
- my_iqsort (n_vertexes_new, id_vertex, perm);+
- +
- /* re-new conectivity */+
- +
- for (i=0; i<n_vertexes_new; i++) iperm [perm[i]] = i;+
- +
- for (i=0; i<n_cells; i++)+
- for (j=0; j<nv; j++)+
- conectivity[i].node[j] = iperm [label[i*nv + j]];+
- +
- g_free(label);+
- g_free(perm);+
- g_free(iperm);+
- g_free(id_vertex);+
- +
- my_alloc(new_coordinate, n_vertexes_new, FttVector, "temporary array for new_coordinate");+
- +
- j=0;+
- for(i=0; i<n; i++) if(flg[i]>0) {+
- new_coordinate[j].x = coordinate[i].x;+
- new_coordinate[j].y = coordinate[i].y;+
-#if FTT_3D+
- new_coordinate[j].z = coordinate[i].z;+
-#endif+
- ++j;+
- }+
- +
- *n_vertex = n_vertexes_new;+
- +
- for(i=0; i<n_vertexes_new; i++) {+
- coordinate[i].x = new_coordinate[i].x;+
- coordinate[i].y = new_coordinate[i].y;+
-#if FTT_3D+
- coordinate[i].z = new_coordinate[i].z;+
-#endif+
- }+
- +
- g_free(flg);+
- g_free(new_coordinate);+
-}+
- +
-/********************************************************************+
- +
- write a tecplot file+
- +
-*********************************************************************/+
- +
-void gfs_write_generic_output (gchar *type, gint n_filed, gchar **filed_name, gint plot_depth, GtsBBox * box, GfsDomain * domain)+
-{+
- gchar file_name[100], temp_string[60];+
- FILE * fp;+
- gint n, n_cells, n_vertex, cell_id, i, j, k, n_real_filed, id_real_filed[20], vort, have_vel_vector[4];+
- gdouble *distance, scene_data[20], *field_data[20];+
- +
- GfsVariable *vorticity;+
- +
- Conectivity *conectivity;+
- FttVector *coordinate;+
- +
- gpointer data[20];+
- +
- /* count number of leaf cells and some useful spatial data */+
- +
- n_cells = 0;+
- +
- scene_data[0] = G_MAXDOUBLE; /* h_min */+
- scene_data[1] = 0.0; /* h_max */+
- +
- scene_data[3] = G_MAXDOUBLE; /* x_min */+
- scene_data[4] = G_MAXDOUBLE; /* y_min */+
- scene_data[5] = G_MAXDOUBLE; /* z_min */+
- +
- scene_data[6] = G_MINDOUBLE; /* x_max */+
- scene_data[7] = G_MINDOUBLE; /* y_max */+
- scene_data[8] = G_MINDOUBLE; /* z_max */+
- +
- data[0] = &n_cells; +
- data[1] = scene_data;+
- +
- if(box != NULL) +
- gfs_domain_cell_traverse_box (domain, box, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, plot_depth,+
- (FttCellTraverseFunc) count_life_cells, data);+
- else +
- gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, plot_depth,+
- (FttCellTraverseFunc) count_life_cells, data);+
- +
- scene_data[9] = scene_data[6] - scene_data[3]; /* lenght_x */+
- scene_data[10] = scene_data[7] - scene_data[4]; /* lenght_y */+
- scene_data[11] = scene_data[8] - scene_data[5]; /* lenght_z */+
- +
- +
- /* allocate temporary memory to save mesh coordinates and its conectivity */+
- +
- n = n_cells * 4*(FTT_DIMENSION - 1);+
- +
- my_alloc(conectivity, n_cells, Conectivity, "temporary mesh conectivity");+
- my_alloc(coordinate, n, FttVector, "temporary mesh coordinates");+
- my_alloc(distance, n, gdouble, "temporary distance of mesh nodes");+
- +
- /* fill coordinates of leaf cells and their distance from a reference point */+
- /* note that in this stage connectivity is priori known */+
- +
- cell_id = 0;+
- +
- data[0] = &cell_id;+
- data[1] = type;+
- data[2] = scene_data;+
- data[3] = coordinate; +
- data[4] = distance;+
- +
- if(box != NULL) +
- gfs_domain_cell_traverse_box (domain, box, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, plot_depth,+
- (FttCellTraverseFunc) fill_mesh_coordinate, data);+
- else +
- gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, plot_depth,+
- (FttCellTraverseFunc) fill_mesh_coordinate, data);+
- +
- /* +
- removing repeated nodes from global connectivity +
- for this: we sort nodes based on the distance from+
- a refrence ponit and then detect repeated ones, remove +
- and renew connectivity. Note that to be compatible with formal +
- post-processors, like Tecplot, VTK, Ensight, we have to convert +
- cell-based data to vertex-based data and use weighted interpolation+
- of cell-wise data to compute data at vertexes. Also precence of +
- repeated vertexes decrease rendering speed of post-processor +
- */+
- +
- remove_repeated_vertex ( n_cells, &n_vertex, scene_data, distance, coordinate, conectivity);+
- +
- /* allocate memory to store vertex interpolation of cell-wise data */+
- k = 5; +
- vort = -1;+
- have_vel_vector[0] = 0;+
- n_real_filed = 0;+
- for(i=0; i<n_filed; i++) {+
- data[k+n_real_filed] = gfs_variable_from_name (domain->variables, filed_name[i]);+
- if(data[k+n_real_filed] == NULL) {+
- if(strcmp( filed_name[i], "Vorticity") == 0 ) vort = i;+
- }+
- else {+
- if(strcmp( filed_name[i], "U") == 0 ) {++have_vel_vector[0]; have_vel_vector[1] = i;}+
- if(strcmp( filed_name[i], "V") == 0 ) {++have_vel_vector[0]; have_vel_vector[2] = i;}+
- if(strcmp( filed_name[i], "W") == 0 ) {++have_vel_vector[0]; have_vel_vector[3] = i;}+
- +
- my_alloc(field_data[n_real_filed], n_vertex, gdouble, "temporary array to store vertex interpolated data");+
- id_real_filed[n_real_filed++] = i;+
- }+
- }+
- +
- if(vort>0) {+
- +
- vorticity = gfs_temporary_variable (domain);+
- +
- data[0] = gfs_domain_velocity (domain);+
- data[1] = vorticity;+
- +
- if(box != NULL) +
- gfs_domain_cell_traverse_box (domain, box, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, (plot_depth == -1 ? (-1):(plot_depth+1)), +
- (FttCellTraverseFunc) compute_cell_vorticity, data);+
- else +
- gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, (plot_depth == -1 ? (-1):(plot_depth+1)),+
- (FttCellTraverseFunc) compute_cell_vorticity, data);+
- +
- data[k+n_real_filed] = vorticity;+
- +
- my_alloc(field_data[n_real_filed], n_vertex, gdouble, "temporary array to store vertex interpolated data");+
- +
- id_real_filed[n_real_filed++] = vort;+
- }+
- +
- if(n_real_filed>0) {+
- for(j=0; j<n_real_filed; j++) +
- for(i=0; i<n_vertex; i++) +
- field_data[j][i] = 0.;+
- cell_id = 0;+
- +
- data[0] = &cell_id;+
- data[1] = &n_real_filed;+
- data[2] = conectivity; +
- data[4] = field_data;+
- +
- if(box != NULL) +
- gfs_domain_cell_traverse_box (domain, box, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, plot_depth, +
- (FttCellTraverseFunc) cell_to_vertex, data);+
- else +
- gfs_domain_cell_traverse (domain, FTT_PRE_ORDER, FTT_TRAVERSE_LEAFS, plot_depth,+
- (FttCellTraverseFunc) cell_to_vertex, data);+
- +
- }+
- +
- +
- /* write data on a file */+
- +
- /* ---------------------------------------------------+
- Tecplot output, it could be used as a general output+
- you could easily modify sort of data based on your+
- specific data format+
- ------------------------------------------------- */+
- +
- if(strcmp(type,"TECPLOT") == 0) {+
- +
- int nn = (FTT_DIMENSION - 1)*4;+
- +
- /* name of output file */+
- +
- ++ tecplot_counter;+
- +
- strcpy(file_name,"tecplot");+
- sprintf(temp_string,"%i", tecplot_counter); +
- strcat(temp_string, ".plt");+
- strcat(file_name,temp_string);+
- +
- fp = fopen (file_name, "w");+
- +
- /* write header of vtk file */+
- +
- fprintf (fp, " TITLE = 'Gerris Flow Solver Output' \n");+
- +
-#if FTT_2D+
- fprintf (fp, " VARIABLES = 'X', 'Y'");+
-#else /* 3D */+
- fprintf (fp, " VARIABLES = 'X', 'Y', 'Z'");+
-#endif+
- +
- for(j=0; j<n_real_filed; j++) +
- fprintf (fp, ", '%s' ",filed_name[id_real_filed[j]]);+
- fprintf (fp, " \n");+
-#if FTT_2D+Another interesting feature is to open in firefox the syntax reference page of the command under the cursor. To do that, you have to add to your .vimrc:
- fprintf (fp, " ZONE N=%i, E=%i, F=FEPOINT, ET=QUADRILATERAL \n", n_vertex, n_cells);+
-#else /* 3D */+
- fprintf (fp, " ZONE N=%i, E=%i, F=FEPOINT, ET=BRICK \n", n_vertex, n_cells);+
-#endif+
- + function! OnlineDoc()
- /* write vertex coordinates and related field data */+ let s:browser = "firefox"
 + let s:wordUnderCursor = expand("<cword>")
 + if s:wordUnderCursor =~ 'Gfs'
 + let s:url = "http://gfs.sourceforge.net/wiki/index.php/".s:wordUnderCursor
 + else
 + let s:url = "http://gfs.sourceforge.net/wiki/index.php/Gfs".s:wordUnderCursor
 + endif
 + let s:cmd = "silent !" . s:browser . " " . s:url
 + "echo s:cmd
 + execute s:cmd
 + redraw!
 + endfunction
- for(i=0; i<n_vertex; i++) {+and then to map this function with some sequence of letters (in my case y use \w )
-#if FTT_2D+ map <Leader>w :call OnlineDoc()<CR>
- fprintf (fp, " %f %f ", coordinate[i].x, coordinate[i].y);+
-#else /* 3D */+
- fprintf (fp, " %f %f %f ", coordinate[i].x, coordinate[i].y, coordinate[i].z);+
-#endif+
- for(j=0; j<n_real_filed; j++) +
- fprintf (fp, " %f ",field_data[j][i]);+
- fprintf (fp, " \n");+
- }+
- /* write conectivity data */+Remember that you can create your own vim plugin for gfs files as
- for(i=0; i<n_cells; i++) {+ au! BufNewFile,BufRead *.gfs set filetype=gfs
- for (j = 0; j < nn; j++)+
- fprintf (fp, " %i ", conectivity[i].node[j]+1);+
- fprintf (fp, " \n");+
- }+
- fclose(fp);+
- fprintf (stderr, "\n %s file %i is dumped n_cells = %i \n", type, tecplot_counter, n_cells);+
- }+
 +and to place all this stuff in .vim/after/ftplugin/gfs.vim (you can also to load the color scheme of c.vim!!!)
- /* ---------------------------------------------------+== Restarting parallel simulations ==
- VTK output, it could be imported and visualized by +
- free visualization package ParaView, due to +
- consistancy with various platform ASCII format is used+
- --------------------------------------------------- */+
- if(strcmp(type,"VTK") == 0) {+Starting from version 2010-07-10 it should be possible to restart large parallel simulations. Just add the standard line:
- int nn = (FTT_DIMENSION - 1)*4;+ [[GfsOutputSimulation|OutputSimulation]] { step = 1 } sim-%g.gfs
- /* name of output file */+and do e.g.
- +
- ++ vtk_counter;+
- +
- strcpy(file_name,"paraview"); +
- sprintf(temp_string,"%i", vtk_counter); +
- strcat(temp_string, ".vtk");+
- strcat(file_name,temp_string);+
- +
- fp = fopen (file_name, "w");+
 + mpirun -np 8 gerris3D sim-10.gfs
- /* write header of vtk file */+with the same number of processors you initially used to create <code>sim-10.gfs</code>.
- fprintf (fp, "# vtk DataFile Version 2.0\n");+Note that if you use this (simple) method each processor will need to read the entire simulation file. If you notice that this phase takes a significant amount of time (e.g. due to filesystem performance when being accessed simultaneously by many processes), you can optimise things using:
- fprintf (fp, "Unstructured Grid Example\n");+
- fprintf (fp, "ASCII\n");+
- fprintf (fp, "DATASET UNSTRUCTURED_GRID\n");+
- fprintf (fp, "\n");+
- /* write spatial coordinates */+ [[GfsOutputSimulation|OutputSimulation]] { step = 1 } sim-%d-%g.gfs
- fprintf (fp, "POINTS %i float\n", n_vertex);+which will generate a separate file for each process. You can then restart the simulation using e.g.
- +
-#if FTT_2D+
- for(i=0; i<n_vertex; i++) +
- fprintf (fp, "%f %f 0 \n", coordinate[i].x, coordinate[i].y);+
-#else /* 3D */+
- for(i=0; i<n_vertex; i++)+
- fprintf (fp, "%f %f %f\n", coordinate[i].x, coordinate[i].y, coordinate[i].z);+
-#endif+
- fprintf (fp, "\n");+
- /* write conectivity data */+ mpirun -np 8 gerris3D sim-%d-10.gfs
- fprintf (fp, "CELLS %i %i \n", n_cells, n_cells*(nn+1));+Note that this requires a version of Gerris later than 2010-07-20.
- for(i=0; i<n_cells; i++) {+== Parallel compilation ==
- fprintf (fp, " %i ",nn);+
- for (j = 0; j < nn; j++)+
- fprintf (fp, " %i ", conectivity[i].node[j]);+
- fprintf (fp, " \n");+
- }+
- fprintf (fp, "\n");+
- fprintf (fp, "CELL_TYPES %i \n",n_cells);+If you have a multi-core system (let's say 4), you can compile Gerris (and other sources) 4 times faster using:
- for(i=0; i<n_cells; i++)+ % make -j4
-#if FTT_2D+
- fprintf (fp, "8 \n"); +
-#else+
- fprintf (fp, "11 \n"); +
-#endif+
- fprintf (fp, "\n");+
- /* write scaler fields */+== Visualizing VTK parallel files using Visit ==
- if(n_real_filed>0) {+[https://wci.llnl.gov/codes/visit/ Visit] is a nice visualization tool that allows us to visualize, among many other formats, VTK files. For parallel runs, one can generate one vtk file per processor. Then, in order to visualize all the simulation domain, we need to generate a text file with extension .visit that gathers the names of the different vtk files at a given timestep. For example, if we have generated file-0-0.vtk and file-0-1.vtk at t=0, we just need to create the following text file:
- fprintf (fp, "POINT_DATA %i\n",n_vertex);+ !NBLOCKS 2
- fprintf (fp, "\n");+ file-0-0.vtk
 + file-0-1.vtk
- for(j=0; j<n_real_filed; j++) {+If you have solutions at different time steps and you want to load them at once, then you will need to do something as follows:
- +
- fprintf (fp, "SCALARS %s float\n",filed_name[id_real_filed[j]]);+
- fprintf (fp, "LOOKUP_TABLE default \n");+
- for(i=0; i<n_vertex; i++) + !NBLOCKS 2
- fprintf(fp, "%f\n",field_data[j][i]);+ file-0-0.vtk
- fprintf (fp, "\n");+ file-0-1.vtk
- }+ file-1-0.vtk
- }+ file-1-1.vtk
- /* write velosity fields as a vector field, +Further information can be found [https://wci.llnl.gov/codes/visit/2.0.0/GettingDataIntoVisIt2.0.0.pdf here].
- if user include all velosity components */+
- if(have_vel_vector[0] == FTT_DIMENSION) {+== How to move in 3D in a 2D space ? ==
- fprintf (fp, "VECTORS Velocity_vector float\n");+If you tried to visualize a Z scalar field in a linear object of GfsView2D you noticed it's difficult to move in 3D
 +because you are in a plane ! ;-)
-#if FTT_2D+May be you tried to modify the q0,q1,q2,q3 parameters of your gfv file by hand, but if you don't like sine or cosine the solution
- for(i=0; i<n_vertex; i++)+is to use the shift key of the keyboard when you move the mouse.
- fprintf (fp, "%f %f 0\n",field_data[have_vel_vector[1]][i], field_data[have_vel_vector[2]][i]);+
-#else+
- for(i=0; i<n_vertex; i++)+
- fprintf (fp, "%f %f %f\n",field_data[have_vel_vector[1]][i], field_data[have_vel_vector[2]][i], field_data[have_vel_vector[3]][i]);+
-#endif+
- }+
- fclose(fp); +
- fprintf (stderr, "\n %s file %i is dumped n_cells = %i \n", type, vtk_counter, n_cells);+
- }+
- g_free(conectivity);+you can see what happens [http://gfs.sourceforge.net/wiki/images/c/c6/Visu_2D3D_before.jpg before] and [http://gfs.sourceforge.net/wiki/images/8/8f/Visu_2D3D_after.jpg after] the magic shift
- g_free(coordinate);+
- g_free(distance);+
- for(j=0; j<n_real_filed; j++) g_free(field_data[j]);+
-}+
-</source>+or look in the [http://gerris.dalembert.upmc.fr/gerris/examples/examples/monai.html#htoc25 Tsunami runup onto a complex three-dimensional beach] example.

Current revision

Contents

Emacs mode for Gerris files

The default Gerris installation comes with an emacs major mode to edit Gerris simulation files. To enable it, you need to find where Gerris is installed. For example do:

% which gerris2D
/usr/bin/gerris2D

Then edit your .emacs file and add the following lines:

;; gerris mode
(add-to-list 'load-path "/usr/share/gerris")
(require 'gfs-mode)

This will setup emacs so that opening files with the ".gfs" extension automatically uses the "gfs-mode" major mode. Alternatively you can manually enable this mode on any file using (within emacs):

M-x gfs-mode

where M- stands for the "emacs Meta key" (usually mapped to "esc" on your keyboard).

If you now restart emacs and open a .gfs file, you will see that emacs highlights Gerris keywords. These keywords are also "clickable". Clicking on a keyword will open the corresponding documentation in your web browser. If you want, you can customize the way links are opened by configuring the 'browse-url' emacs function which is used by gfs-mode.

Indentation

The gfs-mode also knows how to properly indent Gerris simulation files. A line can easily be indented by pressing the 'Tab' key with the cursor positioned anywhere on the line.

To automatically indent a block of text, use Ctrl-space (or the mouse) to select the block and type:

M-x indent-region

(note that you can use the Tab key to autocomplete command names).

Keyword auto-completion

gfs-mode configures dynamic abbrevs to provide auto-completion for Gerris keywords. Dynamic abbrevs are usually enabled by default in emacs. They work with any buffer. To use auto-completion, type the beginning of the Gerris keyword (e.g. "Out") and type

M-/

repeatedly. This will cycle through all the Gerris keywords starting with "Out".

Comments

You can comment out (resp. uncomment) selected blocks of text using

M-x comment-region (resp. M-x uncomment-region)

Generating several movies on-the-fly

Note that using GfsOutputView is simpler and more efficient than the technique described in this section.

While it is fairly simple to use the scripting mode of gfsview and unix pipes to generate a movie on the fly from a running simulation, how does one generate several movies simultaneously?

Using named unix fifos it is fairly easy too. For example if one has three gfsview files called wide.gfv, closeup.gfv and overview.gfv and want to generate the three corresponding movies wide.mpg, closeup.mpg and overview.mpg in one go, one could use the following command:

% gerris3D mysimulation.gfs | gfsview-batch3D

with a simulation file mysimulation.gfs containing the lines:

EventScript { start = 0 } {
  movies="wide closeup overview"
  rm -f $movies
  mkfifo $movies
  for movie in $movies; do
    ppm2mpeg < $movie > $movie.mpg &
  done
}
OutputSimulation { step = 0.01 } stdout
EventScript { step = 0.01 } {
  movies="wide closeup overview"
  for movie in $movies; do
    echo "Clear"
    cat $movie.gfv
    echo "Append $movie { width = 1024 height = 768 }"
  done
}

Compressing simulation files

When it is useful to save simulation results at regular intervals, the size of the files can be reduced by using on-the-fly compression. This can be done like this:

OutputSimulation { istep = 100 } sim-%ld.gfs
EventScript { istep = 100 } { gzip -f -q sim-*.gfs }

GfsView can read compressed GFS files directly.

Adding objects after a simulation has completed

If you want to add and execute Gerris objects after a simulation has completed, you can use something like:

% gerris2D -e "OutputScalarStats stderr { v = P }" end.gfs > /dev/null

Another example would be computing a distance function (i.e. a levelset function) from a VOF interface description for visualisation using isosurfaces in GfsView. This is easily done using something like:

% gerris3D -e "VariableDistance D T" end.gfs | gfsview3D -s iso.gfv

Masking out parts of the domain in GfsView

If you want to display only a subset of the mesh in GfsView you can use the NODATA value in a GfsFunction. For example, setting the scalar field of "Squares" in GfsView to:

(P > 0 ? P : NODATA)

will only display the pressure field where it is positive. This trick also works for GfsOutputPPM and GfsOutputGRD.

VIM features and Gerris

There are plenty of vim features that can be customized to make your life easier when you work with gerris.

Using tab, you can autocomplete the keywords of gerris. To do that, add the following lines to your .vimrc:

function! Tab_Or_Complete()
 if col('.')>1 && strpart( getline('.'), col('.')-2, 3 ) =~ '^\w'
   return "\<C-N>"
 else
   return "\<Tab>"
 endif
endfunction
:inoremap <Tab> <C-R>=Tab_Or_Complete()<CR>

if has("autocmd")
"set complete+=k/etc/dictionaries-common/words isk+=.,(
set complete+=k/usr/share/gerris/gerris.dic isk+=.,(
endif " has("autocmd"

where /usr/share/gerris/gerris.dic may vary according to how gerris was installed on your system.

Another interesting feature is to open in firefox the syntax reference page of the command under the cursor. To do that, you have to add to your .vimrc:

function! OnlineDoc()
       let s:browser = "firefox"
       let s:wordUnderCursor = expand("<cword>")
        if s:wordUnderCursor =~ 'Gfs'
               let s:url = "http://gfs.sourceforge.net/wiki/index.php/".s:wordUnderCursor
       else
               let s:url = "http://gfs.sourceforge.net/wiki/index.php/Gfs".s:wordUnderCursor
       endif
       let s:cmd = "silent !" . s:browser . " " . s:url
       "echo  s:cmd
       execute  s:cmd
       redraw!
endfunction

and then to map this function with some sequence of letters (in my case y use \w )

map <Leader>w :call OnlineDoc()<CR>

Remember that you can create your own vim plugin for gfs files as

au! BufNewFile,BufRead *.gfs set filetype=gfs

and to place all this stuff in .vim/after/ftplugin/gfs.vim (you can also to load the color scheme of c.vim!!!)

Restarting parallel simulations

Starting from version 2010-07-10 it should be possible to restart large parallel simulations. Just add the standard line:

OutputSimulation { step = 1 } sim-%g.gfs

and do e.g.

mpirun -np 8 gerris3D sim-10.gfs

with the same number of processors you initially used to create sim-10.gfs.

Note that if you use this (simple) method each processor will need to read the entire simulation file. If you notice that this phase takes a significant amount of time (e.g. due to filesystem performance when being accessed simultaneously by many processes), you can optimise things using:

OutputSimulation { step = 1 } sim-%d-%g.gfs

which will generate a separate file for each process. You can then restart the simulation using e.g.

mpirun -np 8 gerris3D sim-%d-10.gfs

Note that this requires a version of Gerris later than 2010-07-20.

Parallel compilation

If you have a multi-core system (let's say 4), you can compile Gerris (and other sources) 4 times faster using:

% make -j4

Visualizing VTK parallel files using Visit

Visit is a nice visualization tool that allows us to visualize, among many other formats, VTK files. For parallel runs, one can generate one vtk file per processor. Then, in order to visualize all the simulation domain, we need to generate a text file with extension .visit that gathers the names of the different vtk files at a given timestep. For example, if we have generated file-0-0.vtk and file-0-1.vtk at t=0, we just need to create the following text file:

 !NBLOCKS 2
 file-0-0.vtk
 file-0-1.vtk

If you have solutions at different time steps and you want to load them at once, then you will need to do something as follows:

 !NBLOCKS 2
 file-0-0.vtk
 file-0-1.vtk
 file-1-0.vtk
 file-1-1.vtk

Further information can be found here.

How to move in 3D in a 2D space ?

If you tried to visualize a Z scalar field in a linear object of GfsView2D you noticed it's difficult to move in 3D because you are in a plane ! ;-)

May be you tried to modify the q0,q1,q2,q3 parameters of your gfv file by hand, but if you don't like sine or cosine the solution is to use the shift key of the keyboard when you move the mouse.

you can see what happens before and after the magic shift

or look in the Tsunami runup onto a complex three-dimensional beach example.

Personal tools
communication