Plot in 3d with .vtu and paraview


i have to plot the Following function:

q(x/eps,y/eps,z/eps)=1 if (x/eps,y/eps,z/eps) \in B_{\delta}+ Z^3, else 0.

where B_{delta} is the ball with radius

\delta= ((3/4)(1/(wpi)))^{1/3}

we choise eps=0.1 and w=2000.

The problem is that in .vtu we Don’t ditingue the ball when q=1 to the rest when q=0.

I search in the doc but i Don’t find how we save .vtu in 3d (i choice 3d in paraview). So there is an solution?

And the image after copilation Don’t show the ball.

Thank you in advance to the help. There is the code:

’ ’ ’ freefem
load “UMFPACK64”

load “iovtk”

load “msh3”


func real q (real t, real s, real c,real w) {
real t1 = t - rint(t);
real s1 = s - rint(s);
real c1 = c -rint©;
real delta= ((3./4.)(1./(piw)))^(1/3);
if ( square(s1)+square(t1)+square(c1) <= square(delta))
return 1;
return 0.0;
real eps=0.2, w=2000.;

mesh3 Th=cube(30,30,30);
fespace Wh(Th,P2,periodic=[[1,x,z],[3,x,z],[2,y,z],[4,y,z],[5,x,y],[6,x,y]]);
Wh uplot;
int[int] Order= [1];

savevtk(“q.vtu”,Th,uplot,dataname=“qdelta”, order=Order);

Kind regards