# How to get values of a pde at every position in the space?

Hi,

I may be asking a stupid question, but how can I get the values of the fespace after the problem was solved. Let’s say, I implemented a diffusion pde. I can plot it with the plot function of freefem++. But how can I get the value of the function for each position in the space ? I would say I need to access every finite elements, its geometry and value. Or maybe there is a function to do that.

Mordokkai

HI,

I think, you question is a little be strange because, to get the valeur a point (x0,y0) in 2d of finite element function u
juste du u(x0,y0) (same in 3d with u(x0,y0,z0) )

Event for basis function like (sin(x)) you do not know the value at each point

this peace of code given a idea:
// value a each node of any Lagrange Finite element. here P2

mesh Th=square(10,10);
fespace Vh(Th,P2);

Vh xx=x, yy=x;
Vh u =sin(x)*cos(y);

for(int i=0; i<Vh.ndof; ++i)
cout << xx[][i] << " " << yy[][i] << " u = "<< u[][i] << endl;

I think, you question is a little be strange because, to get the valeur a point (x0,y0) in 2d of finite element function u
juste du u(x0,y0) (same in 3d with u(x0,y0,z0) )

That’s the first thing I tried but I got strange results.
Below you can find the result plot by FreeFem++: plot(u, wait=true, fill=true, dim=2);

Then I try to write everything in a text file:
ofstream out (“out.txt” );
out.fixed;
for(int i=0; i<dims[0]; i++){
for(int j=0; j<dims[1]; j++){
out << u(i,j) << " ";
}
out << endl;
}

And then I try to reopen it with python matplotlib:
import matplotlib.pyplot as plt
import numpy as np
plt.matshow(mat)
plt.show()

And I get the following result:

If I take a basic example:

real[int] dims = [50,50];
border a(t=0, dims[0]){x=t; y=0; label=1;};
border b(t=0, dims[1]){x=dims[0]; y=t; label=2;};
border f(t=0, dims[0]){x=dims[0]-t; y=dims[1]; label=3;};
border d(t=0, dims[1]){x=0; y=dims[1]-t; label=4;};
border g(t=10, dims[0]-10){x=t; y=10; label=5;};
border h(t=10, dims[1]-10){x=dims[0]-10; y=t; label=6;};
border i(t=10, dims[0]-10){x=dims[0]-t; y=dims[1]-10; label=7;};
border j(t=10, dims[1]-10){x=10; y=dims[1]-t; label=8;};

int n=3;
mesh fm = buildmesh(a(n) + b(n) + f(n) + d(n) + g(-n) + h(-n) + i(-n) + j(-n));
plot(fm);

fespace Vh(fm, P1);

func real heatsource(){
return 0;
}

Vh u, v;
Vh myHeat = heatsource();

problem heat(u,v) = int2d(fm)(dx(u)*dx(v) + dy(u)dy(v)) - int2d(fm)(myHeatv)
+ on(1, u=0) + on(2, u=0) + on(3, u=0) + on(4, u=0) + on(5, u=0) + on(6, u=0)
+ on(7, u=10) + on(8, u=0);
heat;
plot(u, wait=true, fill=true, dim=2);

ofstream out (“map.txt” );
out.fixed;
for(int i=0; i<dims[0]; i++){
for(int j=0; j<dims[1]; j++){
out << u(i,j) << " ";
}
out << endl;
}

Then freefem display:

While matplotlib display:

It seems like some values are assigned inside the interior square, while not displayed by freefem

Hello
I often use something like (in the case of magnetic field created by a wire)
for (int i=0;i<Th.nt;i++)
{ for (int j=0; j <3; j++)
ff<<Th[i][j].x << " "<< Th[i][j].y<< " "<<Bx[][Vh(i,j)]<< " "<<By[][Vh(i,j)]<<endl;
ff<<Th[i][0].x << " "<< Th[i][0].y<< " “<<Bx[][Vh(i,0)]<< " “<<By[][Vh(i,0)]<<”\n\n\n”;
}
savemesh(Th,“wire.msh”);

HTH
Denis

Thanks for your help. I think my code in FreeFem++ was ok. It was just the way to display it with opencv or matplotlib that was wrong.
Sorry for that

Dear user
just a small remark,

freefem dowhen the point is outside the domain depend of the finite element.

2 case if the finite element is continuous then the we get a continuous prolongement
and the finite is discontinuous then the prolongement is zero.

In all your case you have a hole in the mesh so the finite element is not clearly difined

Remark: you can use the fonction chi(Th)(xx,yy) to know if the point (xx,yy) is in the mesh Th.