Actual source code: ex10.c
2: static char help[] = "Reads a PETSc matrix and computes the 2 norm of the columns\n\n";
4: /*
5: Include "petscmat.h" so that we can use matrices.
6: automatically includes:
7: petscsys.h - base PETSc routines petscvec.h - vectors
8: petscmat.h - matrices
9: petscis.h - index sets petscviewer.h - viewers
10: */
11: #include <petscmat.h>
13: int main(int argc, char **args)
14: {
15: Mat A; /* matrix */
16: PetscViewer fd; /* viewer */
17: char file[PETSC_MAX_PATH_LEN]; /* input file name */
18: PetscReal *norms;
19: PetscInt n, cstart, cend;
20: PetscBool flg;
21: PetscViewerFormat format;
24: PetscInitialize(&argc, &args, (char *)0, help);
25: /*
26: Determine files from which we read the matrix
27: */
28: PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg);
31: /*
32: Open binary file. Note that we use FILE_MODE_READ to indicate
33: reading from this file.
34: */
35: PetscViewerCreate(PETSC_COMM_WORLD, &fd);
36: PetscViewerSetType(fd, PETSCVIEWERBINARY);
37: PetscViewerSetFromOptions(fd);
38: PetscOptionsGetEnum(NULL, NULL, "-viewer_format", PetscViewerFormats, (PetscEnum *)&format, &flg);
39: if (flg) PetscViewerPushFormat(fd, format);
40: PetscViewerFileSetMode(fd, FILE_MODE_READ);
41: PetscViewerFileSetName(fd, file);
43: /*
44: Load the matrix; then destroy the viewer.
45: Matrix type is set automatically but you can override it by MatSetType() prior to MatLoad().
46: Do that only if you really insist on the given type.
47: */
48: MatCreate(PETSC_COMM_WORLD, &A);
49: MatSetOptionsPrefix(A, "a_");
50: PetscObjectSetName((PetscObject)A, "A");
51: MatSetFromOptions(A);
52: MatLoad(A, fd);
53: PetscViewerDestroy(&fd);
55: MatGetSize(A, NULL, &n);
56: MatGetOwnershipRangeColumn(A, &cstart, &cend);
57: PetscMalloc1(n, &norms);
58: MatGetColumnNorms(A, NORM_2, norms);
59: PetscRealView(cend - cstart, norms + cstart, PETSC_VIEWER_STDOUT_WORLD);
60: PetscFree(norms);
62: PetscObjectPrintClassNamePrefixType((PetscObject)A, PETSC_VIEWER_STDOUT_WORLD);
63: MatViewFromOptions(A, NULL, "-mat_view");
65: MatDestroy(&A);
66: PetscFinalize();
67: return 0;
68: }
70: /*TEST
72: test:
73: suffix: mpiaij
74: nsize: 2
75: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
76: args: -f ${DATAFILESPATH}/matrices/small -a_mat_type mpiaij
77: args: -a_matload_symmetric
79: test:
80: suffix: mpiaij_hdf5
81: nsize: 2
82: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
83: args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat
84: args: -a_matload_symmetric
86: test:
87: suffix: mpiaij_rect_hdf5
88: nsize: 2
89: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
90: args: -f ${DATAFILESPATH}/matrices/matlab/small_rect.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat
92: test:
93: # test for more processes than rows
94: suffix: mpiaij_hdf5_tiny
95: nsize: 8
96: requires: double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
97: args: -f ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system_with_x0.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat
98: args: -a_matload_symmetric
100: test:
101: # test for more processes than rows, complex
102: TODO: not yet implemented for MATLAB complex format
103: suffix: mpiaij_hdf5_tiny_complex
104: nsize: 8
105: requires: double complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
106: args: -f ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system_with_x0_complex.mat -a_mat_type mpiaij -viewer_type hdf5 -viewer_format hdf5_mat
107: args: -a_matload_symmetric
109: test:
110: TODO: mpibaij not supported yet
111: suffix: mpibaij_hdf5
112: nsize: 2
113: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
114: args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type mpibaij -a_mat_block_size 2 -viewer_type hdf5 -viewer_format hdf5_mat
115: args: -a_matload_symmetric
117: test:
118: suffix: mpidense
119: nsize: 2
120: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
121: args: -f ${DATAFILESPATH}/matrices/small -a_mat_type mpidense
122: args: -a_matload_symmetric
124: test:
125: suffix: seqaij
126: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
127: args: -f ${DATAFILESPATH}/matrices/small -a_mat_type seqaij
128: args: -a_matload_symmetric
130: test:
131: suffix: seqaij_hdf5
132: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
133: args: -f ${DATAFILESPATH}/matrices/matlab/small.mat -a_mat_type seqaij -viewer_type hdf5 -viewer_format hdf5_mat
134: args: -a_matload_symmetric
136: test:
137: suffix: seqaij_rect_hdf5
138: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
139: args: -f ${DATAFILESPATH}/matrices/matlab/small_rect.mat -a_mat_type seqaij -viewer_type hdf5 -viewer_format hdf5_mat
141: test:
142: suffix: seqdense
143: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
144: args: -f ${DATAFILESPATH}/matrices/small -a_mat_type seqdense
145: args: -a_matload_symmetric
147: test:
148: suffix: seqdense_hdf5
149: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
150: args: -f ${DATAFILESPATH}/matrices/matlab/small_dense.mat -a_mat_type seqdense -viewer_type hdf5 -viewer_format hdf5_mat
151: args: -a_matload_symmetric
153: test:
154: suffix: seqdense_rect_hdf5
155: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
156: args: -f ${DATAFILESPATH}/matrices/matlab/small_rect_dense.mat -a_mat_type seqdense -viewer_type hdf5 -viewer_format hdf5_mat
158: test:
159: suffix: mpidense_hdf5
160: nsize: 2
161: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
162: args: -f ${DATAFILESPATH}/matrices/matlab/small_dense.mat -a_mat_type mpidense -viewer_type hdf5 -viewer_format hdf5_mat
163: args: -a_matload_symmetric
165: test:
166: suffix: mpidense_rect_hdf5
167: nsize: 2
168: requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) hdf5 defined(PETSC_HDF5_HAVE_ZLIB)
169: args: -f ${DATAFILESPATH}/matrices/matlab/small_rect_dense.mat -a_mat_type mpidense -viewer_type hdf5 -viewer_format hdf5_mat
170: TEST*/