[petsc-users] discontinuous viscosity stokes equation 3D staggered grid

[Dave May]

Yes, the nullspace is important.
There is definitely a pressure null space of constants which needs to be
considered.
I don't believe ONLY using
  -fieldsplit_1_ksp_constant_null_space
is not sufficient.

You need to inform the KSP for the outer syetem (the coupled u-p system)
that there is a null space of constants in the pressure system. This cannot
(as far as I'm aware) be set via command line args. You need to write a
null space removal function which accepts the complete (u,p) vector which
only modifies the pressure dofs.

Take care though when you write this though. Because of the DMDA
formulation you are using, you have introduced dummy/ficticious pressure
dofs on the right/top/front faces of your mesh. Thus your null space
removal function should ignore those dofs when your define the constant
pressure constraint.



Cheers,
  Dave


On 26 July 2013 17:42, Matthew Knepley <knepley at gmail.com> wrote:

> On Fri, Jul 26, 2013 at 10:13 AM, Bishesh Khanal <bisheshkh at gmail.com>wrote:
>
>>
>>
>>
>> On Fri, Jul 26, 2013 at 4:22 PM, Matthew Knepley <knepley at gmail.com>wrote:
>>
>>> On Fri, Jul 26, 2013 at 9:11 AM, Bishesh Khanal <bisheshkh at gmail.com>wrote:
>>>
>>>>
>>>>
>>>>
>>>> On Fri, Jul 26, 2013 at 2:32 PM, Matthew Knepley <knepley at gmail.com>wrote:
>>>>
>>>>> On Fri, Jul 26, 2013 at 7:28 AM, Bishesh Khanal <bisheshkh at gmail.com>wrote:
>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> On Wed, Jul 17, 2013 at 9:48 PM, Jed Brown <jedbrown at mcs.anl.gov>wrote:
>>>>>>
>>>>>>> Bishesh Khanal <bisheshkh at gmail.com> writes:
>>>>>>>
>>>>>>> > Now, I implemented two different approaches, each for both 2D and
>>>>>>> 3D, in
>>>>>>> > MATLAB. It works for the smaller sizes but I have problems solving
>>>>>>> it for
>>>>>>> > the problem size I need (250^3 grid size).
>>>>>>> > I use staggered grid with p on cell centers, and components of v
>>>>>>> on cell
>>>>>>> > faces. Similar split up of K to cell center and faces to account
>>>>>>> for the
>>>>>>> > variable viscosity case)
>>>>>>>
>>>>>>> Okay, you're using a staggered-grid finite difference discretization
>>>>>>> of
>>>>>>> variable-viscosity Stokes.  This is a common problem and I recommend
>>>>>>> starting with PCFieldSplit with Schur complement reduction (make that
>>>>>>> work first, then switch to block preconditioner).
>>>>>>
>>>>>>
>>>>>> Ok, I made my 3D problem work with PCFieldSplit with Schur complement
>>>>>> reduction using the options:
>>>>>> -pc_fieldsplit_type schur -pc_fieldsplit_detect_saddle_point
>>>>>> -fieldsplit_1_ksp_constant_null_space
>>>>>>
>>>>>>
>>>>>> You can use PCLSC or
>>>>>>> (probably better for you), assemble a preconditioning matrix
>>>>>>> containing
>>>>>>> the inverse viscosity in the pressure-pressure block.  This diagonal
>>>>>>> matrix is a spectrally equivalent (or nearly so, depending on
>>>>>>> discretization) approximation of the Schur complement.  The velocity
>>>>>>> block can be solved with algebraic multigrid.  Read the PCFieldSplit
>>>>>>> docs (follow papers as appropriate) and let us know if you get stuck.
>>>>>>>
>>>>>>
>>>>>> Now,  I got a little confused in how exactly to use command line
>>>>>> options to use multigrid for the velocity bock and PCLS for the pressure
>>>>>> block. After going through the manual I tried the following:
>>>>>>
>>>>>
>>>>> You want Algebraic Multigrid
>>>>>
>>>>> -pc_type fieldsplit -pc_fieldsplit_detect_saddle_point
>>>>> -pc_fieldsplit_type schur
>>>>>   -fieldsplit_0_pc_type gamg
>>>>>    -fieldsplit_1_ksp_type fgmres
>>>>>    -fieldsplit_1_ksp_constant_null_space
>>>>>    -fieldsplit_1_ksp_monitor_short
>>>>>    -fieldsplit_1_pc_type lsc
>>>>> -ksp_converged_reason
>>>>>
>>>>> I tried the above set of options but the solution I get seem to be not
>>>> correct. The velocity field I get are quite different than the one I got
>>>> before without using gamg which were the expected one.
>>>> Note: (Also, I had to add one extra option of
>>>> -fieldsplit_1_ksp_gmres_restart 100 , because the fieldsplit_1_ksp residual
>>>> norm did not converge within default 30 iterations before restarting).
>>>>
>>>
>>> These are all iterative solvers. You have to make sure everything
>>> converges.
>>>
>>
>> When I set restart to 100, and do -ksp_monitor, it does converge (for the
>> fieldsplit_1_ksp). Are you saying that in spite of having
>> -ksp_converged_reason in the option and petsc completing the run with the
>> message "Linear solve converged due to CONVERGED_RTOL .." not enough to
>> make sure that everything is converging ? If that is the case what should I
>> do for this particular problem ?
>>
>
> If your outer iteration converges, and you do not like the solution, there
> are usually two possibilities:
>
>   1) Your tolerance is too high, start with it cranked down all the way
> (1e-10) and slowly relax it
>
>   2) You have a null space that you are not accounting for
>
> I have used the MAC scheme with indexing as shown in:  fig 7.5, page 96
>> of:
>> http://books.google.co.uk/books?id=W83gxp165SkC&printsec=frontcover&dq=Introduction+to+Numerical+Geodynamic+Modelling&hl=en&sa=X&ei=v6TmUaP_L4PuOs3agJgE&ved=0CDIQ6AEwAA
>>
>> Thus I have a DM with 4 dof but there are several "ghost values" set to 0.
>> Would this cause any problem when using the multigrid ? (This has worked
>> fine when not using the multigrid.)
>>
>
> I don't know exactly how you have implemented this. These should be rows
> of the identity.
>
>
>>  I do this problem in the tutorial with
>>> a constant viscosity.
>>>
>>
>> Which tutorial are you referring to ? Could you please provide me the
>> link please ?
>>
>
> There are a few on the PETSc Tutorials page, but you can look at this
>
>
> http://www.geodynamics.org/cig/community/workinggroups/short/workshops/cdm2013/presentations/SessionIV_Solvers.pdf
>
> for a step-by-step example of a saddle-point problem at the end.
>
>     Matt
>
>
>>
>>
>>>     Matt
>>>
>>>
>>>>     Matt
>>>>>
>>>>>
>>>>>> but I get the following errror:
>>>>>>
>>>>>> [0]PETSC ERROR: --------------------- Error Message
>>>>>> ------------------------------------
>>>>>> [0]PETSC ERROR: Null argument, when expecting valid pointer!
>>>>>> [0]PETSC ERROR: Null Object: Parameter # 2!
>>>>>> [0]PETSC ERROR:
>>>>>> ------------------------------------------------------------------------
>>>>>> [0]PETSC ERROR: Petsc Release Version 3.4.2, Jul, 02, 2013
>>>>>> [0]PETSC ERROR: See docs/changes/index.html for recent updates.
>>>>>> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
>>>>>> [0]PETSC ERROR: See docs/index.html for manual pages.
>>>>>> [0]PETSC ERROR:
>>>>>> ------------------------------------------------------------------------
>>>>>> [0]PETSC ERROR: src/AdLemMain on a arch-linux2-cxx-debug named
>>>>>> edwards by bkhanal Fri Jul 26 14:23:40 2013
>>>>>> [0]PETSC ERROR: Libraries linked from
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/arch-linux2-cxx-debug/lib
>>>>>> [0]PETSC ERROR: Configure run at Fri Jul 19 14:25:01 2013
>>>>>> [0]PETSC ERROR: Configure options --with-cc=gcc --with-fc=g77
>>>>>> --with-cxx=g++ --download-f-blas-lapack=1 --download-mpich=1
>>>>>> -with-clanguage=cxx --download-hypre=1
>>>>>> [0]PETSC ERROR:
>>>>>> ------------------------------------------------------------------------
>>>>>> [0]PETSC ERROR: MatPtAP() line 8166 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/mat/interface/matrix.c
>>>>>> [0]PETSC ERROR: PCSetUp_MG() line 628 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/impls/mg/mg.c
>>>>>> [0]PETSC ERROR: PCSetUp() line 890 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/interface/precon.c
>>>>>> [0]PETSC ERROR: KSPSetUp() line 278 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c
>>>>>> [0]PETSC ERROR: KSPSolve() line 399 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c
>>>>>> [0]PETSC ERROR: PCApply_FieldSplit_Schur() line 807 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/impls/fieldsplit/fieldsplit.c
>>>>>> [0]PETSC ERROR: PCApply() line 442 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/interface/precon.c
>>>>>> [0]PETSC ERROR: KSP_PCApply() line 227 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/include/petsc-private/kspimpl.h
>>>>>> [0]PETSC ERROR: KSPInitialResidual() line 64 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itres.c
>>>>>> [0]PETSC ERROR: KSPSolve_GMRES() line 239 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/impls/gmres/gmres.c
>>>>>> [0]PETSC ERROR: KSPSolve() line 441 in
>>>>>> /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c
>>>>>>
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> 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
>>>>>
>>>>
>>>>
>>>
>>>
>>> --
>>> 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
>>>
>>
>>
>
>
> --
> 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
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130726/d9241399/attachment-0001.html>