Actual source code: ex43.c
1: static char help[] = "Reads a PETSc matrix from a file and solves a linear system \n\
2: using the aijcusparse class. Input parameters are:\n\
3: -f <input_file> : the file to load\n\n";
5: /*
6: This code can be used to test PETSc interface to other packages.\n\
7: Examples of command line options: \n\
8: ./ex43 -f DATAFILESPATH/matrices/cfd.2.10 -mat_cusparse_mult_storage_format ell \n\
9: ./ex43 -f DATAFILESPATH/matrices/shallow_water1 -ksp_type cg -pc_type icc -mat_cusparse_mult_storage_format ell \n\
10: \n\n";
11: */
13: #include <petscksp.h>
15: int main(int argc, char **argv)
16: {
17: KSP ksp;
18: Mat A;
19: Vec X, B;
20: PetscInt m, its;
21: PetscReal norm;
22: char file[PETSC_MAX_PATH_LEN];
23: PetscBool flg;
24: PetscViewer fd;
26: PetscFunctionBeginUser;
27: PetscCall(PetscInitialize(&argc, &argv, 0, help));
28: /* Load the data from a file */
29: PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
30: PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER_INPUT, "Must indicate binary file with the -f option");
31: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd));
33: /* Build the matrix */
34: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
35: PetscCall(MatSetFromOptions(A));
36: PetscCall(MatLoad(A, fd));
38: /* Build the vectors */
39: PetscCall(MatGetLocalSize(A, &m, NULL));
40: PetscCall(VecCreate(PETSC_COMM_WORLD, &B));
41: PetscCall(VecSetSizes(B, m, PETSC_DECIDE));
42: PetscCall(VecCreate(PETSC_COMM_WORLD, &X));
43: PetscCall(VecSetSizes(X, m, PETSC_DECIDE));
44: PetscCall(VecSetFromOptions(B));
45: PetscCall(VecSetFromOptions(X));
46: PetscCall(VecSet(B, 1.0));
48: /* Build the KSP */
49: PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
50: PetscCall(KSPSetOperators(ksp, A, A));
51: PetscCall(KSPSetType(ksp, KSPGMRES));
52: PetscCall(KSPSetTolerances(ksp, 1.0e-12, PETSC_CURRENT, PETSC_CURRENT, 100));
53: PetscCall(KSPSetFromOptions(ksp));
55: /* Solve */
56: PetscCall(KSPSolve(ksp, B, X));
58: /* print out norm and the number of iterations */
59: PetscCall(KSPGetIterationNumber(ksp, &its));
60: PetscCall(KSPGetResidualNorm(ksp, &norm));
61: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of iterations = %3" PetscInt_FMT "\n", its));
62: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Residual norm %1.5g\n", (double)norm));
64: /* Cleanup */
65: PetscCall(VecDestroy(&X));
66: PetscCall(VecDestroy(&B));
67: PetscCall(MatDestroy(&A));
68: PetscCall(KSPDestroy(&ksp));
69: PetscCall(PetscViewerDestroy(&fd));
70: PetscCall(PetscFinalize());
71: return 0;
72: }
74: /*TEST
76: test:
77: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !CUDA_VERSION_11PLUS
78: args: -f ${DATAFILESPATH}/matrices/cfd.2.10 -mat_type seqaijcusparse -pc_factor_mat_solver_type cusparse -mat_cusparse_storage_format ell -vec_type cuda -pc_type ilu
80: test:
81: suffix: 2
82: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !CUDA_VERSION_11PLUS
83: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type seqaijcusparse -pc_factor_mat_solver_type cusparse -mat_cusparse_storage_format hyb -vec_type cuda -ksp_type cg -pc_type icc
85: testset:
86: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
87: args: -f ${DATAFILESPATH}/matrices/cfd.2.10 -ksp_type bicg -pc_type ilu
89: test:
90: suffix: 3
91: requires: cuda
92: args: -mat_type seqaijcusparse -pc_factor_mat_solver_type cusparse -mat_cusparse_storage_format csr -vec_type cuda
93: test:
94: suffix: 4
95: requires: cuda
96: args: -mat_type seqaijcusparse -pc_factor_mat_solver_type cusparse -mat_cusparse_storage_format csr -vec_type cuda -pc_factor_mat_ordering_type nd
97: test: # Test MatSolveTranspose
98: suffix: 3_kokkos
99: requires: kokkos_kernels
100: args: -mat_type seqaijkokkos -pc_factor_mat_solver_type kokkos -vec_type kokkos
101: output_file: output/ex43_3.out
103: testset:
104: nsize: 2
105: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !CUDA_VERSION_11PLUS
106: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type mpiaijcusparse -mat_cusparse_mult_diag_storage_format hyb -pc_type none -vec_type cuda
107: test:
108: suffix: 5
109: args: -use_gpu_aware_mpi 0
110: test:
111: suffix: 5_gpu_aware_mpi
112: output_file: output/ex43_5.out
114: test:
115: suffix: 6
116: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
117: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type seqaijcusparse -pc_type none -vec_type cuda
119: testset:
120: nsize: 2
121: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
122: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type mpiaijcusparse -pc_type none -vec_type cuda
124: test:
125: suffix: 7
126: args: -use_gpu_aware_mpi 0
127: test:
128: suffix: 7_gpu_aware_mpi
129: output_file: output/ex43_7.out
131: test:
132: suffix: 8
133: requires: viennacl datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
134: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type seqaijviennacl -pc_type none -vec_type viennacl
135: output_file: output/ex43_6.out
137: test:
138: suffix: 9
139: nsize: 2
140: requires: viennacl datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
141: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type mpiaijviennacl -pc_type none -vec_type viennacl
142: output_file: output/ex43_7.out
144: test:
145: suffix: 10
146: nsize: 2
147: requires: kokkos_kernels datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
148: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type aijkokkos -vec_type kokkos
150: TEST*/