[petsc-users] V-cycle multigrid with matrix shells

[Barry Smith]

   Barry


On May 10, 2011, at 9:20 PM, Sylvain Barbot wrote:

> Dear Jed,
> 
> Thank you for your previous comments. I am still trying to design a
> multigrid preconditionner for the Navier's equation of elasticity. I
> have defined shell matrices for the interpolation, restriction,
> smoothing, residual calculation and the fine-level solver (see
> attached my_green.f90 file with stub functions). I thought that all
> these shell matrices would be enough to do the job, but upon running
> KSPSetup(), I get the following error:
> 
> [0]PETSC ERROR: No support for this operation for this object type!
> [0]PETSC ERROR: This matrix type, shell, does not support getting
> diagonal block!
> [0]PETSC ERROR:
> ------------------------------------------------------------------------
> [0]PETSC ERROR: Petsc Release Version 3.1.0, Patch 4, Fri Jul 30
> 14:42:02 CDT 2010
> [0]PETSC ERROR:
> ------------------------------------------------------------------------
> [0]PETSC ERROR: PCSetUp_BJacobi() line 162 in src/ksp/pc/impls/bjacobi/bjacobi.c
> [0]PETSC ERROR: PCSetUp() line 795 in src/ksp/pc/interface/precon.c
> [0]PETSC ERROR: KSPSetUp() line 237 in src/ksp/ksp/interface/itfunc.c
> [0]PETSC ERROR: PCSetUp_MG() line 574 in src/ksp/pc/impls/mg/mg.c
> 
> regarding your previous comments: "For Jacobi, you can have your
> MatShell implement MatGetDiagonal and it will work. (...) If you want
> to use a multiplicative relaxation like SOR, you would have to
> implement it yourself.", I actually want to use SOR with my smoother
> shell matrix. What am I doing wrong to make petsc still believe I want
> to use Jacobi?
> 

 By default, it is using block Jacobi on each level (except the coarsest). You need to tell it to use the relaxation you have coded. For example if your shell matrix provides a SOR/relaxation method then you need to use -mg_levels_pc_type sor  and it will call that on each level instead of block Jacobi. In looking at your code I see the only operation your provide your shell matrices are matrix-vector multiply. Therefore you cannot use SOR/relaxation. The only think you can use is GMRES (or some other Krylov method) with no preconditioner as the smoother. For this you can use -mg_levels_pc_type none -mg_levels_ksp_type gmres 

   In your message  you say  "I actually want to use SOR with my smoother" to do this you need to provide the smoother function and set it with MatShellSetOperations(). 


> another related question: with a direct solver, one needs to assign a
> shell matrix with a KSP object (KSPSetOperators); is it still required
> when using multigrid? if I don't do so, I get the following error:
> 
  Where are you using the direct solver for, on the entire linear system or on the coarse mesh? In either case you NEED to provide an explicit (not shell) sparse matrix otherwise a direct solver cannot be used. From the message below it looks like you did not provide a matrix somewhere. If you don't provide a matrix a direct solver cannot be used.


   Barry

> [0]PETSC ERROR: Null argument, when expecting valid pointer!
> [0]PETSC ERROR: Null Object: Parameter # 1!
> [0]PETSC ERROR: MatGetVecs() line 7265 in src/mat/interface/matrix.c
> [0]PETSC ERROR: KSPGetVecs() line 806 in src/ksp/ksp/interface/iterativ.c
> [0]PETSC ERROR: KSPSetUp_GMRES() line 94 in src/ksp/ksp/impls/gmres/gmres.c
> [0]PETSC ERROR: KSPSetUp() line 199 in src/ksp/ksp/interface/itfunc.c
> 
> however, I don't understand why I would need this matrix in a
> multigrid framework.
> 
> Is there any example out there of using PCMG with shell matrices?
> 
> Thanks in advance for your response.
> Best wishes,
> Sylvain Barbot
> 
> 
> 2010/10/7 Jed Brown <jed at 59a2.org>:
>> On Thu, Oct 7, 2010 at 11:45, Sylvain Barbot <sylbar.vainbot at gmail.com>
>> wrote:
>>> 
>>> questions:
>>> 1) is it at all possible to specify this mode of operation, from the
>>> finest to the coarser level, and back? Any examples out there?
>> 
>> -pc_mg_type multiplicative, or
>>   PCMGSetType(pc,PC_MG_MULTIPLICATIVE);
>> 
>>> 
>>> 2) is it readily possible to use matrix shells with DMMG? I imagine
>>> the Jacobian matrix may simply be provided as a matrix shell. Is there
>>> any examples of multi-grid methods with shell matrices online?
>> 
>> You can do this, but you will have to define a custom smoother.  For Jacobi,
>> you can have your MatShell implement MatGetDiagonal and it will work.  I
>> thought you could implement MatGetDiagonalBlock for PBJacobi, but it's not
>> currently written that way (though it should be and that would be an easy
>> change to make).  If you want to use a multiplicative relaxation like SOR,
>> you would have to implement it yourself.  If you need something like ILU for
>> a smoother, then you will have to pay for a matrix.  Note that one
>> possibility is to assemble a matrix on all but the finest level so you can
>> use stronger smoothers there, then make do with Jacobi on the finest level.
>> 
>>> 
>>> 3) to deal with non-uniform sampling: can I provide the coordinates of
>>> the finest grid with DASetCoordinates, then expect DMMG to provide the
>>> subsampled coordinates at the coarser levels?
>> 
>> Currently no, you have to set them on each level.  Perhaps you could do this
>> by rolling a loop over levels and applying MatRestrict using the restriction
>> matrix from DMMG (this might not be the sampling that you want).
>> Jed
> <my_green.f90>