# Supporting higher order elements in iovtk
**URL:** https://community.freefem.org/t/supporting-higher-order-elements-in-iovtk/1302
**Category:** Feature Request
**Created:** [November 4, 2021, 8:07pm UTC](https://community.freefem.org/t/supporting-higher-order-elements-in-iovtk/1302 "2021-11-04T20:07:29Z")
**Posts on this page:** 2
**Page:** 1
### Author:  [@ivan](https://community.freefem.org/u/ivan)
#### Post date: [November 4, 2021, 8:07pm UTC](https://community.freefem.org/t/supporting-higher-order-elements-in-iovtk/1302/1 "2021-11-04T20:07:29Z")
For some time now Paraview has supported arbitrary order Lagrange finite elements, see
> https://blog.kitware.com/modeling-arbitrary-order-lagrange-finite-elements-in-the-visualization-toolkit/
Currently the `iovtk` module only supports P0 and P1 elements. On the other hand deal.II and mfem already support this:
- deal.ii: [Output and visualization of high-order solutions/meshes](https://github.com/dealii/dealii/issues/6931)
- mfem: [Mesh reader support for VTK high order elements [feature/vtk\_high\_order]](https://github.com/mfem/mfem/pull/975)
It would be nice to see this on the FreeFem roadmap too.
---
### Author:  [@fotios.kasolis](https://community.freefem.org/u/fotios.kasolis)
#### Post date: [August 4, 2022, 9:45am UTC](https://community.freefem.org/t/supporting-higher-order-elements-in-iovtk/1302/2 "2022-08-04T09:45:43Z")
Indeed this would be nice, but on the other hand it is trivial to hack in FreeFEM itself and then save the solution. Whenever better graphics are needed (or when using more exotic elements than Lagrangian), I do an h-refinement and the interpolation in FreeFEM itself, with `trunk` and `u=u`, like in the snippet below. Maybe not a real solution to your feature request, but it actually works, and you can make a macro out of this to get a certain degree of automation.
```auto
load "iovtk"
string tag = "u";
int[int] ord = [1];
border a001(t = 0.0, 2.0*pi){ x = cos(t); y = sin(t);}
mesh T = buildmesh(a001(10));
fespace S(T, P2);
S u, v;
solve Laplace(u, v) = int2d(T)( dx(u)*dx(v) + dy(u)*dy(v) ) + on(a001, u=x*y);
plot(u, fill=true, nbiso=64, wait=true);
savevtk("u-wo.vtu", T, u, dataname=tag, order=ord);
T = trunc(T, true, split=3);
u = u;
plot(u, fill=true, nbiso=64, wait=true);
savevtk("u-w.vtu", T, u, dataname=tag, order=ord);
```
 
Actually I use the same approach when I need to see how vector fields look.
```auto
mesh T = square(1, 1);
fespace S(T, P0);
fespace V(T, RT0Ortho);
V [ux, uy] = [0, 0];
ux[][2] = 1;
T = trunc(T, true, split = 50);
[ux, uy] = [ux, uy];
S unorm = sqrt(ux*ux + uy*uy);
plot(unorm, [ux, uy], fill=true, nbiso=128, wait=true);
```
