/usr/share/doc/freefem++/examples/examples++-mpi/Laplace3d-hypre.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 90 91 92 93 94 95 96 97 98 99 | /*
Warning in before version 3.2-1 the nomarl are interal normal
after the signe is correct and now the noral was exterior normal.
*/
verbosity=2;
load "msh3"
load "hypre_FreeFem"
// for hypre
int[int] iparm(20);
real[int] dparm(6);
for(int iii=0;iii<20;iii++)iparm[iii]=-1;
for(int iii=0;iii<6;iii++) dparm[iii]=-1;
iparm[0]=2; // PCG as krylov method
iparm[1]=0; // AMG as preconditionner 2: if ParaSails
iparm[7]=7; // Interpolation
iparm[9]=6; // AMG Coarsen type
iparm[10]=1; // Measure type
iparm[16]=2; // Additive schwarz as smoother
dparm[0]=1e-13; // Tolerance to convergence
dparm[1]=5e-4; // Threshold
dparm[2]=5e-4; // truncation factor
int nn=10;
mesh Th2=square(nn,nn,region=0);
fespace Vh2(Th2,P1);
Vh2 ux,uz,p2;
int[int] rup=[0,2], rdown=[0,1], rmid=[1,1,2,1,3,1,4,1];
real zmin=0,zmax=1;
mesh3 Th=buildlayers(Th2,nn,
zbound=[zmin,zmax],
labelmid=rmid,
reffaceup = rup,
reffacelow = rdown);
fespace Vh(Th,P1);
func ue = 2*x*x + 3*y*y + 4*z*z + 5*x*y+6*x*z+1;
func uex= 4*x+ 5*y+6*z;
func uey= 6*y + 5*x;
func uez= 8*z +6*x;
func f= -18. ;
Vh uhe = ue; //
cout << " uhe min: " << uhe[].min << " max:" << uhe[].max << endl;
Vh u,v;
macro Grad3(u) [dx(u),dy(u),dz(u)] // EOM
/*
problem Lapl3d(u,v) =
int3d(Th)(Grad3(v)' *Grad3(u)) //') for emacs
+ int2d(Th,2)(u*v)
- int3d(Th)(f*v)
- int2d(Th,2) ( ue*v + (uex*N.x +uey*N.y +uez*N.z)*v )
+ on(1,u=ue);
Lapl3d;
*/
varf Lapl3d(u,v) =
int3d(Th)(Grad3(v)' *Grad3(u)) //') for emacs
+ int2d(Th,2)(u*v)
+ int3d(Th)(f*v)
+ int2d(Th,2) ( ue*v + (uex*N.x +uey*N.y +uez*N.z)*v )
+ on(1,u=ue);
real tgv = -1;
matrix A = Lapl3d(Vh,Vh,tgv=tgv);
cout << A.n << " " << A.m << " " << A.nbcoef << endl;
set(A,solver=sparsesolver,dparams=dparm, lparams=iparm);
cout <<A.n <<" " << A.m << " " << A.nbcoef << endl;
real[int] b = Lapl3d(0,Vh,tgv=tgv);
u[] = A^-1*b;
cout << " b min:: " << b.min << " max: " << b.max << endl;
cout << " u min:: " << u[]. min << " max: " << u[].max << endl;
real err= int3d(Th)( square(u-ue) );
real aa1= int3d(Th,qfV=qfV1)(u) ;
real aa2= int3d(Th,qfV=qfV1lump)(u) ;
cout << " aa1 = " << aa1 << endl;
cout << " aa2 = " << aa2 << endl;
cout << int3d(Th)(1.) << " = " << Th.mesure << endl;
cout << int3d(Th,qfV=qfV1)(1.) << " = " << Th.mesure << endl;
cout << int3d(Th,qfV=qfV1lump)(1.) << " = " << Th.mesure << endl;
Vh d= ue-u;
cout << " err = " << err << " diff l^\intfy = " << d[].linfty << endl;
real aire2=int2d(Th,2)(1.); // bug correct in version 3.0-4
cout << " aire2 = " << aire2 << endl;
func uuu= 2.+x;
cout << uuu(1,1,1) << endl;
assert( abs(aire2-1.) < 1e-6);
plot(u,wait=1);
assert( err < 1e-3);
|