1 #include <petsc.h> 2 3 static char help[] = "Solves a series of linear systems using KSPHPDDM.\n\n"; 4 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; 17 18 PetscFunctionBeginUser; 19 PetscCall(PetscInitialize(&argc,&args,NULL,help)); 20 PetscCall(PetscStrcpy(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) { 34 PetscCall(MatCreateVecs(A,&x,&b)); 35 } 36 PetscCall(PetscSNPrintf(name,sizeof(name),"%s/rhs_%" PetscInt_FMT ".dat",dir,j)); 37 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer)); 38 PetscCall(VecLoad(b,viewer)); 39 PetscCall(PetscViewerDestroy(&viewer)); 40 PetscCall(KSPSetFromOptions(ksp)); 41 PetscCall(KSPSolve(ksp,b,x)); 42 PetscCall(PetscObjectTypeCompare((PetscObject)ksp,KSPHPDDM,&flg)); 43 #if defined(PETSC_HAVE_HPDDM) 44 if (flg && reset) { 45 PetscCall(KSPHPDDMGetDeflationMat(ksp,&U)); 46 PetscCall(KSPReset(ksp)); 47 PetscCall(KSPSetOperators(ksp,A,A)); 48 PetscCall(KSPSetFromOptions(ksp)); 49 PetscCall(KSPSetUp(ksp)); 50 if (U) { 51 PetscCall(KSPHPDDMSetDeflationMat(ksp,U)); 52 PetscCall(MatDestroy(&U)); 53 } 54 } 55 #endif 56 } 57 PetscCall(VecDestroy(&x)); 58 PetscCall(VecDestroy(&b)); 59 PetscCall(MatDestroy(&A)); 60 PetscCall(KSPDestroy(&ksp)); 61 PetscCall(PetscFinalize()); 62 return 0; 63 } 64 65 /*TEST 66 67 test: 68 suffix: 1 69 nsize: 1 70 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 71 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 72 73 test: 74 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 75 suffix: 1_icc 76 nsize: 1 77 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 78 79 testset: 80 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 81 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 82 test: 83 nsize: 1 84 suffix: 2_seq 85 output_file: output/ex75_2.out 86 test: 87 nsize: 2 88 suffix: 2_par 89 output_file: output/ex75_2.out 90 91 testset: 92 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 93 nsize: 1 94 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 95 test: 96 suffix: 2_icc 97 args: 98 test: 99 suffix: 2_icc_atol 100 args: -ksp_atol 1e-12 101 102 test: 103 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES) 104 nsize: 2 105 suffix: symmetric 106 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 107 108 TEST*/ 109