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: PetscInitialize(&argc,&argv,0,help);
27: /* Load the data from a file */
28: PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg);
30: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
32: /* Build the matrix */
33: MatCreate(PETSC_COMM_WORLD,&A);
34: MatSetFromOptions(A);
35: MatLoad(A,fd);
37: /* Build the vectors */
38: MatGetLocalSize(A,&m,NULL);
39: VecCreate(PETSC_COMM_WORLD,&B);
40: VecSetSizes(B,m,PETSC_DECIDE);
41: VecCreate(PETSC_COMM_WORLD,&X);
42: VecSetSizes(X,m,PETSC_DECIDE);
43: VecSetFromOptions(B);
44: VecSetFromOptions(X);
45: VecSet(B,1.0);
47: /* Build the KSP */
48: KSPCreate(PETSC_COMM_WORLD,&ksp);
49: KSPSetOperators(ksp,A,A);
50: KSPSetType(ksp,KSPGMRES);
51: KSPSetTolerances(ksp,1.0e-12,PETSC_DEFAULT,PETSC_DEFAULT,100);
52: KSPSetFromOptions(ksp);
54: /* Solve */
55: KSPSolve(ksp,B,X);
57: /* print out norm and the number of iterations */
58: KSPGetIterationNumber(ksp,&its);
59: KSPGetResidualNorm(ksp,&norm);
60: PetscPrintf(PETSC_COMM_WORLD,"Number of iterations = %3D\n",its);
61: PetscPrintf(PETSC_COMM_WORLD,"Residual norm %1.5g\n",norm);
63: /* Cleanup */
64: VecDestroy(&X);
65: VecDestroy(&B);
66: MatDestroy(&A);
67: KSPDestroy(&ksp);
68: PetscViewerDestroy(&fd);
69: PetscFinalize();
70: return 0;
71: }
73: /*TEST
75: test:
76: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !CUDA_VERSION_11PLUS
77: 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
79: test:
80: suffix: 2
81: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !CUDA_VERSION_11PLUS
82: 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
84: testset:
85: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
86: args: -f ${DATAFILESPATH}/matrices/cfd.2.10 -ksp_type bicg -pc_type ilu
88: test:
89: suffix: 3
90: requires: cuda
91: args: -mat_type seqaijcusparse -pc_factor_mat_solver_type cusparse -mat_cusparse_storage_format csr -vec_type cuda
92: test:
93: suffix: 4
94: requires: cuda
95: args: -mat_type seqaijcusparse -pc_factor_mat_solver_type cusparse -mat_cusparse_storage_format csr -vec_type cuda -pc_factor_mat_ordering_type nd
96: test: # Test MatSolveTranspose
97: suffix: 3_kokkos
98: requires: !sycl kokkos_kernels
99: args: -mat_type seqaijkokkos -pc_factor_mat_solver_type kokkos -vec_type kokkos
100: output_file: output/ex43_3.out
102: testset:
103: nsize: 2
104: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) !CUDA_VERSION_11PLUS
105: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type mpiaijcusparse -mat_cusparse_mult_diag_storage_format hyb -pc_type none -vec_type cuda
106: test:
107: suffix: 5
108: args: -use_gpu_aware_mpi 0
109: test:
110: suffix: 5_gpu_aware_mpi
111: output_file: output/ex43_5.out
113: test:
114: suffix: 6
115: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
116: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type seqaijcusparse -pc_type none -vec_type cuda
118: testset:
119: nsize: 2
120: requires: cuda datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
121: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type mpiaijcusparse -pc_type none -vec_type cuda
123: test:
124: suffix: 7
125: args: -use_gpu_aware_mpi 0
126: test:
127: suffix: 7_gpu_aware_mpi
128: output_file: output/ex43_7.out
130: test:
131: suffix: 8
132: requires: viennacl datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
133: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type seqaijviennacl -pc_type none -vec_type viennacl
134: output_file: output/ex43_6.out
136: test:
137: suffix: 9
138: nsize: 2
139: requires: viennacl datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
140: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type mpiaijviennacl -pc_type none -vec_type viennacl
141: output_file: output/ex43_7.out
143: test:
144: suffix: 10
145: nsize: 2
146: requires: !sycl kokkos_kernels datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
147: args: -f ${DATAFILESPATH}/matrices/shallow_water1 -mat_type aijkokkos -vec_type kokkos
149: TEST*/