Actual source code: ex43.c
petsc-3.5.4 2015-05-23
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>
17: int main(int argc,char **argv)
18: {
19: KSP ksp;
20: Mat A;
21: Vec X,B;
22: PetscInt m, its;
23: PetscReal norm;
24: char file[PETSC_MAX_PATH_LEN];
25: PetscBool flg;
26: PetscViewer fd;
27: PetscErrorCode ierr;
29: PetscInitialize(&argc,&argv,0,help);
31: /* Load the data from a file */
32: PetscOptionsGetString(NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);
33: if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate binary file with the -f option");
34: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
36: /* Build the matrix */
37: MatCreate(PETSC_COMM_WORLD,&A);
38: MatSetFromOptions(A);
39: MatLoad(A,fd);
41: /* Build the vectors */
42: MatGetLocalSize(A,&m,NULL);
43: VecCreate(PETSC_COMM_WORLD,&B);
44: VecSetSizes(B,m,PETSC_DECIDE);
45: VecCreate(PETSC_COMM_WORLD,&X);
46: VecSetSizes(X,m,PETSC_DECIDE);
47: VecSetFromOptions(B);
48: VecSetFromOptions(X);
49: VecSet(B,1.0);
51: /* Build the KSP */
52: KSPCreate(PETSC_COMM_WORLD,&ksp);
53: KSPSetOperators(ksp,A,A);
54: KSPSetType(ksp,KSPGMRES);
55: KSPSetTolerances(ksp,1.0e-12,PETSC_DEFAULT,PETSC_DEFAULT,100);
56: KSPSetFromOptions(ksp);
58: /* Solve */
59: KSPSolve(ksp,B,X);
61: /* print out norm and the number of iterations */
62: KSPGetIterationNumber(ksp,&its);
63: KSPGetResidualNorm(ksp,&norm);
64: PetscPrintf(PETSC_COMM_WORLD,"Number of iterations = %3D\n",its);
65: PetscPrintf(PETSC_COMM_WORLD,"Residual norm %1.5G\n",norm);
67: /* Cleanup */
68: VecDestroy(&X);
69: VecDestroy(&B);
70: MatDestroy(&A);
71: KSPDestroy(&ksp);
72: PetscViewerDestroy(&fd);
73: PetscFinalize();
74: return 0;
75: }