1 #include <petsc.h>
2
3 static char help[] = "Solves a series of linear systems using KSPHPDDM.\n\n";
4
main(int argc,char ** args)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 PetscErrorCode ierr;
18
19 ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
20 ierr = PetscStrcpy(dir,".");CHKERRQ(ierr);
21 ierr = PetscOptionsGetString(NULL,NULL,"-load_dir",dir,sizeof(dir),NULL);CHKERRQ(ierr);
22 ierr = PetscOptionsGetInt(NULL,NULL,"-nmat",&nmat,NULL);CHKERRQ(ierr);
23 ierr = PetscOptionsGetBool(NULL,NULL,"-reset",&reset,NULL);CHKERRQ(ierr);
24 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
25 ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
26 ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr);
27 for (i=0; i<nmat; i++) {
28 j = i+400;
29 ierr = PetscSNPrintf(name,sizeof(name),"%s/A_%d.dat",dir,j);CHKERRQ(ierr);
30 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer);CHKERRQ(ierr);
31 ierr = MatLoad(A,viewer);CHKERRQ(ierr);
32 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
33 if (i == 0) {
34 ierr = MatCreateVecs(A,&x,&b);CHKERRQ(ierr);
35 }
36 ierr = PetscSNPrintf(name,sizeof(name),"%s/rhs_%d.dat",dir,j);CHKERRQ(ierr);
37 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer);CHKERRQ(ierr);
38 ierr = VecLoad(b,viewer);CHKERRQ(ierr);
39 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);
40 ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
41 ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr);
42 ierr = PetscObjectTypeCompare((PetscObject)ksp,KSPHPDDM,&flg);CHKERRQ(ierr);
43 #if defined(PETSC_HAVE_HPDDM)
44 if (flg && reset) {
45 ierr = KSPHPDDMGetDeflationSpace(ksp,&U);CHKERRQ(ierr);
46 ierr = KSPReset(ksp);CHKERRQ(ierr);
47 ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr);
48 ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
49 ierr = KSPSetUp(ksp);CHKERRQ(ierr);
50 if (U) {
51 ierr = KSPHPDDMSetDeflationSpace(ksp,U);CHKERRQ(ierr);
52 ierr = MatDestroy(&U);CHKERRQ(ierr);
53 }
54 }
55 #endif
56 }
57 ierr = VecDestroy(&x);CHKERRQ(ierr);
58 ierr = VecDestroy(&b);CHKERRQ(ierr);
59 ierr = MatDestroy(&A);CHKERRQ(ierr);
60 ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
61 ierr = PetscFinalize();
62 return ierr;
63 }
64
65 /*TEST
66
67 test:
68 suffix: 1
69 nsize: 1
70 requires: hpddm datafilespath double !complex !define(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 !define(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 !define(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 !define(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 !define(PETSC_USE_64BIT_INDICES) slepc
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