Plugins

gsl

The interface with gsl spline is available in FreeFEM, the seven kind of spline are

0. gslinterpcspline: default type of spline

1. gslinterpakima

2. gslinterpsteffen

3. gslinterplinear

4. gslinterppolynomial

5. gslinterpcsplineperiodic

6. gslinterpakimaperiodic

A brief wing example given all the syntax:

 1 load "gsl"
 2 
 3 // Parameters
 4 int n = 10;
 5 real[int, int] dspline(2,n+1); //data points to define the spline
 6 for(int i = 0; i <= n; ++i){ //set data points
 7    real xx = square(real(i)/n);
 8    real yy = sin(xx*pi*2);
 9    dspline(0, i) = xx;
10    dspline(1, i) = yy;
11 }
12 
13 // GSL splines
14 gslspline spline1(gslinterpcspline, dspline); //define the spline1
15 gslspline spline11(dspline); //define the spline11
16 gslspline spline2(gslinterpsteffen, dspline); //define the spline2
17 gslspline spline3(gslinterpcspline, dspline(0, :), dspline(1, :));
18 gslspline spline33(dspline(0, :), dspline(1, :)); //define the spline3
19 spline1 = spline2; //copy spline2 in spline1
20 
21 real t = 1.;
22 real s1 = spline1(t); //evaluate the function spline1 at t
23 cout << "spline1(t) = " << s1 << endl;
24 real ds1 = spline1.d(t); //evaluate the derivative of function spline1 at t
25 cout << "spline1.d(t) = " << ds1 << endl;
26 real dds1 = spline1.dd(t); //evaluate the second derivative of function spline1 at t
27 cout << "spline1.dd(t) = " << dds1 << endl;

This can be usefull to build function from data value.

The list of all gsl functions and the FreeFEM equivalent is available in the Language references (same names without _).

ffrandom

Plugin to linux random functions.

The range of the random generator is from \(0\) to \((2^{31})-1\).

 1 load "ffrandom"
 2 
 3 srandomdev(); //set a true random seed
 4 //warning: under window this command
 5 //change the seed by randinit(random())) so all
 6 //FreeFEM random function are changed
 7 
 8 int maxrang = 2^31 - 1;
 9 cout << " max range " << maxrang << endl;
10 
11 cout << random() << endl;
12 cout << random() << endl;
13 cout << random() << endl;
14 
15 srandom(10);
16 cout << random() << endl;
17 cout << random() << endl;
18 cout << random() << endl;

mmap / semaphore

The idea is just try to use Interprocess communication using POSIX Shared Memory in Linux.

We build a small library libff-mmap-semaphore.c and libff-mmap-semaphore.h to easily interface.

• mmap - allocate memory, or map files or devices into memory

• semaphore - allow processes and threads to synchronize their actions

  A semaphore is an integer whose value is never allowed to fall below zero. Two operations can be performed on semaphores: increment the semaphore value by one (sem_post); and decrement the semaphore value by one (sem_wait).

  If the value of a semaphore is currently zero, then a sem_wait operation will block until the value becomes greater than zero.

The functions of library

First the semaphore interface to make synchronization:

• typedef struct FF_P_sem *ff_Psem; the pointer to data structure

• ff_Psem ffsem_malloc(); malloc an empty data structure

• void ffsem_del(ff_Psem sem); clean and free the pointer

• void ffsem_destroy(ff_Psem sem); clean, close the data structure

• void ffsem_init0(ff_Psem sem); make a correct empty of the data structure

• void ffsem_init(ff_Psem sem,const char *nmm, int crea); create or use a new semaphore

• long ffsem_post(ff_Psem sem); nlocked, the value of the semaphore is incremented, and all threads which are waiting on the semaphore are awakened

• long ffsem_wait(ff_Psem sem); the semaphore referenced by sem is locked. When calling sem_wait(), if the semaphore’s value is zero, the calling thread will block until the lock is acquired or until the call is interrupted by a signal.

  Alternatively, the sem_trywait() function will fail if the semaphore is already locked, rather than blocking on the semaphore

• long ffsem_trywait(ff_Psem p);

Secondly, the mmap functions:

• typedef struct FF_P_mmap *ff_Pmmap; the pointer to data structure

• ff_Psem ffmmap_malloc(); malloc an empty data structure

• void ffmmap_del(ff_Pmmap p); clean and free the pointer

• void ffmmap_destroy(ff_Pmmap p); clean, close the data structure

• void ffmmap_init0(ff_Pmmap p); make a correct empty of the data structure

• long ffmmap_msync(ff_Pmmap p, long off, long ln); call writes modified whole pages back to the filesystem and updates the file modification time. Only those pages containing addr and len-1 succeeding locations will be examined.

• void ffmmap_init(ff_Pmmap p, const char *nmm, long len); allocate memory, or map files or devices into memory.

• long ffmmap_read(ff_Pmmap p, void *t, size_t n, size_t off); read n bytes from the mmap at memory off in pointer t.

• long ffmmap_write(ff_Pmmap p, void *t, size_t n, size_t off); write n bytes to the mmap at memory off in pointer t.

The FreeFEM corresponding functions:

• Pmmap sharedata(filename, 1024); new type to store the mmap informations of name store in string filename with 1024 is the size the sharedata zone and file.

• Psemaphore smff("ff-slave", creat); new type to store the semaphore of name ff-slave where creat is a boolean to create or use a existing semaphore.

• Wait(sem) the semaphore referenced by sem is locked. When calling Wait(sem), if the semaphore’s value is zero, the calling thread will block until the lock is acquired or until the call is interrupted by a signal. Alternatively, the trywait(sem) function will fail if the semaphore is already locked, rather than blocking on the semaphore.

• Post(sem) the semaphore referenced by sem is unlocked, the value of the semaphore is incremented, and all threads which are waiting on the semaphore are awakened.

• Read(sharedata ,offset, data); read the variable data from the place offset in sharedata mmap.

• Write(sharedata, offset, data); write the variable data at the place offset in sharedata mmap.

The full example:

The FFMaster.c file:

 1 #include "libff-mmap-semaphore.h"
 2 #include <unistd.h>
 3 #include <stdlib.h>
 4 #include <stdio.h>
 5 ff_Psem sem_ff, sem_c; //the semaphore for mutex
 6 
 7 int main(int argc, const char ** argv)
 8 {
 9    int debug = 0;
10    ff_Pmmap shd;
11    double cff, rff;
12    long status;
13    int i;
14    if (argc > 1) debug = atoi(argv[1]);
15    ff_mmap_sem_verb = debug;
16 
17    sem_ff = ffsem_malloc();
18    sem_c = ffsem_malloc();
19    shd = ffmmap_malloc();
20 
21    ffsem_init(sem_ff, "ff-slave1", 1);
22    ffsem_init(sem_c, "ff-master1", 1);
23    ffmmap_init(shd, "shared-data", 1024);
24 
25    status = 1;
26    ffmmap_write(shd, &status, sizeof(status), 8);
27    ffmmap_msync(shd, 0, 32);
28 
29    char ff[1024];
30    sprintf(ff, "FreeFem++ FFSlave.edp -nw -ns -v %d&", debug);
31    system(ff); //lauch FF++ in batch no graphics
32    if(debug) printf("cc: before wait\n");
33 
34    if(debug) printf("cc: before wait 0 ff\n");
35    ffsem_wait(sem_ff);
36 
37    for (i = 0; i < 10; ++i){
38       printf(" iter : %d \n", i);
39       cff = 10+i;
40       ffmmap_write(shd, &cff, sizeof(cff), 0);
41       ffsem_post(sem_c);
42 
43       if(debug) printf(" cc: before wait 2\n");
44       ffsem_wait(sem_ff);
45       ffmmap_read(shd, &rff, sizeof(rff), 16);
46       printf(" iter = %d rff= %f\n", i, rff);
47    }
48 
49    status = 0; //end
50    ffmmap_write(shd, &status, sizeof(status), 8);
51    ffsem_post(sem_c);
52    printf("End Master \n");
53    ffsem_wait(sem_ff);
54    ffsem_del(sem_ff);
55    ffsem_del(sem_c);
56    ffmmap_del(shd);
57    return 0;
58 }

The FFSlave.edp file:

 1 load "ff-mmap-semaphore"
 2 
 3 Psemaphore smff("ff-slave1", 0);
 4 Psemaphore smc("ff-master1", 0);
 5 Pmmap sharedata("shared-data", 1024);
 6 if (verbosity < 4) verbosity = 0;
 7 
 8 // Mesh
 9 mesh Th = square(10, 10);
10 int[int] Lab = [1, 2, 3, 4];
11 
12 // Fespace
13 fespace Vh(Th, P1);
14 Vh u, v;
15 
16 // Macro
17 macro grad(u) [dx(u), dy(u)] //
18 
19 int status = 1;
20 cout << " FF status = " << status << endl;
21 real cff, rff;
22 
23 // Problem
24 problem Pb (u, v)
25    = int2d(Th)(
26         grad(u)'*grad(v)
27    )
28    - int2d(Th)(
29         cff*v
30    )
31    + on(Lab, u=0)
32    ;
33 
34 if (verbosity > 9) cout << " FF: before FF post\n";
35 Post(smff); //unlock master end init
36 
37 while (1){
38    if (verbosity > 9) cout << " FF: before FF wait \n";
39    Wait(smc); //wait from cint write ok
40    Read(sharedata, 0, cff);
41    Read(sharedata, 8, status);
42 
43    cout << " After wait .. FF " << cff << " " << status << endl;
44    if(status <= 0) break;
45 
46    // Solve
47    Pb;
48    rff = int2d(Th)(u*u);
49    cout << " ** FF " << cff << " " << rff << endl;
50 
51    // Write
52    Write(sharedata, 16, rff);
53    Post(smff); //unlock cc
54 }
55 
56 Post(smff); //wait from cint
57 cout << " End FreeFEM " << endl;

To test this example of coupling C program and FreeFEM script:

1 cc -c libff-mmap-semaphore.c
2 cc FFMaster.c -o FFMaster libff-mmap-semaphore.o -g -pthread
3 ff-c++ -auto ff-mmap-semaphore.cpp
4 ./FFMaster

The output:

 1 len 1024 size 0
 2 len 1024 size 1024
 3 FF status = 1
 4 iter : 0
 5 After wait .. FF 10 1
 6 ** FF 10 0.161797
 7 iter = 0 rff= 0.161797
 8 iter : 1
 9 After wait .. FF 11 1
10 ** FF 11 0.195774
11 iter = 1 rff= 0.195774
12 iter : 2
13 After wait .. FF 12 1
14 ** FF 12 0.232987
15 iter = 2 rff= 0.232987
16 iter : 3
17 After wait .. FF 13 1
18 ** FF 13 0.273436
19 iter = 3 rff= 0.273436
20 iter : 4
21 After wait .. FF 14 1
22 ** FF 14 0.317121
23 iter = 4 rff= 0.317121
24 iter : 5
25 After wait .. FF 15 1
26 ** FF 15 0.364042
27 iter = 5 rff= 0.364042
28 iter : 6
29 After wait .. FF 16 1
30 ** FF 16 0.414199
31 iter = 6 rff= 0.414199
32 iter : 7
33 After wait .. FF 17 1
34 ** FF 17 0.467592
35 iter = 7 rff= 0.467592
36 iter : 8
37 After wait .. FF 18 1
38 ** FF 18 0.524221
39 iter = 8 rff= 0.524221
40 iter : 9
41 After wait .. FF 19 1
42 ** FF 19 0.584086
43 iter = 9 rff= 0.584086
44 End Master
45 After wait .. FF 19 0