Au début on fera très simple, par exemple "solve(u)". Il faudra donc rajouter l'instruction solve, mais l'argument, ici "u", est une variable et non une expression. Il faudra vérifier que s'il elle n'existe pas elle soit créée. Mais dans un premier temps on peut communiquer les données au solveur entièrement par fichier. Ainsi le solveur, ex dans le cas d'un Laplacien, sera une fonction sans paramètre (ex: void solvepde()) qui lira toutes ses données, la triangulation dans un fichier "test.msh", le second membre dans "f.dta", les conditions de Dirichlet dans "g.dta" et de Neumann dans "h.dta". Cela donne le programme Gfem suivant:
{ pi = 4*atan(1); buildmesh(x,y,ref,1000) { border(t=0,pi*2,50) { ref = 1; x = cos(t); y = sin(t); }; border(t=0,2*pi,20) { ref = 2; x = 0.3+0.3*cos(-t); y = 0.3*sin(-t); }; }; savemesh("test.msh"); fonction2 f(u,v) f = -4; fonction g(w) g = (x*x+y*y)*(ref>=2); fonction h(w) h = 2; plot(g(1)); save("f.dta",f(u,v)); save("g.dta",g(w)); save("h.dta",h(w)); solve(p); }On affichera les résultats à l'intérieur de la fonction Isolve::execute() de sorte que le nom "p",ici, n'intervient pas.
Dans un deuxième temps, et là il y a une difficulté réelle, on ajoutera une structure pour pouvoir utiliser la fonction p dans la suite du programme, par exemple si en fin de programme Gfem on avait plot(2*p-1). Le problème réside dans le fait que la variable p est en fait une fonction mais il lui manque la fonction "eval()".