Actual source code: ex79.c

  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: {
  8:   MatCopy(X,Y,SAME_NONZERO_PATTERN);
  9:   return 0;
 10: }

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

 24:   PetscInitialize(&argc,&args,NULL,help);
 25:   PetscLogDefaultBegin();
 26:   PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
 27:   MatCreateAIJ(PETSC_COMM_WORLD,m,m,PETSC_DECIDE,PETSC_DECIDE,m,NULL,m,NULL,&A);
 28:   MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,m,NULL,&B);
 29:   MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,m,NULL,&X);
 30:   MatSetRandom(A,NULL);
 31:   MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
 32:   MatShift(A,10.0);
 33:   MatSetRandom(B,NULL);
 34:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 35:   KSPSetOperators(ksp,A,A);
 36:   KSPSetFromOptions(ksp);
 37:   KSPGetPC(ksp,&pc);
 38:   PCShellSetMatApply(pc,MatApply);
 39:   KSPMatSolve(ksp,B,X);
 40:   PCMatApply(pc,B,X);
 41:   MatDestroy(&X);
 42:   MatDestroy(&B);
 43:   MatDestroy(&A);
 44:   KSPDestroy(&ksp);
 45:   PetscLogEventRegister("PCApply",PC_CLASSID,&event);
 46:   PetscLogEventGetPerfInfo(PETSC_DETERMINE,event,&info);
 48:   PetscLogEventRegister("PCMatApply",PC_CLASSID,&event);
 49:   PetscLogEventGetPerfInfo(PETSC_DETERMINE,event,&info);
 51:   PetscFinalize();
 52:   return 0;
 53: }

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

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

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

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

137:    testset:
138:       nsize: 1
139:       requires: suitesparse
140:       args: -pc_type qr
141:       test:
142:          suffix: 7
143:          output_file: output/ex77_preonly.out
144:          args: -ksp_type preonly
145:       test:
146:          suffix: 7_hpddm
147:          output_file: output/ex77_preonly.out
148:          requires: hpddm
149:          args: -ksp_type hpddm -ksp_hpddm_type preonly

151: TEST*/