FreeFem Commands (taken from freefem++doc)

What does FreeFem do?

FreeFem is a "high level integrated development environment for numerically solving partial differential equations." It has an automatic mesh generator, is capable of a posteriori mesh adaptation, and has a general purpose elliptic solver interfaced with fast algorithms.

Writing and running FreeFem programs

FreeFem programs should be saved in files with extensions .edp

To run the program poisson.edp, type after the prompt:   FreeFem++ poisson.edp

The symbol // in a FreeFem program is used to write comment statements.

Defining the boundary of the domain

border C(t=0,2*pi) { x=cos(t); y=sin(t); label=1; } // (label optional) defines the boundary to be the unit circle

border Gamma1(t=0,1) {x=t; y=0; label=1;}
border Gamma2(t=0,1) {x=1; y=t; label=1;}
border Gamma3(t=0,1) {x=1-t; y=1; label=1;}
border Gamma4(t=0,1) {x=0; y=1-t; label=1;}
// defines the boundary to be the unit square, traversed counterclockwise
// defining all labels to be the same allows us to specify a single boundary condition on all four sides. Assigning different labels to different sides allows us to specify different boundary conditions on different sides.

Mesh commands

mesh Th=square(20,20);
// generates a 20 x 20 grid in the unit square [0,1]^2.
// Note that with this command, the labels are defined automatically to be 1,2,3,4 on the four sides of the square starting at the bottom and moving counterclockwise.

real x0=1.2,x1=1.8;
real y0=0,y1=1;
int n=5,m=20;
mesh Th=square(n,m,[x0+(x1-x0)*x,y0+(y1-y0)*y]);
// constructs a n by m grid in the rectangle [x0,x1] x [y0,y1]

border C(t=0,2*pi) { x=cos(t); y=sin(t); label=1; } // label optional
mesh Th = buildmesh(C(10));
// constructs a mesh on the unit circle with 10 boundary lines

mesh Th1 = buildmesh (Gamma1(10)+ Gamma2(10)+Gamma3(10)+ Gamma4(10)); // constructs a mesh on the unit square with 11 boundary vertices on each edge.

plot(Th, wait=1);
// plots the mesh and waits for a mouse or keyboard command;
// Clicking the left mouse button moves on.

plot(Th, wait=1, ps = "");
// displays mesh and stores plot in postscript file

savemesh("mesh_sample.msh"); // saves the mesh to a file

mesh th2 = readmesh("mesh_sample.msh"); // reads the mesh into another program

Th = splitmesh(Th,1, split=4);
// gives a uniform mesh refinement

Th = adaptmesh(Th,u);
// gives an adaptive mesh refinement based on the approximate solution

Defining finite element spaces

fespace Vh1(Th,P1); // defines Vh1 as the space of Lagrange P1 elements on the mesh Th
fespace Vh2(Th,P2); // defines Vh2 as the space of Lagrange P2 elements on the mesh Th

Data types

int i,n=20;
real L2error, graderror;
real[int] xx(n), yy(n);
The first statement defines i and n to be integers and sets n=20.
The second statement defines L2error and graderror to be real numbers.
The third statement defines xx and yy to be vectors of length n=20.

Defining functions

func f= x*(1-x)*y*(1-y);
func dxf = (1-2*x)*y*(1-y);
// Function statements defining the functions f and dxf.


The for and while loops are implemented with break and continue keywords. In the for-loop, there are three parameters; the INITIALIZATION of a control variable, the CONDITION to continue, and the CHANGE of the control variable. While CONDITION is true, the for-loop continues.


Below, the sum from 1 to 10 is calculated by a for-loop (the result is in sum).

int sum=0;
for (int i=1; i<=10; i++)
sum += i;

The while-loop

while (CONDITION) { BLOCK of calculations or change of control variables }

is executed repeatedly until CONDITION become false. The sum from 1 to 10 is also calculated by while as follows:

int i=1, sum=0;
while (i<=10) {
sum += i; i++;


To write to the screen, we use statements such as   cout << u << endl;   to print the contents of the variable u, or   cout << "Solution =" << u << endl;   to print the text   Solution =   followed by the contents of the variable u.

To input from the keyboard, we use statements such as

int n; // defines an integer variable n
cin >> n; // the program then waits until a value of n has been input from the keyboard (followed by a press of the RETURN key).

FreeFem Documentation

For more complete documentation on FreeFem, consult freefem++doc.pdf

Sample FreeFem program

Solving - Lap u = f on the unit circle with zero boundary conditions, computing the L2 and H1 errors, and obtaining a plot of the solution

// file: laplace.edp
border C(t=0,2*pi){x=cos(t); y=sin(t); label=1;}; // defines the boundary
mesh Th = buildmesh (C(50)); // the triangulated domain Th is on the left side of its boundary
plot(Th, wait=1, ps = "");
// displays mesh and stores plot in postscript file
fespace Vh(Th,P1); // defines FE space as C^0 piecewise P1 functions
Vh u,v; // the finite element space defined over Th is called Vh
real L2error, graderror;
func f= 8*(1-2*x^2-2*y^2); // definition of a called f function
func uexact = (1- x^2-y^2)^2;
func dxuexact= -4*x*(1- x^2-y^2);
func dyuexact= -4*y*(1- x^2-y^2);
solve Poisson(u,v,solver=LU) = // defines the PDE
int2d(Th)(dx(u)*dx(v) + dy(u)*dy(v)) // bilinear form
- int2d(Th)( f*v) // right hand side
+ on(1,u=0) ; // Dirichlet boundary condition
L2error= sqrt(int2d(Th)((u-uexact)^2));
cout << " L2error = " << L2error << endl;
graderror = sqrt(int2d(Th)((dx(u)-dxuexact)^2 + (dy(u)-dyuexact)^2));
cout << " graderror = "<< graderror << endl;