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;
12: MatLoad(data_mat, inp_viewer);
13: MatAssemblyBegin(data_mat, MAT_FINAL_ASSEMBLY);
14: MatAssemblyEnd(data_mat, MAT_FINAL_ASSEMBLY);
15: MatViewFromOptions(data_mat, NULL, "-view_mat");
17: MatGetSize(data_mat, &M, &N);
18: PetscPrintf(PETSC_COMM_WORLD, "Data matrix size: %D %D\n", M,N);
20: /* compute matrix norms */
21: MatNorm(data_mat, NORM_1, &norms[0]);
22: MatNorm(data_mat, NORM_INFINITY, &norms[1]);
23: MatNorm(data_mat, NORM_FROBENIUS, &norms[2]);
24: PetscPrintf(PETSC_COMM_WORLD, "Data matrix norms: %g %g %g\n", (double)norms[0],(double)norms[1],(double)norms[2]);
26: /* compute autocorrelation matrix */
27: MatMatTransposeMult(data_mat, data_mat, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &corr_mat);
29: /* compute autocorrelation matrix norms */
30: MatNorm(corr_mat, NORM_1, &norms[3]);
31: MatNorm(corr_mat, NORM_INFINITY, &norms[4]);
32: MatNorm(corr_mat, NORM_FROBENIUS, &norms[5]);
33: PetscPrintf(PETSC_COMM_WORLD, "Autocorrelation matrix norms: %g %g %g\n", (double)norms[3],(double)norms[4],(double)norms[5]);
35: MatDestroy(&corr_mat);
36: return(0);
37: }
39: static PetscErrorCode GetReader(MPI_Comm comm, const char option[], PetscViewer *r, PetscViewerFormat *fmt)
40: {
41: PetscBool flg;
45: PetscOptionsGetViewer(PETSC_COMM_SELF, NULL, NULL, option, r, fmt, &flg);
46: if (flg) {
47: PetscFileMode mode;
48: PetscViewerFileGetMode(*r, &mode);
49: flg = (PetscBool) (mode == FILE_MODE_READ);
50: }
51: if (!flg) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Need to specify %s viewer_type:file:format:read", option);
52: return(0);
53: }
55: int main(int argc, char **argv)
56: {
57: PetscErrorCode ierr;
58: PetscInt i;
59: PetscReal norms0[NNORMS], norms1[NNORMS];
60: PetscViewer inp_viewer;
61: PetscViewerFormat fmt;
62: Mat data_mat;
63: char mat_name[PETSC_MAX_PATH_LEN]="dmatrix";
65: PetscInitialize(&argc, &argv, NULL, NULL);if (ierr) return ierr;
66: PetscOptionsGetString(NULL,NULL,"-mat_name",mat_name,sizeof(mat_name),NULL);
68: /* load matrix sequentially */
69: MatCreate(PETSC_COMM_SELF, &data_mat);
70: MatSetType(data_mat,MATDENSE);
71: PetscObjectSetName((PetscObject)data_mat, mat_name);
72: GetReader(PETSC_COMM_SELF, "-serial_reader", &inp_viewer, &fmt);
73: PetscViewerPushFormat(inp_viewer, fmt);
74: MatLoadComputeNorms(data_mat, inp_viewer, norms0);
75: PetscViewerPopFormat(inp_viewer);
76: PetscViewerDestroy(&inp_viewer);
77: MatViewFromOptions(data_mat, NULL, "-view_serial_mat");
78: MatDestroy(&data_mat);
80: /* load matrix in parallel */
81: MatCreate(PETSC_COMM_WORLD, &data_mat);
82: MatSetType(data_mat,MATDENSE);
83: PetscObjectSetName((PetscObject)data_mat, mat_name);
84: GetReader(PETSC_COMM_WORLD, "-parallel_reader", &inp_viewer, &fmt);
85: PetscViewerPushFormat(inp_viewer, fmt);
86: MatLoadComputeNorms(data_mat, inp_viewer, norms1);
87: PetscViewerPopFormat(inp_viewer);
88: PetscViewerDestroy(&inp_viewer);
89: MatViewFromOptions(data_mat, NULL, "-view_parallel_mat");
90: MatDestroy(&data_mat);
92: for (i=0; i<NNORMS; i++) {
93: if (PetscAbs(norms0[i] - norms1[i]) > PETSC_SMALL) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_PLIB, "norm0[%D] = %g != %g = norms1[%D]", i, (double)norms0[i], (double)norms1[i], i);
94: }
96: PetscFinalize();
97: return ierr;
98: }
100: /*TEST
102: test:
103: suffix: 1
104: requires: hdf5 datafilespath complex
105: args: -serial_reader hdf5:${DATAFILESPATH}/matrices/hdf5/sample_data.h5::read -parallel_reader hdf5:${DATAFILESPATH}/matrices/hdf5/sample_data.h5::read
106: nsize: {{1 2 4}}
108: test:
109: requires: hdf5 datafilespath
110: args: -serial_reader hdf5:${DATAFILESPATH}/matrices/hdf5/tiny_rectangular_mat.h5::read -parallel_reader hdf5:${DATAFILESPATH}/matrices/hdf5/tiny_rectangular_mat.h5::read
111: nsize: {{1 2}}
112: test:
113: suffix: 2-complex
114: requires: complex
115: args: -mat_name ComplexMat
116: test:
117: suffix: 2-real
118: requires: !complex
119: args: -mat_name RealMat
121: TEST*/