[petsc-users] Rigid body nullspace for Stokes operator

[Jed Brown]

Yes to 6 null vectors in 3D, 3 null vectors in 2D. The center of mass does not need to be identified because you can algebraically orthogonalize (lines 411-420 here).

https://urldefense.us/v3/__https://gitlab.com/petsc/petsc/-/blob/main/src/dm/impls/plex/plexfem.c?ref_type=heads*L377-425__;Iw!!G_uCfscf7eWS!bc0A2Vc5CGCuDQuOVZRyBUNTEbK3NZvHNmZUECNwUBwJPj0Pvv0CW2rGfEzeSxVXdWkjCE-D5jl2sonSQN8$ 

See also this implementation with raw coordinates. GAMG orthogonalizes within each aggregate (in a later phase of the algorithm) so global orthogonalization is not necessary.

https://urldefense.us/v3/__https://gitlab.com/petsc/petsc/-/blob/main/src/ksp/pc/impls/gamg/agg.c?ref_type=heads*L387__;Iw!!G_uCfscf7eWS!bc0A2Vc5CGCuDQuOVZRyBUNTEbK3NZvHNmZUECNwUBwJPj0Pvv0CW2rGfEzeSxVXdWkjCE-D5jl20dv2S0M$ 

Amneet Bhalla <mail2amneet at gmail.com> writes:

> I think the nullspace for the velocity operator is of the form
>
> vnull = U + ω × r
> in which U is a rigid body velocity and ω is the rigid body rotational
> velocity, and r is the radius vector from the center of mass. I believe I
> need to construct 6 nullspace vectors in 3D and 3 nullspace vectors in 2D.
> Sounds correct? Also does the center of mass coordinates matter when
> defining r?
>
> On Wed, Oct 30, 2024 at 7:53 AM Amneet Bhalla <mail2amneet at gmail.com> wrote:
>
>> @Mark: Is there some document/paper that I can follow to check the algebra
>> of these zero eigenvectors/null space modes?
>>
>> @Jed : We use a projection method preconditioner to solve the coupled
>> velocity pressure system as described here (
>> https://urldefense.us/v3/__https://www.sciencedirect.com/science/article/pii/S0021999123004205__;!!G_uCfscf7eWS!bc0A2Vc5CGCuDQuOVZRyBUNTEbK3NZvHNmZUECNwUBwJPj0Pvv0CW2rGfEzeSxVXdWkjCE-D5jl2BQDymOU$ ). It
>> is an approximation of the Schur complement. As a part of projection
>> preconditioner, we need to solve just the momentum equation separately
>> without considering the pressure part.
>>
>> On Tue, Oct 29, 2024 at 8:03 PM Jed Brown <jed at jedbrown.org> wrote:
>>
>>> And to be clear, we recommend using fieldsplit Schur to separate the
>>> pressure and velocity part (there are many examples). Applying AMG directly
>>> to the saddle point problem will not be a good solver because the
>>> heuristics assume positivity and do not preserve inf-sup stability (nor do
>>> standard smoothers).
>>>
>>> Mark Adams <mfadams at lbl.gov> writes:
>>>
>>> > This is linear elasticity and there are 6 "null" vectors (if you removed
>>> > Dirichlet boundary conditions these are eigenvectors with zero
>>> eigenvalue):
>>> > 3 translations, x, y, z, and three rotatiions xx, yy ,zz.
>>> > x = (1,0,0,1,0,0,1,0 ...)
>>> > and xx is something like (0, z_1, -y_1, 0, z_2, -y_2, ...) where z_1 is
>>> the
>>> > z coordinate of the first vertex, etc.
>>> >
>>> > Mark
>>> >
>>> > On Tue, Oct 29, 2024 at 3:47 PM Amneet Bhalla <mail2amneet at gmail.com>
>>> wrote:
>>> >
>>> >> Hi Mark,
>>> >>
>>> >> Thanks! I am not sure how to construct null space and zero energy modes
>>> >> manually for this operator. Is there some theory or documentation I can
>>> >> follow to figure out what the null space and zero energy modes look
>>> like
>>> >> for this operator? Once I know what these are in symbolic form, I
>>> think I
>>> >> should be able to construct them manually.
>>> >>
>>> >> Best,
>>> >> --Amneet
>>> >>
>>> >> On Tue, Oct 29, 2024 at 7:35 AM Mark Adams <mfadams at lbl.gov> wrote:
>>> >>
>>> >>> Oh my mistake. You are using staggered grids. So you don't have a
>>> block
>>> >>> size that hypre would use for the "nodal" methods.
>>> >>> I'm not sure what you are doing exactly, but try hypre and you could
>>> >>> create the null space, zero energy modes, manually, attach to the
>>> matrix
>>> >>> and try GAMG.
>>> >>> You can run with '-info :pc' and grep on GAMG to see if GAMG is
>>> picking
>>> >>> the null space up (send this output if you can't figure it out).
>>> >>>
>>> >>> Thanks,
>>> >>> Mark
>>> >>>
>>> >>> On Tue, Oct 29, 2024 at 9:28 AM Mark Adams <mfadams at lbl.gov> wrote:
>>> >>>
>>> >>>> This coordinate interface is just a shortcut for vertex based
>>> >>>> discretizations with 3 dof per vertex, etc. (maybe works in 2D).
>>> >>>> You will need to construct the null space vectors manually and
>>> attach it
>>> >>>> to the matrix. Used by GAMG.
>>> >>>>
>>> >>>> Note, for hypre you want to use the "nodal" options and it does not
>>> use
>>> >>>> these null space vectors. That is probably the way you want to go.
>>> >>>> eg: -pc_hypre_boomeramg_nodal_coarsen
>>> >>>>
>>> >>>> I would run with hypre boomerang and -help and grep on nodal to see
>>> all
>>> >>>> the "nodal" options and use them.
>>> >>>>
>>> >>>> Thanks,
>>> >>>> Mark
>>> >>>>
>>> >>>>
>>> >>>> On Mon, Oct 28, 2024 at 8:06 PM Amneet Bhalla <mail2amneet at gmail.com
>>> >
>>> >>>> wrote:
>>> >>>>
>>> >>>>> Hi Folks,
>>> >>>>>
>>> >>>>> I am trying to solve the momentum equation in a projection
>>> >>>>> preconditioner using GAMG or Hypre solver. The equation looks like
>>> for
>>> >>>>> velocity variable *v* looks like:
>>> >>>>>
>>> >>>>>
>>> >>>>>  [image: Screenshot 2024-10-28 at 4.15.17 PM.png]
>>> >>>>>
>>> >>>>> Here, μ is spatially varying dynamic viscosity and λ is spatially
>>> >>>>> varying bulk viscosity. I understand that I need to specify rigid
>>> body
>>> >>>>> nullspace modes to the multigrid solver in order to accelerate its
>>> >>>>> convergence. Looking into this routine
>>> MatNullSpaceCreateRigidBody() (
>>> >>>>>
>>> https://urldefense.us/v3/__https://petsc.org/release/manualpages/Mat/MatNullSpaceCreateRigidBody/__;!!G_uCfscf7eWS!bVR6duCoDqPhZrWS-sm1c5qxsFPjZMhdT86AqLpPzWgVy5qoRhd4_Jue2LJOIS6LRrtV2cHGrqger1Yvb-Y5f-0$
>>> >>>>> <
>>> https://urldefense.us/v3/__https://petsc.org/release/manualpages/Mat/MatNullSpaceCreateRigidBody/__;!!G_uCfscf7eWS!eKqgIJjCdMzIU76f7X65AmGxrU_-lC7W02BMWafJ77DNf_IuQk6O1X3qU1x9Ez8NJ20vZEL-mF6T1yNmDnwv0eWa2w$
>>> >),
>>> >>>>> I see that I need to provide the coordinates of each node. I am
>>> using
>>> >>>>> staggered grid discretization. Do I need to provide coordinates of
>>> >>>>> staggered grid locations?
>>> >>>>>
>>> >>>>> Thanks,
>>> >>>>> --
>>> >>>>> --Amneet
>>> >>>>>
>>> >>>>>
>>> >>>>>
>>> >>>>>
>>> >>
>>> >> --
>>> >> --Amneet
>>> >>
>>> >>
>>> >>
>>> >>
>>>
>>
>>
>> --
>> --Amneet
>>
>>
>>
>>
>
> -- 
> --Amneet