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