Actual source code: ex78.c
2: static char help[] = "Reads in a matrix in ASCII MATLAB format (I,J,A), read in vectors rhs and exact_solu in ASCII format.\n\
3: Writes them using the PETSc sparse format.\n\
4: Note: I and J start at 1, not 0, use -noshift if indices in file start with zero!\n\
5: Input parameters are:\n\
6: -Ain <filename> : input matrix in ascii format\n\
7: -rhs <filename> : input rhs in ascii format\n\
8: -solu <filename> : input true solution in ascii format\n\\n";
10: /*
11: Example: ./ex78 -Ain Ain -rhs rhs -solu solu -noshift -mat_view
12: with the datafiles in the following format:
13: Ain (I and J start at 0):
14: ------------------------
15: 3 3 6
16: 0 0 1.0
17: 0 1 2.0
18: 1 0 3.0
19: 1 1 4.0
20: 1 2 5.0
21: 2 2 6.0
23: rhs
24: ---
25: 0 3.0
26: 1 12.0
27: 2 6.0
29: solu
30: ----
31: 0 1.0
32: 0 1.0
33: 0 1.0
34: */
36: #include <petscmat.h>
38: int main(int argc,char **args)
39: {
40: Mat A = NULL;
41: Vec b,u = NULL,u_tmp;
42: char Ain[PETSC_MAX_PATH_LEN],rhs[PETSC_MAX_PATH_LEN],solu[PETSC_MAX_PATH_LEN];
43: int m,n = 0,nz,dummy; /* these are fscaned so kept as int */
44: PetscInt i,col,row,shift = 1,sizes[3],nsizes;
45: PetscScalar val;
46: PetscReal res_norm;
47: FILE *Afile,*bfile,*ufile;
48: PetscViewer view;
49: PetscBool flg_A,flg_b,flg_u,flg;
50: PetscMPIInt size;
52: PetscInitialize(&argc,&args,(char*)0,help);
53: MPI_Comm_size(PETSC_COMM_WORLD,&size);
56: /* Read in matrix, rhs and exact solution from ascii files */
57: PetscOptionsGetString(NULL,NULL,"-Ain",Ain,sizeof(Ain),&flg_A);
58: PetscOptionsHasName(NULL,NULL,"-noshift",&flg);
59: if (flg) shift = 0;
60: if (flg_A) {
61: PetscPrintf(PETSC_COMM_SELF,"\n Read matrix in ascii format ...\n");
62: PetscFOpen(PETSC_COMM_SELF,Ain,"r",&Afile);
63: nsizes = 3;
64: PetscOptionsGetIntArray(NULL,NULL,"-nosizesinfile",sizes,&nsizes,&flg);
65: if (flg) {
67: m = sizes[0];
68: n = sizes[1];
69: nz = sizes[2];
70: } else {
72: }
73: PetscPrintf(PETSC_COMM_SELF,"m: %d, n: %d, nz: %d \n", m,n,nz);
75: MatCreate(PETSC_COMM_SELF,&A);
76: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);
77: MatSetFromOptions(A);
78: MatSeqAIJSetPreallocation(A,nz/m,NULL);
79: MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
81: for (i=0; i<nz; i++) {
83: row -= shift; col -= shift; /* set index set starts at 0 */
84: MatSetValues(A,1,&row,1,&col,&val,INSERT_VALUES);
85: }
86: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
87: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
88: fclose(Afile);
89: }
91: PetscOptionsGetString(NULL,NULL,"-rhs",rhs,sizeof(rhs),&flg_b);
92: if (flg_b) {
93: VecCreate(PETSC_COMM_SELF,&b);
94: VecSetSizes(b,PETSC_DECIDE,n);
95: VecSetFromOptions(b);
96: PetscPrintf(PETSC_COMM_SELF,"\n Read rhs in ascii format ...\n");
97: PetscFOpen(PETSC_COMM_SELF,rhs,"r",&bfile);
98: for (i=0; i<n; i++) {
100: VecSetValues(b,1,&i,&val,INSERT_VALUES);
101: }
102: VecAssemblyBegin(b);
103: VecAssemblyEnd(b);
104: fclose(bfile);
105: }
107: PetscOptionsGetString(NULL,NULL,"-solu",solu,sizeof(solu),&flg_u);
108: if (flg_u) {
109: VecCreate(PETSC_COMM_SELF,&u);
110: VecSetSizes(u,PETSC_DECIDE,n);
111: VecSetFromOptions(u);
112: PetscPrintf(PETSC_COMM_SELF,"\n Read exact solution in ascii format ...\n");
113: PetscFOpen(PETSC_COMM_SELF,solu,"r",&ufile);
114: for (i=0; i<n; i++) {
116: VecSetValues(u,1,&i,&val,INSERT_VALUES);
117: }
118: VecAssemblyBegin(u);
119: VecAssemblyEnd(u);
120: fclose(ufile);
121: }
123: /* Write matrix, rhs and exact solution in Petsc binary file */
124: PetscPrintf(PETSC_COMM_SELF,"\n Write matrix in binary to 'matrix.dat' ...\n");
125: PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_WRITE,&view);
126: MatView(A,view);
128: if (flg_b) { /* Write rhs in Petsc binary file */
129: PetscPrintf(PETSC_COMM_SELF,"\n Write rhs in binary to 'matrix.dat' ...\n");
130: VecView(b,view);
131: }
132: if (flg_u) {
133: PetscPrintf(PETSC_COMM_SELF,"\n Write exact solution in binary to 'matrix.dat' ...\n");
134: VecView(u,view);
135: }
137: /* Check accuracy of the data */
138: if (flg_A & flg_b & flg_u) {
139: VecDuplicate(u,&u_tmp);
140: MatMult(A,u,u_tmp);
141: VecAXPY(u_tmp,-1.0,b);
142: VecNorm(u_tmp,NORM_2,&res_norm);
143: PetscPrintf(PETSC_COMM_SELF,"\n Accuracy of the reading data: | b - A*u |_2 : %g \n",res_norm);
144: VecDestroy(&u_tmp);
145: }
147: MatDestroy(&A);
148: if (flg_b) VecDestroy(&b);
149: if (flg_u) VecDestroy(&u);
150: PetscViewerDestroy(&view);
151: PetscFinalize();
152: return 0;
153: }
155: /*TEST
157: build:
158: requires: !defined(PETSC_USE_64BIT_INDICES) double !complex
160: test:
161: requires: datafilespath
162: args: -Ain ${DATAFILESPATH}/matrices/indefinite/afiro_A.dat
164: TEST*/