As of writing this, I believe the paraview post files plot the higher-order solutions using their version of cells where the solution points within an element is uniformly distributed. I am not sure, but I think because of this I see some Gibbs phenomenon in my contour plots (it doesn’t go away with clean to grid filter etc). Is there a simple way of plotting the solution at the original solution points (Gauss-Legendre or Gauss-Legendre-Labatto) ?
By default pyfr export will export in the high order .vtu format as you suggest, where the order of the exported polynomials will be the same as the solution polynomial order. You can over rider this order with the -k flag. Alternatively you can export a regular linear .vtu format using the -d flag (details are here User Guide — Documentation). In both cases the solution polynomials are potentially evaluated at new sets of points, but I am not clear why this would lead to Gibbs phenomena?
Equispaced points are not the source of any Gibbs like issues.
Regards, Freddie.