I am new to Freefem and have trouble cycling movemesh. My script is equivalent to (structure is the same, formulas are wrong):
/* iterated rho is the driving parameter */
mesh Th = buildmesh( left(60) + bottom(60) + right(60) + top(60)+ hole(-60));
int nstep=3;
int niter=3;
int iter,istep;
real kappa0 = 83.0;
real kmu=1e4;
fespace Vh(Th,P2);
fespace Vh0(Th,P1);
fespace Vh00(Th,P0);
fespace Eh(Th,HCT);
Vh ux,uy,vx,vy;
Eh [Eixx,Eixxx,Eixxy],[Eiyy,Eiyyx,Eiyyy],[Eixy,Eixyx,Eixyy];
Eh [Fixx,Fixxx,Fixxy],[Fiyy,Fiyyx,Fiyyy],[Fixy,Fixyx,Fixyy];
Vh00 rho, Mises;
rhoini=1e-3;
rho=rhoini;
macro e(ux,uy) [dx(ux),dy(uy),(dx(uy)+dy(ux))/sqrt2]//
macro EE(Eixx,Eiyy,Eixy) [Eixx,Eiyy,sqrt2*Eixy]//
macro tre(ux,uy) (dx(ux)+dy(uy))//
macro trEE(Eixx,Eiyy,Eixy) (Eixx+Eiyy)//
macro Wleft(vx,vy) (-vx)//
macro mut(rho) (kmu/rho)//
/* LOAD STEP LOOP /
for (istep=1;istep<nstep;istep=istep+1)
{
problem incomp([ux,uy,Eixx,Eixxx,Eixxy,Eiyy,Eiyyx,Eiyyy,Eixy,Eixyx,Eixyy],[vx,vy,Fixx,Fixxx,Fixxy,Fiyy,Fiyyx,Fiyyy,Fixy,Fixyx,Fixyy],solver=Crout) =
int2d(Th)((kappa0)(tre(ux,uy)+trEE(Eixx,Eiyy,Eixy))(tre(vx,vy)+trEE(Fixx,Fiyy,Fixy))
+2mut(rho)(e(ux,uy)+EE(Eixx,Eiyy,Eixy))'(e(vx,vy)+EE(Fixx,Fiyy,Fixy)))
-int1d(Th,GammaL)(2*Wleft(vx,vy));
for (iter=1;iter< niter;iter=iter+1)
{
LOOP(istep, iter) for rho readjustment;
Mises=4mut(rho)^2((-1/3)(tre(ux,uy)+trEE(Eixx,Eiyy,Eixy))(tre(ux,uy)+trEE(Eixx,Eiyy,Eixy))+(e(ux,uy)+EE(Eixx,Eiyy,Eixy))'*(e(ux,uy)+EE(Eixx,Eiyy,Eixy)));
}
/* REINTERPOLATE FEs on NEW MESH Th2 > Th */
mesh Th2= movemesh(Th,[x+ux,y+uy]);
plot(Th2,wait=nwait,ps=“displacement-Th2_”+ istep +“.eps”);
/* Finite element spaces */
fespace Vh2(Th2,P2);
fespace Vh002(Th2,P0);
fespace Eh2(Th2,HCT);
Vh2 ux2,uy2,vx2,vy2;
Eh2 [Eixx2,Eixxx2,Eixxy2],[Eiyy2,Eiyyx2,Eiyyy2],[Eixy2,Eixyx2,Eixyy2];
Eh2 [Fixx2,Fixxx2,Fixxy2],[Fiyy2,Fiyyx2,Fiyyy2],[Fixy2,Fixyx2,Fixyy2];
Vh002 Mises2;
// shift to new mesh
ux2=ux(x-ux,y-uy);uy2=uy(x-ux,y-uy);vx2=vx(x-ux,y-uy);vy2=vy(x-ux,y-uy);
[Eixx2,Eixxx2,Eixxy2]=[Eixx(x-ux,y-uy),Eixxx(x-ux,y-uy),Eixxy(x-ux,y-uy)];
[Eiyy2,Eiyyx2,Eiyyy2]=[Eiyy(x-ux,y-uy),Eiyyx(x-ux,y-uy),Eiyyy(x-ux,y-uy)];
[Eixy2,Eixyx2,Eixyy2]=[Eixy(x-ux,y-uy),Eixyx(x-ux,y-uy),Eixyy(x-ux,y-uy)];
[Fixx2,Fixxx2,Fixxy2]=[Fixx(x-ux,y-uy),Fixxx(x-ux,y-uy),Fixxy(x-ux,y-uy)];
[Fiyy2,Fiyyx2,Fiyyy2]=[Fiyy(x-ux,y-uy),Fiyyx(x-ux,y-uy),Fiyyy(x-ux,y-uy)];
[Fixy2,Fixyx2,Fixyy2]=[Fixy(x-ux,y-uy),Fixyx(x-ux,y-uy),Fixyy(x-ux,y-uy)];
Mises2=Mises(x-ux,y-uy);
// rename mesh
Th=Th2;
// rename FE functions to initial names
ux=ux2;uy=uy2;vx=vx2;vy=vy2;
[Eixx,Eixxx,Eixxy]=[Eixx2,Eixxx2,Eixxy2];
[Eiyy,Eiyyx,Eiyyy]=[Eiyy2,Eiyyx2,Eiyyy2];
[Eixy,Eixyx,Eixyy]=[Eixy2,Eixyx2,Eixyy2];
[Fixx,Fixxx,Fixxy]=[Fixx2,Fixxx2,Fixxy2];
[Fiyy,Fiyyx,Fiyyy]=[Fiyy2,Fiyyx2,Fiyyy2];
[Fixy,Fixyx,Fixyy]=[Fixy2,Fixyx2,Fixyy2];
Mises=Mises2;
// repeat macros
macro Wleft(vx,vy) (-vx)//
macro e(ux,uy) [dx(ux),dy(uy),(dx(uy)+dy(ux))/sqrt2]//
macro EE(Eixx,Eiyy,Eixy) [Eixx,Eiyy,sqrt2*Eixy]//
macro tre(ux,uy) (dx(ux)+dy(uy))//
macro trEE(Eixx,Eiyy,Eixy) (Eixx+Eiyy)//
}
It does not work but I cannot see why. Would you please help?