/usr/share/doc/freefem++/examples/examples++-3d/Laplace-Adapt-aniso-3d.edp is in freefem++-doc 3.19.1-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 | // in test wrong skip
if(0)
{
load "msh3"
load "tetgen"
load "mshmet"
load "mmg3d-v4.0"
load "freeyams"
load "medit"
searchMethod=1; // more safe seach algo ..
macro det2(a,b,c,d) (a*d-c*b) //
macro det3s(m11,m21,m31,m22,m32,m33)
( m11*det2(m22,m32,m32,m33)
-m21*det2(m21,m31,m31,m33)
+m31*det2(m11,m21,m21,m22) ) //
int nn = 6;
int[int] l1111=[1,1,1,1],l01=[0,1],l11=[1,1];// label numbering
mesh3 Th3=buildlayers(square(nn,nn,region=0,label=l1111),
nn, zbound=[0,1], labelmid=l11, labelup = l01, labeldown = l01);
// remove the $]0.5,1[^3 cube$
Th3 = trunc(Th3,(x<0.5) | (y < 0.5) | (z < 0.5) ,label=1);
fespace Vh(Th3,P1);
fespace Mh(Th3,[P1,P1,P1,P1,P1,P1]);
Vh u,v,usol,h3;
Mh [m11,m21,m22,m31,m32,m33];
macro Grad(u) [dx(u),dy(u),dz(u)] // EOM
problem Poisson(u,v,solver=CG) = int3d(Th3)( Grad(u)'*Grad(v) ) // ') for emacs
-int3d(Th3)( 1*v ) + on(1,u=0);
int nbregu=1;
int[int] MSHloptions=[1,0,0,0,verbosity,nbregu,0];// mesh aniso ..
real[int] MSHdoptions=[1e-3,0.2,0.01,0];
real lerr=0.01;
verbosity=4;
for(int ii=0; ii<2; ii++) // BUG trap in interation 3
{
Poisson;
plot(u,wait=1);
h3=0;
[m11,m21,m22,m31,m32,m33]=[0,0,0,0,0,0];
cout <<" u min, max = " << u[].min << " "<< u[].max << endl;
real cc=(u[].max-u[].min);// rescale coefficiant
// real[in t] met=mshmet(Th3,u,loptions=MSHloptions,doptions=MSHdoptions);
real[int] met=mshmet(Th3,u,hmin=1e-3,hmax=0.2,err=lerr,aniso=1);
// cot << met.n << " " << met2.n << endl;
m11[]=met;
// savemesh(Th3,"Th3.mesh");
// savesol("Th3.sol",Th3, [m11,m21,m22,m31,m32,m33]);
medit("met",Th3,[m11,m21,m22,m31,m32,m33]);
cout << " met size "<< met.n << " " << m11[].n << endl;
verbosity=1;
mesh3 Thb=freeyams(Th3,metric=m11[],hmin=0.0001,hmax=0.3,gradation=2.,verbosity=-10);
h3=det3s(m11,m21,m31,m22,m32,m33);
real[int] rltetg = [0.1,0.1,.1,1,lerr*5];
Th3=tetg(Thb,switch="aAQpYY",regionlist=rltetg);
[m11,m21,m22,m31,m32,m33]=[m11,m21,m22,m31,m32,m33];
h3 = 1/sqrt(det3s(m11,m21,m31,m22,m32,m33));
cout <<" h3 min, max = " << h3[].min << " "<< h3[].max << endl;
Th3=tetgreconstruction(Th3,switch="raAQYY",sizeofvolume=h3/6.);
[m11,m21,m22,m31,m32,m33]=[m11,m21,m22,m31,m32,m33];
medit("met3",Th3,[m11,m21,m22,m31,m32,m33]);
int[int] mopt=[1,0,64,0,0,-10];
fespace Ph(Th3,P0);
Ph vol= volume;
medit("Vol",Th3,vol,order=0);
Th3=mmg3d(Th3,metric=m11[]);
[m11,m21,m22,m31,m32,m33]=[m11,m21,m22,m31,m32,m33];
Th3=mmg3d(Th3,metric=m11[]);
medit("Th3-"+ii,Th3,wait=1);
lerr *= 0.6;// change the level of error
cout << " Th3" << Th3.nv < " " << Th3.nt << endl;
medit("U-adap-iso-"+ii,Th3,u,wait=1);
}
cout <<"end Laplace Adapt aniso 3d. edp " <<endl;
}
|