# Singular Linear System - Pressure Poisson Equation

 Register Blogs Members List Search Today's Posts Mark Forums Read

 September 10, 2008, 14:06 Singular Linear System - Pressure Poisson Equation #1 Allan Guest   Posts: n/a Hello, I'm trying to solve the pressure Poisson equation with FULL Neumann boundary condition. Unfortunately, it leads to the sparse linear system equation Ax = b, where A is singular (because the Neumann BC in all boundaries). Well, I've tried Gauss-Seidel and it works very well, but it is very very slow. Now, I'm trying to solve this equation with GMRES or BiCGStab with the library gmm++, but it seems that all Krylov methods get blocked with singular matrixes and do not converge. I'd like to know how can I solve the pressure Poisson equation with GMRES? Is there any library which implements these methods with possibility to solve singular linear systems?

 September 10, 2008, 15:57 Re: Singular Linear System - Pressure Poisson Equa #2 Jed Guest   Posts: n/a You need to remove the null space. The best way is to tell the Krylov solver about this (see KSPSetNullSpace(...) in PETSc) but you can also fix the pressure at one point. That is, pick any row of your matrix and set it to a row of the identity and put 0 on the right hand side. This should also make the matrix nonsingular, but it disrupts the spectrum so may effect convergence rates (highly preconditioner dependent).

 September 10, 2008, 16:10 Re: Singular Linear System - Pressure Poisson Equa #3 Allan Guest   Posts: n/a What can you talk about the difference of PETSc and AztecOO in Trilinos package? It seems that PETSc has more options when solving a linear system.

 September 10, 2008, 18:07 Re: Singular Linear System - Pressure Poisson Equa #4 Jed Guest   Posts: n/a I'm not very familiar with AztecOO but there should be a way to remove the null space. I appreciate the PETSc philosophy of making everything available from the command line. Avoid at all costs writing code that you have to recompile to change solvers or preconditioners. I think PETSc is somewhat more comprehensive and it provides many third-party preconditioners (including ML which is arguably the best part of Trilinos, and Hypre from LLNL) and direct solvers all of which are tunable from the command line. Also PETSc cleanly separates the preconditioning matrix from the Jacobian. This is especially useful for me because I have a fast way to apply the Jacobian matrix-free, but use different (sparser) matrices for preconditioning. This is possible with AztecOO, but less clean in my opinion. The matrix and vector interface is certainly much richer in PETSc which makes it possible to write solvers which are more generic. my \$.02

 September 10, 2008, 20:17 Re: Singular Linear System - Pressure Poisson Equa #5 Q Guest   Posts: n/a Allan I would image that the PPE formulation would have neumann BC for closed cell problems. but what sorts of problems are you designing this code for where you anticipate all the pressure boundaries to be gradient (Neumann)? Also, I haven't experience with the CG solver you describe but shouldn't it solve for N BC if the Neumann condition is upheld? curious. Regards, D

 September 10, 2008, 23:00 Re: Singular Linear System - Pressure Poisson Equa #6 Allan Guest   Posts: n/a The case is that in the development of the numerical procedure for solving the incompressible Navier-Stokes equation, you end up with the Neumann BC for all boundaries in a Poisson equation for a variable phi (which is not the pressure, but is related to it). laplacian(phi) = divergence(V*)/dt, interior nabla(phi) = 0, boundary If you wish, I can cite to you the paper where I got it.

 September 10, 2008, 23:06 Re: Singular Linear System - Pressure Poisson Equa #7 Allan Guest   Posts: n/a I've found it in the PETSc's manual: 4.5 Solving Singular Systems Sometimes one is required to solver linear systems that are singular. That is systems with the matrix has a null space. For example, the discretization of the Laplacian operator with Neumann boundary conditions as a null space of the constant functions. PETSc has tools to help solve these systems. First, one must know what the null space is and store it using an orthonormal basis in an array of PETSc Vecs. (The constant functions can be handled separately, since they are such a common case). Create a MatNullSpace object with the command MatNullSpaceCreate(MPI Comm,PetscTruth hasconstants,int dim,Vec *basis,MatNullSpace *nsp); Here dim is the number of vectors in basis and hasconstants indicates if the null space contains the constant functions. (If the null space contains the constant functions you do not need to include it in the basis vectors you provide). One then tells the KSP object you are using what the null space is with the call KSPSetNullSpace(KSP ksp,MatNullSpace nsp); The PETSc solvers will now handle the null space during the solution process. But if one chooses a direct solver (or an incomplete factorization) it may still detect a zero pivot. You can run with the additional options -pc_factor_shift_nonzero or -pc_ factor_shift_nonzero to prevent the zero pivot. A good choice for the damping factor is 1.e-10.

 September 11, 2008, 04:41 Re: Singular Linear System - Pressure Poisson Equa #8 Q Guest   Posts: n/a Hi Allan But if you have an outflow boundary, you would normally set the PPE boundary as fixed Dirichlet. Gradient BC on vixed velocity faces derive P,n = mu Delta(V.n). Thus you would automatically satisfy the Neumann condition. Essentially, by analogy, the diffusion of pressure (via the PPE equation) has access to this fixed boundary thereby eliminating the non-uniqueness I think you are talking about.

 September 11, 2008, 04:46 Re: Singular Linear System - Pressure Poisson Equa #9 Tom Guest   Posts: n/a Actually BiCGstab should still work without fixing the unknown constant in the pressure - in practise it actually converges faster when the constant isn't fixed.

 September 11, 2008, 09:17 Re: Singular Linear System - Pressure Poisson Equa #10 Allan Guest   Posts: n/a Hello Q, In my case, if there is an outflow BC, I still have to set up the Neumann BC in all boundaries for the Pressure Poisson equation. This is because of an enforcement of the BC on the Momentum equation (without pressure). The treatment of the outflow goes in other place. 1º) Calculate u* and v* (without pressure [p] in the momentum eq.) with Crank-Nicholson for viscous term and Adams-Bashforth for advective term. The boundary conditions are u^(n+1) and v^(n+1). 2º) Calculate phi^(n+1) with laplacian(phi) = divergence(u*, v*) / dt. Boundary conditions must be Neumann in all boundaries. 3º) Calculate u^(n+1) and v^(n+1) with the projection eq. V^(n+1) = V* - dt*Grad(phi) You can find all that I've written here in: http://citeseerx.ist.psu.edu/viewdoc...=10.1.1.7.5949 This is an excellent paper.

 August 12, 2010, 08:10 #11 New Member   Shyam Join Date: Apr 2009 Posts: 29 Rep Power: 10 Hi Allan, I am wondering if you were able to solve the problem you had mentioned in the earlier posts. I am trying to use GMRES based solver for solution of pressure poisson equation. Please let me know if you managed to solve the problem. Thanks and Regards, Shyam

 March 17, 2011, 10:39 #12 New Member   Carlos Osorio Join Date: Mar 2011 Posts: 1 Rep Power: 0 The problem maybe is that the "pressure" that enters is not the same as in the boundaries gets out. If you have for example a Laplace equation (laplacian(P) = 0) it says that "not any pressure inside is being generated". The same as a thermic wall (L=1) with no generation (the heat that enters in one side SHOULD BE THE SAME as in the other side). And if you have 0 flux in the left side and 1 in the right, the only Poisson 1D equation that can be solved is d^2 P /dx^2 = 1 (because we have to generate 1 "heat" to receive 1 heat in the right side). Whe i was simulating that "singular system" in 2D i always got Planes wherever the f (right side) was. And the GMRES converged only when the generation and the Neumman BC's was in "equilibrium".

 Thread Tools Display Modes Linear Mode

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post TheBoyce Main CFD Forum 29 November 11, 2011 07:43 haha Main CFD Forum 1 January 27, 2006 16:39 Maciej Matyka Main CFD Forum 9 November 10, 2004 12:30 Dan Moskal Main CFD Forum 0 October 24, 2002 22:02 DS & HB Main CFD Forum 0 January 8, 2000 16:00

All times are GMT -4. The time now is 09:21.