Actual source code: ex75.c

  1: #include <petsc.h>

  3: static char help[] = "Solves a series of linear systems using KSPHPDDM.\n\n";

  5: int main(int argc, char **args)
  6: {
  7:   Vec x, b; /* computed solution and RHS */
  8:   Mat A;    /* linear system matrix */
  9:   KSP ksp;  /* linear solver context */
 10: #if defined(PETSC_HAVE_HPDDM)
 11:   Mat U; /* deflation space */
 12: #endif
 13:   PetscInt    i, j, nmat = 10;
 14:   PetscViewer viewer;
 15:   char        dir[PETSC_MAX_PATH_LEN], name[256];
 16:   PetscBool   flg, reset = PETSC_FALSE;

 18:   PetscFunctionBeginUser;
 19:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 20:   PetscCall(PetscStrncpy(dir, ".", sizeof(dir)));
 21:   PetscCall(PetscOptionsGetString(NULL, NULL, "-load_dir", dir, sizeof(dir), NULL));
 22:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nmat", &nmat, NULL));
 23:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-reset", &reset, NULL));
 24:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 25:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
 26:   PetscCall(KSPSetOperators(ksp, A, A));
 27:   for (i = 0; i < nmat; i++) {
 28:     j = i + 400;
 29:     PetscCall(PetscSNPrintf(name, sizeof(name), "%s/A_%" PetscInt_FMT ".dat", dir, j));
 30:     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, &viewer));
 31:     PetscCall(MatLoad(A, viewer));
 32:     PetscCall(PetscViewerDestroy(&viewer));
 33:     if (i == 0) PetscCall(MatCreateVecs(A, &x, &b));
 34:     PetscCall(PetscSNPrintf(name, sizeof(name), "%s/rhs_%" PetscInt_FMT ".dat", dir, j));
 35:     PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, name, FILE_MODE_READ, &viewer));
 36:     PetscCall(VecLoad(b, viewer));
 37:     PetscCall(PetscViewerDestroy(&viewer));
 38:     PetscCall(KSPSetFromOptions(ksp));
 39:     PetscCall(KSPSolve(ksp, b, x));
 40:     PetscCall(PetscObjectTypeCompare((PetscObject)ksp, KSPHPDDM, &flg));
 41: #if defined(PETSC_HAVE_HPDDM)
 42:     if (flg && reset) {
 43:       PetscCall(KSPHPDDMGetDeflationMat(ksp, &U));
 44:       PetscCall(KSPReset(ksp));
 45:       PetscCall(KSPSetOperators(ksp, A, A));
 46:       PetscCall(KSPSetFromOptions(ksp));
 47:       PetscCall(KSPSetUp(ksp));
 48:       if (U) {
 49:         PetscCall(KSPHPDDMSetDeflationMat(ksp, U));
 50:         PetscCall(MatDestroy(&U));
 51:       }
 52:     }
 53: #endif
 54:   }
 55:   PetscCall(VecDestroy(&x));
 56:   PetscCall(VecDestroy(&b));
 57:   PetscCall(MatDestroy(&A));
 58:   PetscCall(KSPDestroy(&ksp));
 59:   PetscCall(PetscFinalize());
 60:   return 0;
 61: }

 63: /*TEST

 65:    test:
 66:       suffix: 1
 67:       nsize: 1
 68:       requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
 69:       args: -nmat 1 -pc_type none -ksp_converged_reason -ksp_type {{gmres hpddm}shared output} -ksp_max_it 1000 -ksp_gmres_restart 1000 -ksp_rtol 1e-10 -ksp_hpddm_type {{gmres bgmres}shared output} -options_left no -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR

 71:    test:
 72:       requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
 73:       suffix: 1_icc
 74:       nsize: 1
 75:       args: -nmat 1 -pc_type icc -ksp_converged_reason -ksp_type {{gmres hpddm}shared output} -ksp_max_it 1000 -ksp_gmres_restart 1000 -ksp_rtol 1e-10 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR

 77:    testset:
 78:       requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
 79:       args: -nmat 3 -pc_type none -ksp_converged_reason -ksp_type hpddm -ksp_max_it 1000 -ksp_gmres_restart 40 -ksp_rtol 1e-10 -ksp_hpddm_type {{gcrodr bgcrodr}shared output} -ksp_hpddm_recycle 20 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
 80:       test:
 81:         nsize: 1
 82:         suffix: 2_seq
 83:         output_file: output/ex75_2.out
 84:       test:
 85:         nsize: 2
 86:         suffix: 2_par
 87:         output_file: output/ex75_2.out

 89:    testset:
 90:       requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
 91:       nsize: 1
 92:       args: -nmat 3 -pc_type icc -ksp_converged_reason -ksp_type hpddm -ksp_max_it 1000 -ksp_gmres_restart 40 -ksp_rtol 1e-10 -ksp_hpddm_type {{gcrodr bgcrodr}shared output} -ksp_hpddm_recycle 20 -reset {{false true}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
 93:       test:
 94:         suffix: 2_icc
 95:         args:
 96:       test:
 97:         suffix: 2_icc_atol
 98:         args: -ksp_atol 1e-12

100:    test:
101:       requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
102:       nsize: 2
103:       suffix: symmetric
104:       args: -nmat 3 -pc_type jacobi -ksp_converged_reason -ksp_type hpddm -ksp_max_it 1000 -ksp_gmres_restart 40 -ksp_atol 1e-11 -ksp_hpddm_type bgcrodr -ksp_hpddm_recycle 20 -reset {{false true}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR -ksp_hpddm_recycle_symmetric true

106: TEST*/