Actual source code: ex84.c
1: #include <petscmat.h>
3: #define NNORMS 6
5: static PetscErrorCode MatLoadComputeNorms(Mat data_mat, PetscViewer inp_viewer, PetscReal norms[])
6: {
7: Mat corr_mat;
8: PetscInt M,N;
10: MatLoad(data_mat, inp_viewer);
11: MatAssemblyBegin(data_mat, MAT_FINAL_ASSEMBLY);
12: MatAssemblyEnd(data_mat, MAT_FINAL_ASSEMBLY);
13: MatViewFromOptions(data_mat, NULL, "-view_mat");
15: MatGetSize(data_mat, &M, &N);
16: PetscPrintf(PETSC_COMM_WORLD, "Data matrix size: %" PetscInt_FMT " %" PetscInt_FMT "\n", M,N);
18: /* compute matrix norms */
19: MatNorm(data_mat, NORM_1, &norms[0]);
20: MatNorm(data_mat, NORM_INFINITY, &norms[1]);
21: MatNorm(data_mat, NORM_FROBENIUS, &norms[2]);
22: PetscPrintf(PETSC_COMM_WORLD, "Data matrix norms: %g %g %g\n", (double)norms[0],(double)norms[1],(double)norms[2]);
24: /* compute autocorrelation matrix */
25: MatMatTransposeMult(data_mat, data_mat, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &corr_mat);
27: /* compute autocorrelation matrix norms */
28: MatNorm(corr_mat, NORM_1, &norms[3]);
29: MatNorm(corr_mat, NORM_INFINITY, &norms[4]);
30: MatNorm(corr_mat, NORM_FROBENIUS, &norms[5]);
31: PetscPrintf(PETSC_COMM_WORLD, "Autocorrelation matrix norms: %g %g %g\n", (double)norms[3],(double)norms[4],(double)norms[5]);
33: MatDestroy(&corr_mat);
34: return 0;
35: }
37: static PetscErrorCode GetReader(MPI_Comm comm, const char option[], PetscViewer *r, PetscViewerFormat *fmt)
38: {
39: PetscBool flg;
41: PetscOptionsGetViewer(PETSC_COMM_SELF, NULL, NULL, option, r, fmt, &flg);
42: if (flg) {
43: PetscFileMode mode;
44: PetscViewerFileGetMode(*r, &mode);
45: flg = (PetscBool) (mode == FILE_MODE_READ);
46: }
48: return 0;
49: }
51: int main(int argc, char **argv)
52: {
53: PetscInt i;
54: PetscReal norms0[NNORMS], norms1[NNORMS];
55: PetscViewer inp_viewer;
56: PetscViewerFormat fmt;
57: Mat data_mat;
58: char mat_name[PETSC_MAX_PATH_LEN]="dmatrix";
60: PetscInitialize(&argc, &argv, NULL, NULL);
61: PetscOptionsGetString(NULL,NULL,"-mat_name",mat_name,sizeof(mat_name),NULL);
63: /* load matrix sequentially */
64: MatCreate(PETSC_COMM_SELF, &data_mat);
65: MatSetType(data_mat,MATDENSE);
66: PetscObjectSetName((PetscObject)data_mat, mat_name);
67: GetReader(PETSC_COMM_SELF, "-serial_reader", &inp_viewer, &fmt);
68: PetscViewerPushFormat(inp_viewer, fmt);
69: MatLoadComputeNorms(data_mat, inp_viewer, norms0);
70: PetscViewerPopFormat(inp_viewer);
71: PetscViewerDestroy(&inp_viewer);
72: MatViewFromOptions(data_mat, NULL, "-view_serial_mat");
73: MatDestroy(&data_mat);
75: /* load matrix in parallel */
76: MatCreate(PETSC_COMM_WORLD, &data_mat);
77: MatSetType(data_mat,MATDENSE);
78: PetscObjectSetName((PetscObject)data_mat, mat_name);
79: GetReader(PETSC_COMM_WORLD, "-parallel_reader", &inp_viewer, &fmt);
80: PetscViewerPushFormat(inp_viewer, fmt);
81: MatLoadComputeNorms(data_mat, inp_viewer, norms1);
82: PetscViewerPopFormat(inp_viewer);
83: PetscViewerDestroy(&inp_viewer);
84: MatViewFromOptions(data_mat, NULL, "-view_parallel_mat");
85: MatDestroy(&data_mat);
87: for (i=0; i<NNORMS; i++) {
89: }
91: PetscFinalize();
92: return 0;
93: }
95: /*TEST
97: test:
98: suffix: 1
99: requires: hdf5 datafilespath complex
100: args: -serial_reader hdf5:${DATAFILESPATH}/matrices/hdf5/sample_data.h5::read -parallel_reader hdf5:${DATAFILESPATH}/matrices/hdf5/sample_data.h5::read
101: nsize: {{1 2 4}}
103: test:
104: requires: hdf5 datafilespath
105: args: -serial_reader hdf5:${DATAFILESPATH}/matrices/hdf5/tiny_rectangular_mat.h5::read -parallel_reader hdf5:${DATAFILESPATH}/matrices/hdf5/tiny_rectangular_mat.h5::read
106: nsize: {{1 2}}
107: test:
108: suffix: 2-complex
109: requires: complex
110: args: -mat_name ComplexMat
111: test:
112: suffix: 2-real
113: requires: !complex
114: args: -mat_name RealMat
116: TEST*/