Actual source code: ex79.c

petsc-master 2020-08-25
Report Typos and Errors
  1:  #include <petsc.h>

  3: static char help[] = "Solves a linear system with a block of right-hand sides, apply a preconditioner to the same block.\n\n";

  5: PetscErrorCode MatApply(PC pc, Mat X, Mat Y)
  6: {

 10:   MatCopy(X,Y,SAME_NONZERO_PATTERN);
 11:   return(0);
 12: }

 14: int main(int argc,char **args)
 15: {
 16:   Mat                X,B;         /* computed solutions and RHS */
 17:   Mat                A;           /* linear system matrix */
 18:   KSP                ksp;         /* linear solver context */
 19:   PC                 pc;          /* preconditioner context */
 20:   PetscInt           m = 10;
 21:   PetscLogEvent      event;
 22:   PetscEventPerfInfo info;
 23:   PetscErrorCode     ierr;

 25:   PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
 26:   PetscLogDefaultBegin();
 27:   PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
 28:   MatCreateAIJ(PETSC_COMM_WORLD,m,m,PETSC_DECIDE,PETSC_DECIDE,m,NULL,m,NULL,&A);
 29:   MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,m,NULL,&B);
 30:   MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,m,NULL,&X);
 31:   MatSetRandom(A,NULL);
 32:   MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
 33:   MatShift(A,10.0);
 34:   MatSetRandom(B,NULL);
 35:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 36:   KSPSetOperators(ksp,A,A);
 37:   KSPSetFromOptions(ksp);
 38:   KSPGetPC(ksp,&pc);
 39:   PCShellSetMatApply(pc,MatApply);
 40:   KSPMatSolve(ksp,B,X);
 41:   PCMatApply(pc,B,X);
 42:   MatDestroy(&X);
 43:   MatDestroy(&B);
 44:   MatDestroy(&A);
 45:   KSPDestroy(&ksp);
 46:   PetscLogEventRegister("PCApply",PC_CLASSID,&event);
 47:   PetscLogEventGetPerfInfo(PETSC_DETERMINE,event,&info);
 48:   if (PetscDefined(USE_LOG) && m > 1 && info.count) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_PLIB,"PCApply() called %d times",info.count);
 49:   PetscLogEventRegister("PCMatApply",PC_CLASSID,&event);
 50:   PetscLogEventGetPerfInfo(PETSC_DETERMINE,event,&info);
 51:   if (PetscDefined(USE_LOG) && m > 1 && !info.count) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_PLIB,"PCMatApply() never called");
 52:   PetscFinalize();
 53:   return ierr;
 54: }

 56: /*TEST
 57:    # KSPHPDDM does either pseudo-blocking or "true" blocking, all tests should succeed with other -ksp_hpddm_type
 58:    testset:
 59:       nsize: 1
 60:       args: -pc_type {{bjacobi lu ilu mat cholesky icc none shell asm gasm}shared output}
 61:       test:
 62:          suffix: 1
 63:          output_file: output/ex77_preonly.out
 64:          args: -ksp_type preonly
 65:       test:
 66:          suffix: 1_hpddm
 67:          output_file: output/ex77_preonly.out
 68:          requires: hpddm
 69:          args: -ksp_type hpddm -ksp_hpddm_type preonly

 71:    testset:
 72:       nsize: 1
 73:       args: -pc_type ksp
 74:       test:
 75:          suffix: 2
 76:          output_file: output/ex77_preonly.out
 77:          args: -ksp_ksp_type preonly -ksp_type preonly
 78:       test:
 79:          suffix: 2_hpddm
 80:          output_file: output/ex77_preonly.out
 81:          requires: hpddm
 82:          args: -ksp_ksp_type hpddm -ksp_type hpddm -ksp_hpddm_type preonly -ksp_ksp_hpddm_type preonly

 84:    testset:
 85:       nsize: 1
 86:       requires: hara
 87:       args: -pc_type hara
 88:       test:
 89:          suffix: 3
 90:          output_file: output/ex77_preonly.out
 91:          args: -ksp_type preonly
 92:       test:
 93:          suffix: 3_hpddm
 94:          output_file: output/ex77_preonly.out
 95:          requires: hpddm
 96:          args: -ksp_type hpddm -ksp_hpddm_type preonly

 98:    testset:
 99:       nsize: 1
100:       requires: spai
101:       args: -pc_type spai
102:       test:
103:          suffix: 4
104:          output_file: output/ex77_preonly.out
105:          args: -ksp_type preonly
106:       test:
107:          suffix: 4_hpddm
108:          output_file: output/ex77_preonly.out
109:          requires: hpddm
110:          args: -ksp_type hpddm -ksp_hpddm_type preonly
111:    # special code path in PCApplyMat() for PCBJACOBI when a block is shared by multiple processes
112:    testset:
113:       nsize: 2
114:       args: -pc_type bjacobi -pc_bjacobi_blocks 1 -sub_pc_type none
115:       test:
116:          suffix: 5
117:          output_file: output/ex77_preonly.out
118:          args: -ksp_type preonly -sub_ksp_type preonly
119:       test:
120:          suffix: 5_hpddm
121:          output_file: output/ex77_preonly.out
122:          requires: hpddm
123:          args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm
124:    # special code path in PCApplyMat() for PCGASM when a block is shared by multiple processes
125:    testset:
126:       nsize: 2
127:       args: -pc_type gasm -pc_gasm_total_subdomains 1 -sub_pc_type none
128:       test:
129:          suffix: 6
130:          output_file: output/ex77_preonly.out
131:          args: -ksp_type preonly -sub_ksp_type preonly
132:       test:
133:          suffix: 6_hpddm
134:          output_file: output/ex77_preonly.out
135:          requires: hpddm
136:          args: -ksp_type hpddm -ksp_hpddm_type preonly -sub_ksp_type hpddm


139: TEST*/