// Parameters int nbsolids=2; // Number of solids real ratio=0.5; // Minimal distances (relative to the edges) - as to be lower that 1. real alpha=0.05; // Uzawa step real deltamin=1.e-8; // Precision // to plot the current deformation mesh[int] Shdef(nbsolids); macro plotdef(waitornot,comment) for (int i=0;i0); n1(i1*nbvertices+j1,i2*nbvertices+j2)= alpha*ty ; n2(i1*nbvertices+j1,i2*nbvertices+j2)=-alpha*tx; } else { if ( sg >= 0) { t1=(phi1[i2][][a0(j2)]-phi1[i1][][a0(j1)])^2+(phi2[i2][][a0(j2)]-phi2[i1][][a0(j1)])^2; t1=sqrt(t1); n1(i1*nbvertices+j1,i2*nbvertices+j2)=-(phi1[i2][][a0(j2)]-phi1[i1][][a0(j1)])/t1; n2(i1*nbvertices+j1,i2*nbvertices+j2)=-(phi2[i2][][a0(j2)]-phi2[i1][][a0(j1)])/t1;} else { t1=(phi1[i2][][a0(j2)]-phi1[i1][][a1(j1)])^2+(phi2[i2][][a0(j2)]-phi2[i1][][a1(j1)])^2; t1=sqrt(t1); n1(i1*nbvertices+j1,i2*nbvertices+j2)=-(phi1[i2][][a0(j2)]-phi1[i1][][a1(j1)])/t1; n2(i1*nbvertices+j1,i2*nbvertices+j2)=-(phi2[i2][][a0(j2)]-phi2[i1][][a1(j1)])/t1; }}} else { n1(i1*nbvertices+j1,i2*nbvertices+j2)=0; n2(i1*nbvertices+j1,i2*nbvertices+j2)=0;}} // Intialization of the Lagrange Multipliers lambda0=0.;lambda1=0.; //for (int iterconvex=0;iterconvex<10;iterconvex++){ int iterconvex=0; int iterUzawa=0; Vh1 phiprec1,phiprec2; // The deformations while (iterUzawa!=1) { iterUzawa=0; // Computation of the initial normals computen(); { fespace Wh1(Shdef[0],P1); Wh1 Nt1,Nt2; Nt1=0;Nt2=0; int p=3; for (int i=0;i"<"<0)+(C1(i1*nbvertices+j1,i2*nbvertices+j2)>0); if (!((C0(i1*nbvertices+j1,i2*nbvertices+j2)"<"<