[petsc-users] Use block Jacobi preconditioner with SNES

[Smith, Barry F.]

   Actually you do not need my new branch to achieve what you desired. All you need in your main program is something like

  ierr = SNESCreate(PETSC_COMM_WORLD,&snes);CHKERRQ(ierr);
  ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr);
  ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
  ierr = PCSetType(pc,PCBJACOBI);CHKERRQ(ierr);
  ierr = PCBJacobiSetTotalBlocks(pc,3,lens);CHKERRQ(ierr);  /* here you set your block sizes to whatever you need */

Then simply do not call PCBJacobiGetSubKSP() but use the options database to set the inner solver with -sub_pc_type lu -sub_ksp_type preonly 

    I have updated the branch to move the PCBJacobiSetTotalBlocks() to the main program but left the callback in there for setting the inner solver types (though as I just said you don't need to use the callback since you can control the solver from the options database). The callback is needed, if, for example, you wished to use a different solver on different blocks (which is not your case).

   Barry

  PETSc developers - do you think we should put the callback functionality into PETSc? It allows doing things that are otherwise not doable but is rather ugly (perhaps too specialized)?



> On Aug 27, 2018, at 4:01 PM, Smith, Barry F. <bsmith at mcs.anl.gov> wrote:
> 
> 
>    I have added new functionality within SNES to allow one to do as you desire in the branch barry/feature-snessetkspsetupcallback. Please see src/snes/examples/tests/ex3.c 
> 
>    Please let us know if it works for you and if you have any questions or problems.
> 
>    Barry
> 
>   Note that it has to be handled by a callback called from within the SNES solver because that is the first time the matrix exists in a form that the block sizes may be set.
> 
>    Sorry for the runaround with so many emails. 
> 
> 
>> On Aug 24, 2018, at 3:48 PM, Ali Reza Khaz'ali <arkhazali at cc.iut.ac.ir> wrote:
>> 
>> Hi,
>> 
>> I am trying to use block Jacobi preconditioner in SNES (SNESNEWTONLS). However, PCBJacobiGetSubKSP function returns an error stating "Object is in wrong state, Must call KSPSetUp() or PCSetUp() first". When I add KSPSetUp, I got and error from them as: "Must call DMShellSetGlobalVector() or DMShellSetCreateGlobalVector()", and if PCSetUp is added, "Object is in wrong state, Matrix must be set first" error is printed.
>> 
>> Below is a part of my code. It is run serially. Any help is much appreciated.
>> 
>>    ierr = SNESGetKSP(snes, &Petsc_ksp);
>>    CHKERRQ(ierr);
>>    ierr = KSPGetPC(Petsc_ksp, &Petsc_pc);
>>    CHKERRQ(ierr);
>>    ierr = KSPSetTolerances(Petsc_ksp, 1.e-3, 1e-3, PETSC_DEFAULT, 20000);
>>    CHKERRQ(ierr);
>>    ierr = SNESSetTolerances(snes, 1e-1, 1e-1, 1e-1, 2000, 2000);
>>    CHKERRQ(ierr);
>>    ierr = SNESSetType(snes, SNESNEWTONLS);
>>    CHKERRQ(ierr);
>>    ierr = KSPSetType(Petsc_ksp, KSPGMRES);
>>    CHKERRQ(ierr);
>>    ierr = PCSetType(Petsc_pc, PCBJACOBI);
>>    CHKERRQ(ierr);
>>    ierr = PCSetType(Petsc_pc, PCBJACOBI);
>>    CHKERRQ(ierr);
>>    ierr = PCBJacobiSetTotalBlocks(Petsc_pc, 2*Nx*Ny*Nz, SadeqSize);
>>    CHKERRQ(ierr);
>> 
>>    SNESSetUp(snes);
>>    CHKERRQ(ierr);
>>    ierr = PCBJacobiGetSubKSP(Petsc_pc, &nLocal, &Firstly, &subKSP);
>>    CHKERRQ(ierr);
>> 
>>    for (i = 0; i < nLocal; i++) {
>>        ierr = KSPGetPC(subKSP[i], &SubPc);
>>        CHKERRQ(ierr);
>>        ierr = PCSetType(SubPc, PCLU);
>>        CHKERRQ(ierr);
>>        ierr = PCFactorSetMatSolverPackage(SubPc, "mkl_pardiso");
>>        CHKERRQ(ierr);
>>        ierr = KSPSetType(subKSP[i], KSPPREONLY);
>>        CHKERRQ(ierr);
>>        ierr = KSPSetTolerances(subKSP[i], 1.e-6, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT);
>>        CHKERRQ(ierr);
>>    }
>>    ierr = SNESSolve(snes, NULL, Petsc_X);
>>    CHKERRQ(ierr);
>> 
>> 
>> 
>> -- 
>> Ali Reza Khaz’ali
>> Assistant Professor of Petroleum Engineering,
>> Department of Chemical Engineering
>> Isfahan University of Technology
>> Isfahan, Iran
>> 
>