Actual source code: ex75f.F90
petsc-3.13.1 2020-05-02
1: !
2: ! Description: Solves a series of linear systems using KSPHPDDM.
3: !
5: program main
6: #include <petsc/finclude/petscksp.h>
7: use petscksp
8: implicit none
9: Vec x,b
10: Mat A
11: #if defined(PETSC_HAVE_HPDDM)
12: Mat U
13: #endif
14: KSP ksp
15: PetscInt i,j,nmat
16: PetscViewer viewer
17: character*(128) name
18: character*(128) dir
19: character*(8) fmt
20: character(3) cmat
21: PetscBool flg,reset
22: PetscErrorCode ierr
24: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
25: if (ierr .ne. 0) then
26: print *,'Unable to initialize PETSc'
27: stop
28: endif
29: dir = '.'
30: call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-load_dir',dir,flg,ierr);CHKERRA(ierr)
31: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-nmat',nmat,flg,ierr);CHKERRA(ierr)
32: reset = PETSC_FALSE
33: call PetscOptionsGetBool(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-reset',reset,flg,ierr);CHKERRA(ierr)
34: call MatCreate(PETSC_COMM_WORLD,A,ierr);CHKERRA(ierr)
35: call KSPCreate(PETSC_COMM_WORLD,ksp,ierr);CHKERRA(ierr)
36: call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
37: do 50 i=0,nmat-1
38: j = i+400
39: fmt = '(I3)'
40: write (cmat,fmt) j
41: write (name,'(a)')trim(dir)//'/A_'//cmat//'.dat'
42: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
43: call MatLoad(A,viewer,ierr);CHKERRA(ierr)
44: call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
45: if (i .eq. 0) then
46: call MatCreateVecs(A,x,b,ierr);CHKERRA(ierr)
47: endif
48: write (name,'(a)')trim(dir)//'/rhs_'//cmat//'.dat'
49: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
50: call VecLoad(b,viewer,ierr);CHKERRA(ierr)
51: call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
52: call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
53: call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
54: call KSPSolve(ksp,b,x,ierr);CHKERRA(ierr)
55: #if defined(PETSC_HAVE_HPDDM)
56: call PetscObjectTypeCompare(ksp,KSPHPDDM,flg,ierr);CHKERRA(ierr)
57: if (flg .and. reset) then
58: call KSPHPDDMGetDeflationSpace(ksp,U,ierr);CHKERRA(ierr)
59: call KSPReset(ksp,ierr);CHKERRA(ierr)
60: call KSPSetOperators(ksp,A,A,ierr);CHKERRA(ierr)
61: call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
62: call KSPSetUp(ksp,ierr);CHKERRA(ierr)
63: if (U .ne. PETSC_NULL_MAT) then
64: call KSPHPDDMSetDeflationSpace(ksp,U,ierr);CHKERRA(ierr)
65: call MatDestroy(U,ierr);CHKERRA(ierr)
66: endif
67: endif
68: #endif
69: 50 continue
70: call VecDestroy(x,ierr);CHKERRA(ierr)
71: call VecDestroy(b,ierr);CHKERRA(ierr)
72: call MatDestroy(A,ierr);CHKERRA(ierr)
73: call KSPDestroy(ksp,ierr);CHKERRA(ierr)
74: call PetscFinalize(ierr)
75: end
77: !/*TEST
78: !
79: ! test:
80: ! suffix: 1
81: ! output_file: output/ex75_1.out
82: ! nsize: 1
83: ! requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
84: ! 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 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
85: !
86: ! test:
87: ! requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
88: ! suffix: 1_icc
89: ! output_file: output/ex75_1_icc.out
90: ! nsize: 1
91: ! 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
92: !
93: ! testset:
94: ! requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
95: ! 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 -ksp_hpddm_recycle 20 -load_dir ${DATAFILESPATH}/matrices/hpddm/GCRODR
96: ! test:
97: ! nsize: 1
98: ! suffix: 2_seq
99: ! output_file: output/ex75_2.out
100: ! test:
101: ! nsize: 2
102: ! suffix: 2_par
103: ! output_file: output/ex75_2.out
104: !
105: ! test:
106: ! requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
107: ! suffix: 2_icc
108: ! output_file: output/ex75_2_icc.out
109: ! nsize: 1
110: ! 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
111: !
112: !TEST*/