[petsc-users] MatSetNullSpace results in nan with MUMPS

[Matthew Knepley]

On Sat, Apr 2, 2022 at 8:59 PM Bhargav Subramanya <
bhargav.subramanya at kaust.edu.sa> wrote:

> Dear All,
>
> I am trying to solve Ax = b in parallel using MUMPS, where x is composed
> of velocity, pressure, and temperature. There is a null space due to the
> homogeneous Neumann pressure boundary conditions.
>
> 1. Without explicitly removing the null space, the solution is bounded and
> robust, and I am able to iterate in time, although I need to check if the
> solution is physically correct. I am not sure if MUMPS is able to
> detect the null space in this case.
>

MUMPS should fail with a 0 pivot. I do not understand this.


> 2. However, when I explicitly use MatSetNullSpace as shown below, the
> solution blows up in 2-3 temporal iterations, resulting in nan. Could you
> please share your thoughts regarding this problem?
>

Did you check that this is a nullspace vector using
https://petsc.org/main/docs/manualpages/Mat/MatNullSpaceTest.html ?

  Thanks,

     Matt


> Thanks,
> Bhargav
>
> ---------------
> Output for -pc_type svd -pc_svd_monitor
> SVD: condition number 3.159257564649e+16, 1 of 1080 singular values are
> (nearly) zero
>     SVD: smallest singular values: 3.693890048955e-17 2.902795363876e-10
> 2.902795486255e-10 2.902796682780e-10 7.614221123356e-10
>     SVD: largest singular values : 1.166981958972e+00 1.166981958972e+00
> 1.166995008014e+00 1.166995008014e+00 1.166995008014e+00
> -----------------
>
> PetscErrorCode *SetNullVector*(Mat MatGlobalLinear,Mesh mesh,Equations
> eq){
>
> Vec NullVector;
>
> PetscInt numNullVectors,disp,i,j;
>
> PetscScalar constantValue = 1.0,nrm;
>
> MatNullSpace NullSp = NULL;
>
> PetscErrorCode ierr;
>
>
> PetscFunctionBeginUser;
>
> numNullVectors = 1;
>
> disp = mesh.numLayers*mesh.numEdges + (eq.pressure-1)*mesh.numCells;
>
> ierr = *MatCreateVecs*(MatGlobalLinear,&NullVector,NULL);CHKERRQ(ierr);
>
> *for*(i=0; i<mesh.numLayers; i++){
>
> *for*(j=0; j<mesh.numCells; j++){
>
> ierr = VecSetValue(NullVector,j + i*mesh.numCells + disp,constantValue,
> *INSERT_VALUES*);CHKERRQ(ierr);
>
> }
>
> }
>
> ierr = *VecAssemblyBegin*(NullVector);CHKERRQ(ierr);
>
> ierr = *VecAssemblyEnd*(NullVector);CHKERRQ(ierr);
>
> ierr = *VecNorm*(NullVector,*NORM_2*,&nrm);CHKERRQ(ierr);
>
> ierr = *VecScale*(NullVector,1.0/nrm);CHKERRQ(ierr);
>
> /* set null space */
>
> ierr = *MatNullSpaceCreate*(MPI_COMM_WORLD,*PETSC_FALSE*
> ,numNullVectors,&NullVector,&NullSp);CHKERRQ(ierr);
>
> ierr = *MatSetNullSpace*(MatGlobalLinear,NullSp);CHKERRQ(ierr);
>
> /* clean up */
>
> ierr = *MatNullSpaceDestroy*(&NullSp);CHKERRQ(ierr);
>
> ierr = *VecDestroy*(&NullVector);CHKERRQ(ierr);
>
> PetscFunctionReturn(0);
>
> }
> -----------------
>
>
> ------------------------------
> This message and its contents, including attachments are intended solely
> for the original recipient. If you are not the intended recipient or have
> received this message in error, please notify me immediately and delete
> this message from your computer system. Any unauthorized use or
> distribution is prohibited. Please consider the environment before printing
> this email.



-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20220403/4834a7ec/attachment-0001.html>