FlowSolver: wrong results in 1-dim pipe flow problem

Clearly defined bug reports and their fixes
Post Reply
carsten_broese
Posts: 4
Joined: 15 May 2022, 17:05
Antispam: Yes

FlowSolver: wrong results in 1-dim pipe flow problem

Post by carsten_broese »

Dear ElmerTeam,

I want to model a 1-dimensional pipe flow with FlowSolver and get wrong results for the output velocity.
Model Options for FlowSolver are: Stokes flow, incompressible, vanishing viscosity, no body forces.

The pipe is of rectangular cross section and modeled in length (x-direction) by one 8-node hexahedron (808) element.
The side at x=0 is Input and the side at x=2 is Output. The flow through all other sides is forbidden in normal direction.
At the Input side, the water flows with velocity v0 into the pipe.
At the Output side, Elmer calculates a flow with velocity v0 also INTO the pipe (see e.g. the velocity vector in output file flow/meshout.ep).
Physically, this is wrong, the flow at x=2 should be OUT of the pipe.
Therefore, the calculated results by the FlowSolver are wrong (divergence of velocity is not 0).

I use ElmerSolver Version: 9.0 (Rev: 5bc675ac, Compiled: 2022-05-12).
My Gnu C, C++ and Fortran Compilers are version 9.4.0
Operating System is UBUNTU 20.04

Input and Mesh Files are attached. Output Files are flow/meshout.*

Thanks for help,
Carsten Broese
Attachments
problem_description.pdf
Describes Problem graphically
(134.15 KiB) Downloaded 33 times
flow.zip
Mesh, sif and result files
(3.63 KiB) Downloaded 37 times
carsten_broese
Posts: 4
Joined: 15 May 2022, 17:05
Antispam: Yes

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by carsten_broese »

Hi Community again,

An Addendum:
Maybe the problem is due to inconsistency of normal vectors of the surfaces of my mesh.

I made the mesh with GMSH and exported it to Elmer via ElmerGrid:

gmsh -3 -format msh22 flow.geo
ElmerGrid 14 2 flow.msh -removelowdim -bulkorder -autoclean

Thanks,
Carsten Broese
Attachments
flow.geo
GMSH geometry file
(755 Bytes) Downloaded 31 times
kevinarden
Posts: 1417
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by kevinarden »

It is a condition of normals, but not the way the mesh was created. In Elmer the positive normal of a boundary is into the solid it is attached to. So a negative on the right side is out of the solid. This is accounted for in the post processing files, I use vtu and paraview. In paraview it displays and lists positive, since it has been converted to global coordinates.
carsten_broese
Posts: 4
Joined: 15 May 2022, 17:05
Antispam: Yes

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by carsten_broese »

Thank you for replay and for the information about the convention in Elmer, that positive normals shows into the solid.
Today, i checked the calculated normal vectors in Elmer (Module: NavierStokes.F90, Subroutine NavierStokesBoundary) and they where calculated correctly.

Therefore, the problem is not the normal vectors (as i thought yesterday). But it still exists.
I attached the calculated result in Elemer *.ep (flowout.ep).
The last 8 lines and 4 columns are the 3 velocity and 1 pressure component of all 8 nodes of the element.
The last 4 lines in column 1 are the component 1 of outflow velocity and should have opposite sign according to the global coordinate system. This is my problem.
Attachments
flowout.ep
Result file
(1.89 KiB) Downloaded 29 times
kevinarden
Posts: 1417
Joined: 25 Jan 2019, 01:28
Antispam: Yes

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by kevinarden »

I don't know if ep uses the global coordinate system or not. Very few people, if anyone uses ep anymore, I don't. Most users use paraview. I don't know how to use ep to check if it is global coordinate or not.
raback
Site Admin
Posts: 4345
Joined: 22 Aug 2009, 11:57
Antispam: Yes
Location: Espoo, Finland
Contact:

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by raback »

Hi

Did not run the case, just looked at the sif, so this is just guessing...

Now the "GradP Discretization" means that the continuity equation is integrated by parts which allows setting the Dirichlet condition for pressure. Otherwise this cannot be done. Now the downside is that there is a new boundary integral term popping up that should be given: the normal velocity over the outlet. The ideal use of this alternative is when we have either no-flux or pressure conditions. It does not handle well this outlet condition by default because the natural BC for continuity equation is zero flux.

What I would do is to use the standard discretization and just omit trying to set the pressure. The natural BCs does define the pressure implicitely as the outlet is assumed to have zero normal traction.

You could also define Dirichlet conditions for pressure both at the inlet and outlet with the gradp discretization and letting normal velocity follow.

In the n-t system the normal is indeed computed such that it points out of the material. In 3D there are additionally 2 tangent vectors. The 1st one is always in (x,y) plane and the 3rd one is defined by cross product the two. Then "Velocity 1" refers to normal direction and components "2" and "3" to the tangent directions.

If my guess is correct there is per se nothing wrong with the existing code. However, there should definately be some warning or documentation of the dangers of choosing this "gradp discretization".

-Peter
Rich_B
Posts: 136
Joined: 24 Aug 2009, 20:18

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by Rich_B »

Hello,

Attached is an archive with three example folders.

1) flow-vin-pin-gradp contains the a slightly modified version of the original example. BCs are Velocity at the inlet and Pressure at the inlet, with gradp turned on. The geometry is a single 8 node brick/hexahedron, so there are 8 nodes. The post file 'flow-vin-pin-gradp.vtk' shows the results for the eight nodes, velocity: 3, pressure: 1, as follows:

1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
-1.01197022E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
-9.70325151E-001 0.00000000E+000 0.00000000E+000 1.00000000E+000
-1.02905333E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
-1.00327803E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000

One would expect the nodal velocity, Vx to be positive for all eight nodes. Note that the above are nodal values, not normal or tangential components.

2) flow-pin-pout-gradp contains a modified version, with Pressure at the inlet and outlet, without setting the inlet Vx, with gradp turned on. The post file 'flow-pin-pout-gradp.vtk' shows the results for the eight nodes, velocity: 3, pressure: 1, as follows:

2.45650889E+030 0.00000000E+000 0.00000000E+000 1.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 1.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 1.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 1.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 0.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 0.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 0.00000000E+000
2.45650889E+030 0.00000000E+000 0.00000000E+000 0.00000000E+000

One can see the inlet and outlet pressures are as expected, and all 8 Vx are positive.

3) flow-vin-pin-no-gradp contains a third modified version. BCs are Velocity at the inlet and Pressure at the inlet, with gradp turned off. The post file 'flow-vin-pin-no-gradp.vtk' shows the results for the eight nodes, velocity: 3, pressure: 1, as follows:

1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
1.00000000E+000 0.00000000E+000 0.00000000E+000 1.00000000E+000
4.00000000E+013 0.00000000E+000 0.00000000E+000 -9.20000000E-001
4.00000000E+013 0.00000000E+000 0.00000000E+000 -9.20000000E-001
4.00000000E+013 0.00000000E+000 0.00000000E+000 -9.20000000E-001
4.00000000E+013 0.00000000E+000 0.00000000E+000 -9.20000000E-001

One can see that all 8 Vx are positive, and the pressure decreases in the flow direction.

Based on the above, it seems to confirm Peter's observations, as well as confirming the original post.

Rich.
Attachments
1D pipe flow.zip
(18.41 KiB) Downloaded 28 times
Rich_B
Posts: 136
Joined: 24 Aug 2009, 20:18

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by Rich_B »

Hello,

Another observation, the original sif set 'Max Output Level = 3'. Changing the max output level for all three examples, from 3 to 5, and then to 32, did not seem to show any warnings or notices about setting Vin, Pin, with gradp turned on.

Rich.
Rich_B
Posts: 136
Joined: 24 Aug 2009, 20:18

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by Rich_B »

Hello,

For completeness, I ran the fourth combination, setting Pin, Pout, but without turning on gradp.

Velocity: 3, Pressure: 1
-1.66039629E+014 0.00000000E+000 0.00000000E+000 1.00000000E+000
-1.66039629E+014 0.00000000E+000 0.00000000E+000 1.00000000E+000
-1.66039629E+014 0.00000000E+000 0.00000000E+000 1.00000000E+000
-1.66039629E+014 0.00000000E+000 0.00000000E+000 1.00000000E+000
3.33960371E+014 0.00000000E+000 0.00000000E+000 0.00000000E+000
3.33960371E+014 0.00000000E+000 0.00000000E+000 0.00000000E+000
3.33960371E+014 0.00000000E+000 0.00000000E+000 0.00000000E+000
3.33960371E+014 0.00000000E+000 0.00000000E+000 0.00000000E+000

So the pressure is as expected, but the 8 values of Vx are not all positive.

Rich.
carsten_broese
Posts: 4
Joined: 15 May 2022, 17:05
Antispam: Yes

Re: FlowSolver: wrong results in 1-dim pipe flow problem

Post by carsten_broese »

@Peter and @Rich, thank you all for your comments.
I will try out the examples tomorrow and think about Peters observation.

-Carsten
Post Reply