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 PetscBool flg; 13 #endif 14 PetscInt i,j,nmat = 10; 15 PetscViewer viewer; 16 char dir[PETSC_MAX_PATH_LEN],name[256]; 17 PetscBool reset = PETSC_FALSE; 18 PetscErrorCode ierr; 19 20 ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr; 21 ierr = PetscStrcpy(dir,".");CHKERRQ(ierr); 22 ierr = PetscOptionsGetString(NULL,NULL,"-load_dir",dir,sizeof(dir),NULL);CHKERRQ(ierr); 23 ierr = PetscOptionsGetInt(NULL,NULL,"-nmat",&nmat,NULL);CHKERRQ(ierr); 24 ierr = PetscOptionsGetBool(NULL,NULL,"-reset",&reset,NULL);CHKERRQ(ierr); 25 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 26 ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr); 27 ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr); 28 for (i=0; i<nmat; i++) { 29 j = i+400; 30 ierr = PetscSNPrintf(name,sizeof(name),"%s/A_%d.dat",dir,j);CHKERRQ(ierr); 31 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer);CHKERRQ(ierr); 32 ierr = MatLoad(A,viewer);CHKERRQ(ierr); 33 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 34 if (i == 0) { 35 ierr = MatCreateVecs(A,&x,&b);CHKERRQ(ierr); 36 } 37 ierr = PetscSNPrintf(name,sizeof(name),"%s/rhs_%d.dat",dir,j);CHKERRQ(ierr); 38 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer);CHKERRQ(ierr); 39 ierr = VecLoad(b,viewer);CHKERRQ(ierr); 40 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 41 ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr); 42 ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr); 43 #if defined(PETSC_HAVE_HPDDM) 44 ierr = PetscObjectTypeCompare((PetscObject)ksp,KSPHPDDM,&flg);CHKERRQ(ierr); 45 if (flg && reset) { 46 ierr = KSPHPDDMGetDeflationSpace(ksp,&U);CHKERRQ(ierr); 47 ierr = KSPReset(ksp);CHKERRQ(ierr); 48 ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr); 49 ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr); 50 ierr = KSPSetUp(ksp);CHKERRQ(ierr); 51 if (U) { 52 ierr = KSPHPDDMSetDeflationSpace(ksp,U);CHKERRQ(ierr); 53 ierr = MatDestroy(&U);CHKERRQ(ierr); 54 } 55 } 56 #endif 57 } 58 ierr = VecDestroy(&x);CHKERRQ(ierr); 59 ierr = VecDestroy(&b);CHKERRQ(ierr); 60 ierr = MatDestroy(&A);CHKERRQ(ierr); 61 ierr = KSPDestroy(&ksp);CHKERRQ(ierr); 62 ierr = PetscFinalize(); 63 return ierr; 64 } 65 66 /*TEST 67 68 test: 69 suffix: 1 70 nsize: 1 71 requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 72 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 73 74 test: 75 requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 76 suffix: 1_icc 77 nsize: 1 78 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 79 80 testset: 81 requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 82 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 83 test: 84 nsize: 1 85 suffix: 2_seq 86 output_file: output/ex75_2.out 87 test: 88 nsize: 2 89 suffix: 2_par 90 output_file: output/ex75_2.out 91 92 test: 93 requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) 94 suffix: 2_icc 95 nsize: 1 96 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 -ksp_hpddm_recycle 20 -reset {{false true}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR 97 98 TEST*/ 99