[petsc-users] Problem in PCShellSetApply()

[Barry Smith]

On Apr 18, 2011, at 9:42 PM, Rongliang Chen wrote:

> Hi,
> 
> I faced a problem when I used the function PCShellSetApply to set a composite PC in Petsc 3.1. There was no problem when I compiled my code but I got the following message when I run the code.
> ------------------------------------------------------------------------------------------------------------
> [0]PETSC ERROR: [1]PETSC ERROR: --------------------- Error Message ------------------------------------
> [1]PETSC ERROR: No support for this operation for this object type!
> [1]PETSC ERROR: Matrix format mpiaij does not have a built-in PETSc direct solver!
> [1]PETSC ERROR: ------------------------------------------------------------------------
> [1]PETSC ERROR: Petsc Release Version 3.1.0, Patch 8, Thu Mar 17 13:37:48 CDT 2011
> ...............................
> [7]PETSC ERROR: MatGetFactor() line 3644 in src/mat/interface/matrix.c
> [7]PETSC ERROR: PCSetUp_LU() line 133 in src/ksp/pc/impls/factor/lu/lu.c
> [7]PETSC ERROR: PCSetUp() line 795 in src/ksp/pc/interface/precon.c
> [7]PETSC ERROR: PCApply() line 353 in src/ksp/pc/interface/precon.c
> [7]PETSC ERROR: PCApply_Composite_Additive() line 102 in src/ksp/pc/impls/composite/composite.c
> [7]PETSC ERROR: PCApply() line 357 in src/ksp/pc/interface/precon.c
> [7]PETSC ERROR: PCApplyBAorAB() line 582 in src/ksp/pc/interface/precon.c
> [7]PETSC ERROR: GMREScycle() line 161 in src/ksp/ksp/impls/gmres/gmres.c
> [7]PETSC ERROR: KSPSolve_GMRES() line 241 in src/ksp/ksp/impls/gmres/gmres.c
> [7]PETSC ERROR: KSPSolve() line 396 in src/ksp/ksp/interface/itfunc.c
> [1]PETSC ERROR: SNES_KSPSolve() line 2944 in src/snes/interface/snes.c
> [1]PETSC ERROR: SNESSolve_LS() line 191 in src/snes/impls/ls/ls.c
> [1]PETSC ERROR: SNESSolve() line 2255 in src/snes/interface/snes.c
> ..............................
> --------------------------------------------------------------------------------------------------------------
> 
> I found that the function PCShellSetApply did not call the user defined function CoarseSolvePCApply (I called PCShellSetApply like this: ierr = PCShellSetApply(coarsesolve, CoarseSolvePCApply);CHKERRQ(ierr); ). Can anyone tell me what may be the problem? The related code attached.
> 
> -----------------------The related code-----------------------------------
> PetscErrorCode SetupPreconditioner(JoabCtx* ctx, SNES snes)
> {
>   JoabParameters *parameters = &ctx->parameters;
>   JoabView *view = &ctx->view;
>   JoabGrid *grid = &ctx->grid;
>   JoabGrid *coarsegrid = &ctx->coarsegrid;
>   JoabAlgebra *algebra = &ctx->algebra;
>   JoabAlgebra *coarsealgebra = &ctx->coarsealgebra;
> 
>   KSP fineksp;
>   PC finepc, coarsepc;
>   PC asmpc;
>   PC coarsesolve;
>   Vec fineones;
>   PetscInt veclength;
>   PetscScalar *values;
>   int i;
> 
>   PetscErrorCode ierr;
>   
>   PetscFunctionBegin;
>   ierr = SNESGetKSP(snes,&fineksp);CHKERRQ(ierr);
>   ierr = KSPGetPC(fineksp,&finepc);CHKERRQ(ierr);
>   ierr = KSPCreate(PETSC_COMM_WORLD,&ctx->coarseksp);CHKERRQ(ierr);
> 
>   ierr = SetMyKSPDefaults(ctx->coarseksp);CHKERRQ(ierr); 
>   ierr = KSPSetFromOptions(ctx->coarseksp);CHKERRQ(ierr);
>   ierr = KSPSetOperators(ctx->coarseksp, coarsealgebra->H, coarsealgebra->H,
>                SAME_NONZERO_PATTERN);CHKERRQ(ierr);
>   if (parameters->geometric_asm) {
>     ierr = KSPGetPC(ctx->coarseksp, &coarsepc);CHKERRQ(ierr);
>     ierr = PCASMSetOverlap(coarsepc,0);CHKERRQ(ierr);
>     ierr = PCASMSetLocalSubdomains(coarsepc,1,&coarsegrid->df_global_asm, PETSC_NULL);CHKERRQ(ierr);
>   }
> 
>   ierr = PCSetType(finepc,PCCOMPOSITE);CHKERRQ(ierr);
>   ierr = PCCompositeAddPC(finepc,PCSHELL);CHKERRQ(ierr);
>   ierr = PCCompositeAddPC(finepc,PCASM);CHKERRQ(ierr); 
> 
>   /* set up asm (fine) part of two-level preconditioner */
>   ierr = PCCompositeGetPC(finepc,1,&asmpc);CHKERRQ(ierr);
>   if (parameters->geometric_asm) {

   Is this flag set so it actually sets the type to PCASM

>     ierr = PCSetType(asmpc,PCASM);CHKERRQ(ierr);
>     ierr = PCASMSetOverlap(asmpc,0);CHKERRQ(ierr);
>     ierr = PCASMSetLocalSubdomains(asmpc,1,&grid->df_global_asm, PETSC_NULL);CHKERRQ(ierr);
>   }
>   ierr = SetMyPCDefaults(asmpc);CHKERRQ(ierr); 

   What is this setting the solver to?

   One of the two composite PC's is being set to LU and it is parallel so it fails. My guess is that SetMyPCDefaults() is setting the the PCType to LU. If not you'll need to track through your code, perhaps put a break point in PCSetType() and see where the type is being set to LU.

   Barry

>   ierr = PCSetFromOptions(asmpc);CHKERRQ(ierr);
> 
>   /* set up coarse solve part of two-level preconditioner */
>   ierr = PCCompositeGetPC(finepc,0,&coarsesolve);CHKERRQ(ierr);
>   ierr = PCShellSetContext(coarsesolve,ctx);CHKERRQ(ierr);
>   ierr = PCShellSetApply(coarsesolve, CoarseSolvePCApply);CHKERRQ(ierr);
>   
>   PetscFunctionReturn(0);
> }
> 
> PetscErrorCode CoarseSolvePCApply(PC pc, Vec xin, Vec xout)
> {
>   JoabCtx *ctx;
>   PetscErrorCode ierr;
>   JoabParameters *parameters;
>   JoabView *view;
>   JoabGrid *finegrid;
>   JoabGrid *coarsegrid;
>   JoabAlgebra *finealgebra;
>   JoabAlgebra *coarsealgebra;
>   PetscInt its;
>   PetscLogDouble t1,t2,v1,v2;
>   KSPConvergedReason reason;
> 
>   PetscFunctionBegin;
>   ierr = PetscPrintf(PETSC_COMM_WORLD,"Setup coarse level preconditioner.....\n");CHKERRQ(ierr);
>   ierr = PCShellGetContext(pc,(void**)&ctx);CHKERRQ(ierr);
>   parameters = &ctx->parameters;
>   view = &ctx->view;
>   finegrid = &ctx->grid;
>   coarsegrid = &ctx->coarsegrid;
>   finealgebra = &ctx->algebra;
>   coarsealgebra = &ctx->coarsealgebra;
> 
>   ierr = PetscGetTime(&t1);CHKERRQ(ierr);
> 
>   parameters->whichlevel = COARSE_GRID;
> 
>   /* restrict fine grid to coarse grid */
>   ierr = PetscGetTime(&v1);CHKERRQ(ierr);
>   ierr = VecSet(coarsealgebra->predictedShape_opt,0.0);CHKERRQ(ierr);CHKERRQ(ierr);
>   ierr = ApplyRestriction(ctx,xin,coarsealgebra->predictedShape_opt);CHKERRQ(ierr); 
>   
>   ierr = VecSet(coarsealgebra->solutionShape_opt,0.0);CHKERRQ(ierr);
>   ierr = KSPSetTolerances(ctx->coarseksp,1e-6,1e-14,PETSC_DEFAULT,1000);CHKERRQ(ierr);
>   ierr = KSPSolve(ctx->coarseksp, coarsealgebra->predictedShape_opt, coarsealgebra->solutionShape_opt);CHKERRQ(ierr);
> 
>   /* interpolate coarse grid to fine grid */
>   ierr = MatInterpolate(ctx->Interp,coarsealgebra->solutionShape_opt,xout);CHKERRQ(ierr); 
> 
>   PetscFunctionReturn(0);
> }
> 
> 
> Regards,
> Rongliang
> 
>