Hi everyone,
today I was trying one of the examples that were worked out during the basic tutorial of FreeFem Days 2019.
I was trying to visualise the mesh in Paraview.
The code is correct and medit gives the right output. however, when I try to load the vtk file onto paraview, the console display a message saying that the file cannot be opened.
I am using ArchLinux.
Could you please help me understand if I am doing something wrong? If I plot in medit the mesh is displayed correctly.
Best regards,
Paolo
PS. Since I do not have the permissions to upload the file, I type the code here below:
load "msh3"
load "medit"
load "iovtk"
real R =1, theta= 20*pi/180, hcap =0.05;
real hsize=hcap/3;
real y2=R*sin(theta);
real x2=R*cos(theta);
real theta1 = asin(y2/(R-hcap));
real y1=(R-hcap)*sin(theta1);
real x1 =(R-hcap)*cos(theta1);
assert( abs(y1-y2)< 1e-10);
border cext(t=theta,pi/2){ x = R*cos(t); y = R*sin(t); label =1;}
border cint(t=theta1,pi/2){ x = (R-hcap)*cos(t); y = (R-hcap)*sin(t); label =2;}
border vseg (t= R-hcap, R) { x=0; y=t; label = 3;}
border vhor (t= x1,x2) { x=t; y=y1; label = 4;}
int nn = (pi/2-theta)*R/hsize, mm=hcap/hsize;
plot( cext(nn)+cint(-nn)+vseg(-mm)+vhor(mm));
mesh Th2=buildmesh(cext(nn)+cint(-nn)+vseg(-mm)+vhor(mm));
Th2=movemesh(Th2,[y,x]);
plot(Th2,wait=1);
int MaxLayersT=(int(2*pi*R/hsize)/4)*4;
func zminT = 0;
func zmaxT = 2*pi;
func fx= y*cos(z);// / max( abs(cos(z) ), abs(sin(z)));
func fy= y*sin(z);// / max( abs(cos(z) ), abs(sin(z)));
func fz= x;
int[int] r1T=[0,0], r2T=[0,0,2,2];
int[int] r4T=[0,2];
mesh3 Th=buildlayers(Th2,coef= max(.01,y/max(x,0.4) ), MaxLayersT,zbound=[zminT,zmaxT],transfo=[fx,fy,fz],facemerge=1, region=r1T, labelmid=r2T);
medit("cap",Th,wait=1);
// FFCS: testing 3d plots
plot(Th,cmm="cap");
savevtk("my_shell.vtk",Th);