[petsc-users] Implementation of Multilevel preconditioner

Barry Smith bsmith at mcs.anl.gov
Thu Apr 5 10:38:37 CDT 2012


On Apr 5, 2012, at 10:22 AM, Abdul Hanan Sheikh wrote:

> Hi users and developers, 
> Summer greetings from Delft. 
> 
> 
> For the system Au = f, I want to apply a multilevel preconditioner with a KSP (say FGMRES) . The preconditioner
> reads as
> 
> Prec = C +  M^-1 (I - A*C) , where
> 
> 	•  C reads coarse grid correction operator i.e. C =  P*A_2h\R  [ R  restriction, A_2h coarse operator, P interpolation]
> 	• M is say some sparse matrix (resembling to A) 
> What makes it multilevel is I have to approximate coarse operator "A_2h" with few KSP iterations preconditioned by the above defined
> preconditioner "Prec" , but off course at the level 2h and hence at every coarse level until coarsest. 
> 
> I only know, it can be implemented with PCMG. I am little afraid of DMMG. 
> Explanatory suggestions to customize PCMG for what is required are very appreciated. 


    C +  M^-1 (I - A*C)   is exactly a multiplicative Schwarz method with two preconditioners. Or in another way to put it it is exactly a two level multigrid so you can use PCMG in a straightforward way. 

    Do NOT use DMMG it is deprecated and being removed from PETSc.

    To use PCMG you must provide the P and R with PCMGSetInterpolation and PCMGSetRestriction. You would us
 PCMGGetCoarseSolve(pc,&user.ksp_coarse);CHKERRQ(ierr);
   ierr = KSPSetOptionsPrefix(user.ksp_coarse,"coarse_");CHKERRQ(ierr);
  ierr = KSPSetFromOptions(user.ksp_coarse);CHKERRQ(ierr);
  ierr = KSPSetOperators(user.ksp_coarse,user.coarse.J,user.coarse.J,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);

   to provide the coarse A_c and depending what your "smoother" M is you would possibly pass a different matrix into the second location of KSPSetOperators() for each level to represent the "M is say some sparse matrix resembling A."

    src/ksp/ksp/examples/tests/ex19.c is a simple example showing how all the values are set for two levels.

   Barry

> 
> Thanks in advance. 
> 
> 
> regards, Abdul 
> 
> 



More information about the petsc-users mailing list