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