[matdehaast]

Good day Ladies and Gentlemen

I am having an issue with the A Matrix created through the pressure poisson equation. I am finding that the condition number of this matrix can be really large and for large aspect ratios on quad meshes becomes way to large to solve. For instance if I have elements that have volume ranging from O(1e-9) to O(1e-3) I cant solve the system.

Has anyone got a better idea on improving the condition number of this system. I have tried adding values to the diagonal, which does help in some instances but not the above mentioned problem.

Thanks
Matt

[FMDenaro]

Originally Posted by matdehaast

Good day Ladies and Gentlemen

I am having an issue with the A Matrix created through the pressure poisson equation. I am finding that the condition number of this matrix can be really large and for large aspect ratios on quad meshes becomes way to large to solve. For instance if I have elements that have volume ranging from O(1e-9) to O(1e-3) I cant solve the system.

Has anyone got a better idea on improving the condition number of this system. I have tried adding values to the diagonal, which does help in some instances but not the above mentioned problem.

Thanks
Matt

Theoretically, the pressure equation with Neumann BC.s is singular

[matdehaast]

Thank you I understand that but it has been and is been used in this context among many codes. How do they resolve this issue?

[FMDenaro]

you can search for similar posts...
the matrix is singular but if you fulfill the compatibility condition with your Neumann BC.s, the system has infinite solutions (that means a solution apart a function of time)

[mprinkey]

As Prof Denaro stated, the Neumann BCs on the pressure equation for incompressible flow is formally singular and admits solutions of the form P + C where P is a pressure field and C is an arbitrary scalar uniform constant. So, you have to lock that value of C down somehow.

OpenFOAM solvers generally set one cell pressure value if there are no pressure inlets or outlets. You can see something like this in almost every PEqn.H file:

pEqn.setReference(pRefCell, pRefValue);

That makes the diagonal entry in A corresponding to rRefCell equal to 1.0, all off-diagonal entries on that row equal to 0.0, and the b vector for that row equal to pRefValue.

Some relaxation schemes and non-stationary solvers can get away without doing this...they naturally relax to a solution that doesn't drift with a gauge value, but you cannot necessarily count on that.