Dear all,

I’m a newbie at Julia. Probably I’m doing something wrong and/or inefficiently: I was trying to replicate one of the examples to generate my own linear elasticity solver. Starting with a very simple example (following incompressible elasticity example ).

```
function create_cook_grid(nx, ny)
corners = [Vec{2}((0.0, 0.0)),
Vec{2}((48.0, 44.0)),
Vec{2}((48.0, 60.0)),
Vec{2}((0.0, 44.0))]
grid = generate_grid(Triangle, (nx, ny), corners);
# facesets for boundary conditions
addfaceset!(grid, "clamped", x -> norm(x[1]) ≈ 0.0);
addfaceset!(grid, "traction", x -> norm(x[1]) ≈ 48.0);
return grid
end;
```

and copying assembly functions, I can found a solution using P1 and P2 finite elements. In VTK looks great. However, when I try to get nodal coordinates and its connectivity:

```
grid.nodes; # has 99-element Array{Node{2,Float64},1}
grid.cells; # has 160-element Array{Cell{2,3,3},1}
```

while the length of the 2d vector solution, called `u`

, is 5202. The latter coincides with the info provided by Paraview: num of cells = 5000, num of points = 2601.

Trying to collect the nodal coordinates (certainly this looks sub-optimal):

```
nodes = grid.nodes
elems = grid.cells;
nnodes = length(nodes); # << 99-element Array{Node{2,Float64},1}
nelems = length(elems); # << 160-element Array{Cell{2,3,3},1}
coords = Array{Float64,2}(undef,nnodes,2)
cells = Array{Int,2}(undef,nelems,3)
for i in 1:nnodes
coords[i,:] = [nodes[i].x[1],nodes[i].x[2]];
end
for i in 1:nelems
cells[i,:] = [elems[i].nodes[1],elems[i].nodes[2],elems[i].nodes[3]]
end
```

and then

```
using PyPlot
PyPlot.plot_trisurf(coords[:,1],coords[:,2],triangles=cells)
```

shows a nice mesh, but coarser than that shown in Paraview.

I’m completely lost. How can I get the *real* nodal coordinates and connectivity?

Another question: Had someone import the VTK library in Julia? I’m thinking of having something similar to Mayavi in Julia.

Thanks.

J