Actual source code: ex77f.F90
1: !
2: ! Description: Solves a linear system with a block of right-hand sides using KSPHPDDM.
3: !
5: program main
6: #include <petsc/finclude/petscksp.h>
7: use petscksp
8: implicit none
9: Mat X,B
10: Mat A
11: KSP ksp
12: PC pc
13: Mat F
14: PetscScalar alpha
15: PetscReal norm
16: PetscInt m,K
17: PetscViewer viewer
18: character*(PETSC_MAX_PATH_LEN) name
19: PetscBool flg
20: PetscErrorCode ierr
22: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
23: if (ierr .ne. 0) then
24: print *,'Unable to initialize PETSc'
25: stop
26: endif
27: call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-f',name,flg,ierr);CHKERRA(ierr)
28: if (flg .eqv. PETSC_FALSE) then
29: SETERRA(PETSC_COMM_WORLD,PETSC_ERR_SUP,'Must provide a binary file for the matrix')
30: endif
31: K = 5
32: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',K,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: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,viewer,ierr);CHKERRA(ierr)
37: call MatLoad(A,viewer,ierr);CHKERRA(ierr)
38: call PetscViewerDestroy(viewer,ierr);CHKERRA(ierr)
39: call MatGetLocalSize(A,m,PETSC_NULL_INTEGER,ierr);CHKERRA(ierr)
40: call MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,K,PETSC_NULL_SCALAR,B,ierr);CHKERRA(ierr)
41: call MatCreateDense(PETSC_COMM_WORLD,m,PETSC_DECIDE,PETSC_DECIDE,K,PETSC_NULL_SCALAR,X,ierr);CHKERRA(ierr)
42: call MatSetRandom(B,PETSC_NULL_RANDOM,ierr);CHKERRA(ierr)
43: call KSPSetFromOptions(ksp,ierr);CHKERRA(ierr)
44: call KSPSetUp(ksp,ierr);CHKERRA(ierr)
45: call KSPMatSolve(ksp,B,X,ierr);CHKERRA(ierr)
46: call KSPGetMatSolveBatchSize(ksp,M,ierr);CHKERRA(ierr)
47: if (M .ne. PETSC_DECIDE) then
48: call KSPSetMatSolveBatchSize(ksp,PETSC_DECIDE,ierr);CHKERRA(ierr)
49: call MatZeroEntries(X,ierr);CHKERRA(ierr)
50: call KSPMatSolve(ksp,B,X,ierr);CHKERRA(ierr)
51: endif
52: call KSPGetPC(ksp,pc,ierr);CHKERRA(ierr)
53: call PetscObjectTypeCompare(pc,PCLU,flg,ierr);CHKERRA(ierr)
54: if (flg) then
55: call PCFactorGetMatrix(pc,F,ierr);CHKERRA(ierr)
56: call MatMatSolve(F,B,B,ierr);CHKERRA(ierr)
57: alpha = -1.0
58: call MatAYPX(B,alpha,X,SAME_NONZERO_PATTERN,ierr);CHKERRA(ierr)
59: call MatNorm(B,NORM_INFINITY,norm,ierr);CHKERRA(ierr)
60: if (norm > 100*PETSC_MACHINE_EPSILON) then
61: SETERRA(PETSC_COMM_WORLD,PETSC_ERR_PLIB,'KSPMatSolve() and MatMatSolve() difference has nonzero norm')
62: endif
63: endif
64: call MatDestroy(X,ierr);CHKERRA(ierr)
65: call MatDestroy(B,ierr);CHKERRA(ierr)
66: call MatDestroy(A,ierr);CHKERRA(ierr)
67: call KSPDestroy(ksp,ierr);CHKERRA(ierr)
68: call PetscFinalize(ierr)
69: end
71: !/*TEST
72: !
73: ! testset:
74: ! nsize: 2
75: ! requires: datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
76: ! args: -ksp_converged_reason -ksp_max_it 1000 -f ${DATAFILESPATH}/matrices/hpddm/GCRODR/A_400.dat
77: ! test:
78: ! suffix: 1
79: ! output_file: output/ex77_1.out
80: ! args:
81: ! test:
82: ! suffix: 2a
83: ! requires: hpddm
84: ! output_file: output/ex77_2_ksp_hpddm_type-gmres.out
85: ! args: -ksp_type hpddm -pc_type asm -ksp_hpddm_type gmres
86: ! test:
87: ! suffix: 2b
88: ! requires: hpddm
89: ! output_file: output/ex77_2_ksp_hpddm_type-bgmres.out
90: ! args: -ksp_type hpddm -pc_type asm -ksp_hpddm_type bgmres
91: ! test:
92: ! suffix: 3a
93: ! requires: hpddm
94: ! output_file: output/ex77_3_ksp_hpddm_type-gcrodr.out
95: ! args: -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type gcrodr
96: ! test:
97: ! suffix: 3b
98: ! requires: hpddm
99: ! output_file: output/ex77_3_ksp_hpddm_type-bgcrodr.out
100: ! args: -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type bgcrodr
101: ! test:
102: ! nsize: 4
103: ! suffix: 4
104: ! requires: hpddm
105: ! output_file: output/ex77_4.out
106: ! args: -ksp_rtol 1e-4 -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type bgcrodr -ksp_view_final_residual -N 12 -ksp_matsolve_batch_size 5
107: ! test:
108: ! nsize: 1
109: ! suffix: preonly
110: ! requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
111: ! output_file: output/ex77_preonly.out
112: ! args: -N 6 -f ${DATAFILESPATH}/matrices/hpddm/GCRODR/A_400.dat -pc_type lu -ksp_type hpddm -ksp_hpddm_type preonly
113: ! test:
114: ! nsize: 4
115: ! suffix: 4_slepc
116: ! requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES) slepc define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
117: ! output_file: output/ex77_4.out
118: ! filter: sed "/^ksp_hpddm_recycle_ Linear eigensolve converged/d"
119: ! args: -ksp_converged_reason -ksp_max_it 500 -f ${DATAFILESPATH}/matrices/hpddm/GCRODR/A_400.dat -ksp_rtol 1e-4 -ksp_type hpddm -ksp_hpddm_recycle 5 -ksp_hpddm_type bgcrodr -ksp_view_final_residual -N 12 -ksp_matsolve_batch_size 5 -ksp_hpddm_recycle_redistribute 2 -ksp_hpddm_recycle_mat_type dense -ksp_hpddm_recycle_eps_converged_reason -ksp_hpddm_recycle_st_pc_type redundant
120: !
121: !TEST*/