Strange issue: no vector display when plotting

Hello everyone !

I haven’t been able to find a topic that deals with this weird issue, so here it is:
I’m trying to solve a 3D Navier-Stokes problem using IPCS ( Incremental Pressure Correction Scheme ). It is the same core principle as in Fenics tutorial ft07.

Here is my method:
-First, I create a STEP file with Onshape (here, a cylinder)
-I mesh the STEP file and tag all the boundaries in GMSH and then convert it to XDMF.
-Finally I import the XDMF file in Fenics, compute Navier-Stokes solution and plot the results with Matplotlib and Paraview.

Here is my issue:
I can see no arrows to represent the velocity field in matplotlib, only points. When plotting in Paraview arrows seems to be unit vectors on the x direction, which is very unrealistic.
However, the norm values seems quite relevant: zero at the cylinder’s wall and maximum velocity near the center.

Any idea on why this is happening ?

(dowload the github repository)
Minimal working example:

Script named that is needed to import the create_mesh function:

.msh file needed:

Since you are only applying constant pressure gradient in x-direction, I think it is not that strange that the velocity vector only has x component. One thing I am curious is that I don’t see any arrows inside the mesh but only on the boundary parts of the mesh. How do you plot in Paraview ?

Hi @Kei_Yamamoto !

Thanks for your answer.Tell me if i’m wrong, but the pressure gradient is only along the z-axis, with a pression of 18.6 at the inflow and 0 at the outflow. That’s why I don’t understand why the velocity field is along the x-axis.
For plotting in paraview I just open the file velocity_field.xdmf with Xdmf3-ReaderS then select 3D glyphs and coloring : f_28

here is a screenshot:

Okay, I misunderstood your Paraview plot and it is indeed strange. I have done similar simulation and the variational formulation looks no problem to me. One difference is that I have never used keep_diagonal=True and A1.ident_zeros(). I don’t really know what they are doing but could you tell me why they are needed ??

Also, maybe you already know it but
this might be helpful.

before building this code I made a 2D one that wasn’t working without those two lines :
keep_diagonal=True and A1.ident_zeros()

Following what you said I just tried without A1.ident_zeros() and turning keep_diagonal to False. It makes no difference.

Actually i knew about this tutorial (Solving PDEs in Python - <br> The FEniCS Tutorial Volume I) and used it to implement the IPCS scheme.

Actually what is really weird is that the two-dimensionnal simulation and plotting works fine. I think i did something wrong when switching from 2D to 3D

Hi gpr,

I think I figured out why it does not work.

When you create triangle_mesh on line 38,39 , you do not set prune_z=true. This needs to be done to create 2D mesh that consists of boundary part.

You should modify your code as follows.

triangle_mesh = create_mesh(mesh3D_from_msh, “triangle”, prune_z=True)
meshio.write("facet_mesh3D_from_msh.xdmf", triangle_mesh)
1 Like

Hi @Kei_Yamamoto ! I don’t think the z coordinate needs to be pruned. We still need it for the surfaces as they have a z-coordinate (they aren’t in the x,y plane). I tried setting prune_z to true and it doesn’t make any difference. Still the idea was relevant and i thank you for your suggestion :slight_smile:

Eventually i found the solution. The problem lied in how i applied the filters in Paraview . It didn’t had anything to do with the Fenics code.


Out of curiosity, what was the filter problem & resolution (in case someone modelling ft07 runs into the same issue)?

At the beginning I only selected “glyph” in the top menu without applying any filters.
It indeed displayed glyphs, but the orientation of the vector field wasn’t taken into account.
In order to display properly a vector field, one must select a glyph filter and enable scaling by magnitude. Then you have to play a bit with the scale factor for the arrows to appear.

This might seem trivial for a regular paraview user, but i was a complete beginner so i got a bit confused :sweat_smile: