Re: [Getfem-users] using matplotlib

[jeff witz]

Hello,

You can find below the discussion with Luis Saavedra about plotting results with matplotlib.
To summarize, the main issue is for extracting datas that are compatible the the pid i.e the mesh connectivity table. A delaunay triangulation can't handle the hole in a plate for example, and once the geometry is lost it's not possible to handle it correctly. 
It must be possible to sort "u = md.variable('u')" to get a vector field compatible with mesh.
Does any one know a way to sort it out ?
if it's not possible I think it could be a good idea to add a methos that allows this kind of transformation in the model class.

Furthermore, I observed that the triangulation is not the same for the slice use in case of export and even if I can understand it in case of quadratique shape functions, it's difficult to justify for linear ones.

Best Regards.




>at the moment the only solution that I can think is to remove leftover triangles as does the example of matplotlib...

2011/3/17 jeff witz <address@hidden>

Thanks,

I will send a post during the day.

Doy you find any solution for having good connectivity for the plot I gave you ?

My question is how to rearrrange the solution vector to the pid order ?

If you don't know I will post on the list

On 17 March 2011 07:55, Luis Saavedra <address@hidden> wrote:

Hi jeff,

you're right, the example that you sent me finishes with a segment fault, I still don't find the reason for this, please send it to the list,

regards,
Luis

2011/3/7 jeff witz <address@hidden>

Thanks for this response.

First even If completly do agree about the KISS principle, I think there is a problem.
If you look to the geometry file I give you, it's not possible to recover the mesh using Delaunay triangulation. The arc circle in the 4 corners leads to unwanted triangles. So I think this question is not solved yet. Could you please try with the plate geometry I give you ?

I've completed the script you gave me with one of the most simple example I can find. It's what I give to the student to do manually.
So I send you this script.

I've noticed that there is an error (sementation fault) when I try to change from Dirichlet to Neuman (line 60) conditions on the right side of the mesh, I don't understand why, so If you can reproduce  this behavior I will send on the support list.

So I send you an advanced example of plot, in order to see that for basics curves this method is quite interesting.

 

On 7 March 2011 06:21, Luis Saavedra <address@hidden> wrote:

Dear jeff,

sorry the delay in response, I've attached a simple example where I use matplotlib.tri to plot a solution obtained with getfem++,

regards,
Luis

PD: http://en.wikipedia.org/wiki/KISS_principle

2011/3/4 jeff witz <address@hidden>

Hi,
Could you please just help me in sorting the solution vector I think I just need this to finish my course ...
 
Regards.

On 4 March 2011 09:01, jeff witz <address@hidden> wrote:

Hello,

On 4 March 2011 03:42, Luis Saavedra <address@hidden> wrote:

Hi jeff,

three questions:

 * are using the svn version?

 

For now I don't use the svn version but I can easily switch from 4.4.1 to svn.
 

 * if you had the correct values, how they would be plotted?

First, you need at least the version 1.0.1 of matplotlib. There are very instructive example on the gallery :

http://matplotlib.sourceforge.net/examples/pylab_examples/tripcolor_demo.html

Beware if you already have a older version of matplotlb it's quite complex to handle (I remove the packaged version)
I regret that all the linux distributions only handle the 0.99 version.

For now I just want for my course to use a simple delaunay triangulation, but there is much more that can be done.

 I'm asking you because if I know what you want, it might be a good idea to integrate this into the interface without having to write the code yourself,

It would be great !

But I think that for teaching purposes building connectivity, inverse connectivity and node table is interesting.

I send you the source code that allow to plot the displacement on the mesh, but for now with wrong values.
I still believe it's quite strange that the displacements are not ordinate as the coordinate  pid, is this the normal behaviour ? 
You need to have gmsh installed.
If you are interested in I also write a little script in python that allows to get the mesh and results and work with. But all the faces are independent in this format so it's not fem behavior ....

You can find below what I think to be the good way to integrate plot in python using matplotlib.

If you look to this example :
http://matplotlib.sourceforge.net/examples/pylab_examples/shading_example.html
you can see that ray-tracing is available in matplotlib.
My point is that it is possible to have the perfect visualization tool imshow function. The main issue with complex element is that there is a need to create a lot of submesh (as X-FEM) to correctly display the picture, I think this strategy (using in paraview) is not the good one. It is much more versatile to decompose the plot in two parts.
First the the FEM value on a selected griddata (projection of the delaunay triangulation to a regular grid) that allows you to get a pixel description,
Then to plot the mesh, the main difficulty  lie on the fact that there is a need to create a Z-buffer to plot only the visible element of the mesh.

By this way, it can be quite simple to produce light but very professionnal plot (perfect for article), indeed the color variation can be stored in png (the most efficient format for this purpose) and the visible mesh and axes in vectorial (svg,pdf,eps) in order to have good results.


 

regards,
Luis

2011/3/3 jeff witz <address@hidden>

Hello,

I almost succeed in using matplotlib for plotting 2 dimensionnal linear triangle !

So I have a two dimensional NonLinearelasticity problem, so I do the computation and extract the solution

u = md.variable('u')

In order to be able to plot solution in matplotlib 1.0.1 I need node table and element table. So I build this table :


NbNod=m.pts().shape[1]
NbElem = m.cvid().shape[0]
NbNodElem=m.pid_from_cvid(0)[0].shape[0]
ElemTab = np.zeros((NbElem,NbNodElem))
for i in range(NbElem):
    for j in range(NbNodElem):
        ElemTab[i,j]=m.pid_in_cvids(i)[j]

NodTab=np.transpose(m.pts())

The main issue remain in the fact that I don't understand how the u vector is sorted.
So when I plot the results the values don't correspond to the computed solution. Could you explain me how I can sort it.

Logically this vector is 324 : 2* NbNod

I try to use mf.basic_dof_nodes(), but I don't really understand the way output vector

mf.basic_dof_nodes().shape
Out[100]: (2, 324)

which is equal to (2,2*NbNod)

Could you help me to solve this out ?

I think it would be great to be able to treat the data without leaving ipython.

I just want to plot a simple tripcolor (pcolor in matlab) with a quiver over it. It will be quite simple once I have the good u vector

I hope I can solve this before the Course I have to give Monday to show the student that everything can be done simply in python.

Best regards.




--
Jean-François WITZ

_______________________________________________
Getfem-users mailing list
address@hidden
https://mail.gna.org/listinfo/getfem-users

--
Jean-François WITZ

--
Jean-François WITZ

--
Jean-François WITZ

--
Jean-François WITZ

Example.zip
Description: Zip archive