/usr/share/doc/freefem++/examples/examples++-load/testFE-PkEdge.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 | load "Element_PkEdge"
// a macro the compute numerical derivative
macro DD(f,hx,hy) ( (f(x1+hx,y1+hy)-f(x1-hx,y1-hy))/(2*(hx+hy))) //
mesh Th=square(1,1,[10*(x+y/3),10*(y-x/3)]);
// Th=square(1,1);
//savemesh(Th,"Th.msh");
real x1=0.7,y1=0.9, h=1e-7;
int it1=Th(x1,y1).nuTriangle;
macro Check(PkEdge)
{
cout << " Test Finite Element " <<endl << endl;
fespace Vh(Th,PkEdge);
Vh a1,b1,c1;
c1=x+y;
real[int] viso(100);
real v0= -100;
real dv = 200./(viso.n+1);
for(int i=0;i<viso.n;++i)
viso[i]=v0+i*dv;
plot(c1,fill=1,wait=1,viso=viso);
for (int i=0;i<Vh.ndofK;++i)
cout << i << " " << Vh(0,i) << endl;
for (int i=0;i<Vh.ndofK;++i)
{
a1[]=0;
int j=Vh(it1,i);
a1[][j]=1;
plot(a1, wait=1);
b1=a1;/* do the interpolation */
c1[] = a1[] - b1[];
cout << " ---------" << i << " " << c1[].max << " " << c1[].min << endl;
cout << " a = " << a1[] <<endl;
cout << " b = " << b1[] <<endl;
assert(c1[].max < 1e-9 && c1[].min > -1e-9); /* check if the interpolation is correct */
/* check the derivative and numerical derivative */
/* no derivative in this element. */
/*
cout << " dx(a1)(x1,y1) = " << dx(a1)(x1,y1) << " == " << DD(a1,h,0) << endl;
assert( abs(dx(a1)(x1,y1)-DD(a1,h,0) ) < 1e-5);
assert( abs(dx(a2)(x1,y1)-DD(a2,h,0) ) < 1e-5);
assert( abs(dy(a1)(x1,y1)-DD(a1,0,h) ) < 1e-5);
assert( abs(dy(a2)(x1,y1)-DD(a2,0,h) ) < 1e-5);
*/
}
} //EOM
Check(P1edge)
Check(P2edge)
Check(P3edge)
Check(P4edge)
Check(P5edge)
|