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 PetscCall(PetscInitialize(&argc,&args,NULL,help)); 19 PetscCall(PetscStrcpy(dir,".")); 20 PetscCall(PetscOptionsGetString(NULL,NULL,"-load_dir",dir,sizeof(dir),NULL)); 21 PetscCall(PetscOptionsGetInt(NULL,NULL,"-nmat",&nmat,NULL)); 22 PetscCall(PetscOptionsGetBool(NULL,NULL,"-reset",&reset,NULL)); 23 PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 24 PetscCall(KSPCreate(PETSC_COMM_WORLD,&ksp)); 25 PetscCall(KSPSetOperators(ksp,A,A)); 26 for (i=0; i<nmat; i++) { 27 j = i+400; 28 PetscCall(PetscSNPrintf(name,sizeof(name),"%s/A_%d.dat",dir,j)); 29 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer)); 30 PetscCall(MatLoad(A,viewer)); 31 PetscCall(PetscViewerDestroy(&viewer)); 32 if (i == 0) { 33 PetscCall(MatCreateVecs(A,&x,&b)); 34 } 35 PetscCall(PetscSNPrintf(name,sizeof(name),"%s/rhs_%d.dat",dir,j)); 36 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer)); 37 PetscCall(VecLoad(b,viewer)); 38 PetscCall(PetscViewerDestroy(&viewer)); 39 PetscCall(KSPSetFromOptions(ksp)); 40 PetscCall(KSPSolve(ksp,b,x)); 41 PetscCall(PetscObjectTypeCompare((PetscObject)ksp,KSPHPDDM,&flg)); 42 #if defined(PETSC_HAVE_HPDDM) 43 if (flg && reset) { 44 PetscCall(KSPHPDDMGetDeflationSpace(ksp,&U)); 45 PetscCall(KSPReset(ksp)); 46 PetscCall(KSPSetOperators(ksp,A,A)); 47 PetscCall(KSPSetFromOptions(ksp)); 48 PetscCall(KSPSetUp(ksp)); 49 if (U) { 50 PetscCall(KSPHPDDMSetDeflationSpace(ksp,U)); 51 PetscCall(MatDestroy(&U)); 52 } 53 } 54 #endif 55 } 56 PetscCall(VecDestroy(&x)); 57 PetscCall(VecDestroy(&b)); 58 PetscCall(MatDestroy(&A)); 59 PetscCall(KSPDestroy(&ksp)); 60 PetscCall(PetscFinalize()); 61 return 0; 62 } 63 64 /*TEST 65 66 test: 67 suffix: 1 68 nsize: 1 69 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 70 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 72 test: 73 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 74 suffix: 1_icc 75 nsize: 1 76 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 78 testset: 79 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 80 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 81 test: 82 nsize: 1 83 suffix: 2_seq 84 output_file: output/ex75_2.out 85 test: 86 nsize: 2 87 suffix: 2_par 88 output_file: output/ex75_2.out 89 90 testset: 91 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 92 nsize: 1 93 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 94 test: 95 suffix: 2_icc 96 args: 97 test: 98 suffix: 2_icc_atol 99 args: -ksp_atol 1e-12 100 101 test: 102 requires: hpddm datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES) 103 nsize: 2 104 suffix: symmetric 105 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 107 TEST*/ 108