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";

 61:   PetscInitialize(&argc, &argv, NULL, NULL);
 62:   PetscOptionsGetString(NULL, NULL, "-mat_name", mat_name, sizeof(mat_name), NULL);

 64:   /* load matrix sequentially */
 65:   MatCreate(PETSC_COMM_SELF, &data_mat);
 66:   MatSetType(data_mat, MATDENSE);
 67:   PetscObjectSetName((PetscObject)data_mat, mat_name);
 68:   GetReader(PETSC_COMM_SELF, "-serial_reader", &inp_viewer, &fmt);
 69:   PetscViewerPushFormat(inp_viewer, fmt);
 70:   MatLoadComputeNorms(data_mat, inp_viewer, norms0);
 71:   PetscViewerPopFormat(inp_viewer);
 72:   PetscViewerDestroy(&inp_viewer);
 73:   MatViewFromOptions(data_mat, NULL, "-view_serial_mat");
 74:   MatDestroy(&data_mat);

 76:   /* load matrix in parallel */
 77:   MatCreate(PETSC_COMM_WORLD, &data_mat);
 78:   MatSetType(data_mat, MATDENSE);
 79:   PetscObjectSetName((PetscObject)data_mat, mat_name);
 80:   GetReader(PETSC_COMM_WORLD, "-parallel_reader", &inp_viewer, &fmt);
 81:   PetscViewerPushFormat(inp_viewer, fmt);
 82:   MatLoadComputeNorms(data_mat, inp_viewer, norms1);
 83:   PetscViewerPopFormat(inp_viewer);
 84:   PetscViewerDestroy(&inp_viewer);
 85:   MatViewFromOptions(data_mat, NULL, "-view_parallel_mat");
 86:   MatDestroy(&data_mat);


 90:   PetscFinalize();
 91:   return 0;
 92: }

 94: /*TEST

 96:   test:
 97:     suffix: 1
 98:     requires: hdf5 datafilespath complex
 99:     args:  -serial_reader hdf5:${DATAFILESPATH}/matrices/hdf5/sample_data.h5::read -parallel_reader hdf5:${DATAFILESPATH}/matrices/hdf5/sample_data.h5::read
100:     nsize: {{1 2 4}}

102:   test:
103:     requires: hdf5 datafilespath
104:     args:  -serial_reader hdf5:${DATAFILESPATH}/matrices/hdf5/tiny_rectangular_mat.h5::read -parallel_reader hdf5:${DATAFILESPATH}/matrices/hdf5/tiny_rectangular_mat.h5::read
105:     nsize: {{1 2}}
106:     test:
107:       suffix: 2-complex
108:       requires: complex
109:       args: -mat_name ComplexMat
110:     test:
111:       suffix: 2-real
112:       requires: !complex
113:       args: -mat_name RealMat

115: TEST*/