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.
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.
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 = "laplacemeshplot.ps");
// displays mesh and stores plot in postscript file laplacemeshplot.ps
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
for (INITIALIZATION; CONDITION; CHANGE) { BLOCK of calculations }
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 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).
// 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 = "laplacemeshplot.ps");
// displays mesh and stores plot in postscript file laplacemeshplot.ps
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;