An Example of Use of the Gerris Object System

From Gerris

Revision as of 18:36, 8 November 2007; view current revision
←Older revision | Newer revision→
Jump to: navigation, search


A new class to initialize periodic velocity fields

In this session, we are going to attempt to create a new class that will initialize the velocity field in the form of an array of vortices. While it is possible to do this with a little C code section in the simulation file, we will do it here as an educational exercise.

A template for new object classes

GTS comes with a simple script called gtstemplate which will generate a C code template for the new class we want to create. For a summary of its syntax just type:

% gtstemplate
Usage: gtstemplate [OPTIONS] Class ParentClass

As an example we are going to try and create a template for a new initial condition class called InitPeriodic. Just type:

% gtstemplate --overload=read --overload=write --overload=event \
InitPeriodic GfsGenericInit > init_periodic.c

We have just created a template for a new class called InitPeriodic derived from GfsGenericInit and where the read, write and event methods are overloaded. Fire up your favourite editor and have a look at the file generated: init_periodic.c. The first thing you see is that the file is divided in two sections starting with

/* InitPeriodic: Header */


/* InitPeriodic: Object */

As in the previous examples these sections correspond respectively to the declaration of the structures and functions (header part) and to the corresponding definitions (object part). If we were to use these functions and structures in a library, these two parts would be in separate files (a .h and a corresponding .c file). For what we are interested in they are fine staying in the same file. There is a function we will not need: init_periodic_new, just remove the lines declaring and defining it in the header and object sections. We are left with only one function: init_periodic_class. This function essentially registers our new class and its associated attributes, which are things like:

  • the name of the class: "InitPeriodic",
  • the sizes of the object and class structures: sizeof (InitPeriodic) and sizeof (InitPeriodicClass),
  • function to call to initialize the class: init_periodic_class_init,
  • function to call to initialize a new object: init_periodic_init,

If we now look at what the init_periodic_class function does, we see that it indeed overloads the event, read and write methods of our class with locally defined functions. The GFS_EVENT_CLASS and GTS_OBJECT_CLASS macro calls are casting operators which change the type of our new class (InitPeriodicClass) to the types of one of its parent class. GtsObjectClass is the ancestor of all object classes and we see that the read and write methods are thus defined for all objects. The parent class GfsGenericInitClass is also an alias of the GfsEvent class (as we have seen before it is a special type of event occurring once at the start of the simulation) which has an associated event method.

The goal: initial conditions

With this preliminary overall understanding of what the template does, we are ready to define our new initialization class. What we want to do is initialize the 2D velocity field using the following functions:

u(x,y) = − cos(2mπ x)sin(2mπ y) (1)

v(x,y) = sin(2mπ x)cos(2mπ y) (2)

where m is a parameter. This is an exact stationary solution of the 2D Euler equations with periodic boundary conditions (each vortex is a kind of flywheel).


We need to store the value of m somewhere. A clean way to do that is to add m as a data associated with our InitPeriodic object. To do this just modify the definition of the structure in the header part of the file, like this:

struct _InitPeriodic {
/*< private >*/
GfsGenericInit parent;
/*< public >*/
gdouble m;

Note again that the GfsGenericInit parent; a declaration of member data is the first member in the structure. The whole object-inheritance mechanism in C relies on this data alignment.


We now need to specify the value of m. The first value we want to assign is the default value of m for an object newly created. We can do that using the initialisation function init_periodic_init which is called everytime a new InitPeriodic object is created (it is the constructor function of the InitPeriodicClass).

Init method

Just add:

static void init_periodic_init (InitPeriodic * object)
object->m = 1.;

Our default value for m is now one. What we really want is being able to specify the value of m in the parameter file used by gerris. We can do that by using the read and write methods.

Write method

We will start with the write method: init_periodic_write. This function first calls the write method of the parent class of our object:

/* call write method of parent */
if (GTS_OBJECT_CLASS (init_periodic_class ())->parent_class->write)
(* GTS_OBJECT_CLASS (init_periodic_class ())->parent_class->write)
(o, fp);

The parent class is given by the parent_class field of any GtsObjectClass. The write method might not be defined for the parent class, so it is always a good (safe) idea to first test that it is indeed defined.

Now that we have written the data associated with the parent class, we can write our own data like this:

static void init_periodic_write (GtsObject * o, FILE * fp)
/* call write method of parent */
if (GTS_OBJECT_CLASS (init_periodic_class ())->parent_class->write)
(* GTS_OBJECT_CLASS (init_periodic_class ())->parent_class->write)
(o, fp);
/* do object specific write here */
fprintf (fp, " %g", INIT_PERIODIC (o)->m);

Note that we used a space to separate this new data from the data of the parent class and that we did not add anything after our own data (no space or newline). This is so that we can eventually extend (through inheritance) this class by adding more data in exactly the same way.

Read Method

We now need to define a “symmetrical” read method which will read the value of m from the parameter file. In the init_periodic_read method, we first read the parameters associated with the parent class. The line:

if (fp->type == GTS_ERROR)

just checks if an error occurred while reading parameters for the parent class, in which case the method returns prematurely. To read our parameter we are going to use functions associated with the GtsFile data type (have a look at the GTS reference manual for details):

static void init_periodic_read (GtsObject ** o, GtsFile * fp)
/* call read method of parent */
if (GTS_OBJECT_CLASS (init_periodic_class ())->parent_class->read)
(* GTS_OBJECT_CLASS (init_periodic_class ())->parent_class->read)
(o, fp);
if (fp->type == GTS_ERROR)
if (fp->type != GTS_INT && fp->type != GTS_FLOAT) {
gts_file_error (fp, "expecting a number (m)");
INIT_PERIODIC (*o)->m = atof (fp->token->str);
/* do not forget to prepare for next read */
gts_file_next_token (fp);

We first check that the current token is an integer (GTS_INT) or a floating point number (GTS_FLOAT). If it is not we set an error message describing the problem and return prematurely. If it is, we convert the string describing the current token (fp->token->str) into a floating point number (using the standard C atof function) and assign it to m. Note that we need to use the type conversion macro INIT_PERIODIC because the argument passed to init_periodic_read is a generic GtsObject. As we are dealing with a method of our InitPeriodicClass we know that this GtsObject is also a InitPeriodic object and the type conversion is then legitimate.

Event method

We can now create, read and write our new object. However we are not doing anything with it yet. Like for any other GfsEvent, the action performed is controlled by the event method. If we look at the init_periodic_event function, we see that it returns a gboolean. If TRUE this return value means that the event took place. We then call the event method of the parent class and check its return value. If TRUE we do something and return TRUE. We know that the parent class is GfsInit, an event which takes place once at the start of the simulation. What we want to do is initialise the velocity field using our formula. We can do that like this:

typedef struct {
gdouble * m;
GfsVariable * u;
GfsVariable * v;
} MyData;
static void init_velocity (FttCell * cell, MyData * data)
FttVector pos;
gdouble * m = data->m;
GfsVariable * v = data->v;
GfsVariable * u = data->u;
ftt_cell_pos (cell, &pos);
GFS_VARIABLE (cell, u->i) = - cos (2.*(*m)*M_PI*pos.x)*sin (2.*(*m)*M_PI*pos.y);
GFS_VARIABLE (cell, v->i) = sin (2.*(*m)*M_PI*pos.x)*cos (2.*(*m)*M_PI*pos.y);
static gboolean init_periodic_event (GfsEvent * event, GfsSimulation * sim)
InitPeriodic * inputdata = INIT_PERIODIC(event);
if ((* GFS_EVENT_CLASS (GTS_OBJECT_CLASS (init_periodic_class ())\
->parent_class)->event) (event, sim)) {
/* do object-specific event here */
MyData data;
GfsDomain * domain = GFS_DOMAIN (sim);
data.m = &(inputdata->m);
data.u = gfs_variable_from_name (domain->variables, "U");
data.v = gfs_variable_from_name (domain->variables, "V");
gfs_domain_cell_traverse (GFS_DOMAIN (sim),
(FttCellTraverseFunc) init_velocity,
return TRUE;
return FALSE;

The gfs_domain_cell_traverse function traverses the cell tree and calls the init_velocity function for each leaf cell (FTT_TRAVERSE_LEAFS). A GfsSimulation is an object derived from GfsDomain which justifies the GFS_DOMAIN (sim) type casting. We shall have a look at these structures and functions later in the course.

We also pass an extra parameter to init_velocity: &data. This extra argument cojntains some data obtained from our object through the object InitPeriodic * init .

What happens in init_velocity is straightforward. We first get the position of the center of the cell using ftt_cell_pos and then assign the values of the two components of the velocities using our formula and the m parameter defined in init.

The full source of our new object can be found here

Discussion of the example

Cast to the class

The example shows how one can extract information from an object using a cast. This applies if we want information in a parent. For example, if we want to use the information of the sim variable which came from the parent GfsDomain class, we can have access to it by means of the following variable "domain":

GfsDomain * domain = GFS_DOMAIN (sim);

In this way, we load the information of the "sim" object coming from its parent GfsDomain into this variable and we can either use it or manipulate if we want. For example domain->variables will give us the variables of the simulation sim.


This macro is used everywhere in the code to access the value of a GFS_VARIABLE in a cell. An example of its use is found in the init_periodic object above:

GFS_VARIABLE (cell, v->i) = sin (2.*(*m)*M_PI*pos.x)*cos (2.*(*m)*M_PI*pos.y);

To use this macro, there are two points of view.

  • In the "magical" approach we just now that GFS_VARIABLE(somecell,somevariable->i)allows us to write to the variable somevariable in the cell somecell. For us dummies a question may arise: what value should the index i have ? Is it the cell number ? Well i is not an index: it is a member of each GfsVariable. We do not have to set it: it is already there in the variable.
  • In the rational approach, we would like to read the code and understand the definition of GFS_VARIABLE. However the definition of the GFS_VARIABLE macro is a bit baffling. It can be found in fluid.h.
#define GFS_VARIABLE(cell, index)     ((&GFS_STATE (cell)->place_holder)[index])

Notice that we have a special member of the cell state that tells us where the variable is stored. The macro GFS_STATE is defined by

#define GFS_STATE(cell)               ((GfsStateVector *) (cell)->data)

In some other places of the code a different construction is used

GFS_STATE (cell)->v = sin (2.*init->m*M_PI*pos.x)*cos (2.*init->m*M_PI*pos.y);

In both cases we need to access a member of the cell object which belongs to the FttCell . We thus turn to the description of this class:

The FttCell class

This class is defined as follows:

struct _FttCell {
/*< public >*/
guint flags;
gpointer data;
/*< private >*/
struct _FttOct * parent, * children;

The data member of this class is a gpointer. Remember from last session that it is an untyped pointer, so there is no difficulty to cast to a GfsStateVector by using

struct _GfsStateVector {
/* temporary face variables */
GfsFaceStateVector f[FTT_NEIGHBORS];
/* solid boundaries */
GfsSolidVector * solid;
gdouble place_holder;

Tip: Navigate using the TAGS file

In the above example, we have seen that we can gain some information about the macros, classes, structures and function of the code by navigating in the code from one definition to the next. This navigation is not always immediately enlightening: we must resort to it only when we have enough general understanding of the code. To perform this navigation in a more efficient way than grep, we can use the following trick:

  • Create the TAGS file :
% cd src
% make tags

You can then use the tags in emacs, nedit or vim.

  1. In emacs:
    1. open any *.c or *.h file with emacs.
    2. Position the cursor on a function or variable.
    3. do ESC . to find its definition. (or M-. using the emacs Meta key convention )
    4. do ESC * to return to the previous location (s) .
  2. In vim, nedit ... : vim, nedit users please contribute ...

You can now try to use this trick on your own editor. It is a good idea to try the following: find any place in the code where the GFS_VARIABLE macro is used ---that should be easy--- and then click on it. Type M-. in emacs ---or the equivalent in vim--- and get to the definition, then repeat the operation to display the declarations and definitions above.

Testing our new class

Creating a module

We have now almost all that is needed for our new object. How do we use this new piece of code with gerris? What we want to do is to create a dynamically loadable module. First of all we need to check that we can compile the code. There are a few missing headers we need to add (at the top of the file):

#include <math.h>
#include <stdlib.h>
#include <gfs.h>

The gfs.h header contains all the function declarations for the gerris library. We can now try to compile using for example:

% cc `gfs-config --2D --cflags` -Wall -g -O2 -c init_periodic.c

where gfs-config –2D –cflags (quoted using inverted quotes) defines the compilation flags needed to use the 2D version of the gerris library.

To make a proper module we also need to add the following at the end of the file:

/* Initialize module */
const gchar * g_module_check_init (void);
const gchar * g_module_check_init (void)
init_periodic_class ();
return NULL;

This just tells gerris how to initialise the module after it has been loaded. The g_module_check_init function will be called. In our case, it does only one thing which is to instantiate our new class: this is necessary so that gerris registers how to handle it.

We are now all set to create our new module. Sadly, dynamically loadable module creation is not a standardised process and the command-line arguments vary from compiler to compiler and from system to system. In the following, I will assume that you use gcc on a linux box. If you are using another system supporting dynamically loadable modules, you will need to read your local manual. On a linux box:

% cc `gfs-config --2D --cflags` -Wall -g -O2 -c -fPIC init_periodic.c
% cc -shared init_periodic.o -lc -o

should work. Note that we add the 2D extension to indicate that this module uses the 2D gerris library. We could have built both a 2D and a 3D version of the same module. At runtime the gerris executable uses this extension to check which version to load.

Installing the module

Our module is now ready to use. We just need to install it in a directory where it will be found by gerris. If you installed gerris in /home/joe/local just type:

% cp /home/joe/local/lib/gerris

Using the module

We can now use it directly in a parameter file, for example:

1 2 GfsSimulation GfsBox GfsGEdge {} {
GfsTime { end = 50 }
GfsRefine 6
GModule init_periodic
InitPeriodic {} 2
GfsOutputTime { istep = 10 } stdout
GfsOutputProjectionStats { istep = 10 } stdout
GfsOutputPPM { step = 0.1 } vorticity-%4.2f.ppm { v = Vorticity }
GfsBox {}
1 1 right
1 1 top

The first part of the object definition: InitPeriodic {} is the generic GfsEvent definition, the second part: 2 is the parameter m we added. What if we do not specify the right parameter? Just try to replace InitPeriodic {} 2 with InitPeriodic {} a. You should get a message like

gerris: file 'periodic.gfs' is not a valid simulation file
periodic.gfs:5:18: expecting a number (m)

which tells you where the error occurred (in file periodic.gfs, line 5, character 18) together with the error message which we specified in the read method of our class.

Go to course top level. Go to top of page. Go to next Session.

Personal tools