Actual source code: ex72.c
2: #if !defined(PETSC_USE_COMPLEX)
4: static char help[] = "Reads in a Symmetric matrix in MatrixMarket format. Writes\n\
5: it using the PETSc sparse format. It also adds a Vector set to random values to the\n\
6: output file. Input parameters are:\n\
7: -fin <filename> : input file\n\
8: -fout <filename> : output file\n\n";
10: #include petscmat.h
14: int main(int argc,char **args)
15: {
16: Mat A;
17: Vec b;
18: char filein[PETSC_MAX_PATH_LEN],fileout[PETSC_MAX_PATH_LEN],buf[PETSC_MAX_PATH_LEN];
19: PetscInt i,m,n,nnz,col,row;
21: PetscMPIInt size;
22: PetscScalar val;
23: FILE* file;
24: PetscViewer view;
25: PetscRandom r;
27: PetscInitialize(&argc,&args,(char *)0,help);
29: MPI_Comm_size(PETSC_COMM_WORLD,&size);
30: if (size > 1) SETERRQ(1,"Uniprocessor Example only\n");
32: /* Read in matrix and RHS */
33: PetscOptionsGetString(PETSC_NULL,"-fin",filein,PETSC_MAX_PATH_LEN-1,PETSC_NULL);
34: PetscFOpen(PETSC_COMM_SELF,filein,"r",&file);
36: /* Ignore the first line */
37: /* while (getc(file) != '\n'); */
38: fgets(buf,PETSC_MAX_PATH_LEN-1,file);
39: printf("%s",buf);
40: fscanf(file,"%d %d %d\n",&m,&n,&nnz);
41: printf ("m = %d, n = %d, nnz = %d\n",m,n,nnz);
43: MatCreateSeqAIJ(PETSC_COMM_WORLD,m,n,20,0,&A);
44: VecCreate(PETSC_COMM_WORLD,&b);
45: VecSetSizes(b,PETSC_DECIDE,n);
46: VecSetFromOptions(b);
47: PetscRandomCreate(PETSC_COMM_SELF,RANDOM_DEFAULT,&r);
48: VecSetRandom(b,r);
50: for (i=0; i<nnz; i++) {
51: fscanf(file,"%d %d %le\n",&row,&col,(double*)&val);
52: row = row-1; col = col-1 ;
53: MatSetValues(A,1,&row,1,&col,&val,INSERT_VALUES);
54: if (row != col) {
55: MatSetValues(A,1,&col,1,&row,&val,INSERT_VALUES);
56: }
57: }
58: fclose(file);
60: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
61: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
63: PetscPrintf(PETSC_COMM_SELF,"Reading matrix completes.\n");
64: PetscOptionsGetString(PETSC_NULL,"-fout",fileout,PETSC_MAX_PATH_LEN-1,PETSC_NULL);
65: PetscViewerBinaryOpen(PETSC_COMM_WORLD,fileout,FILE_MODE_WRITE,&view);
66: MatView(A,view);
67: VecView(b,view);
68: PetscViewerDestroy(view);
70: VecDestroy(b);
71: MatDestroy(A);
72: PetscRandomDestroy(r);
74: PetscFinalize();
75: return 0;
76: }
77: #else
78: #include <stdio.h>
79: int main(int argc,char **args)
80: {
81: fprintf(stdout,"This example does not work for complex numbers.\n");
82: return 0;
83: }
84: #endif