Actual source code: ex49.c


  2: static char help[] = "Tests MatTranspose(), MatNorm(), and MatAXPY().\n\n";

  4: #include <petscmat.h>

  6: int main(int argc, char **argv)
  7: {
  8:   Mat         mat, tmat = 0;
  9:   PetscInt    m = 4, n, i, j;
 10:   PetscMPIInt size, rank;
 11:   PetscInt    rstart, rend, rect = 0;
 12:   PetscBool   flg;
 13:   PetscScalar v;
 14:   PetscReal   normf, normi, norm1;
 15:   MatInfo     info;

 18:   PetscInitialize(&argc, &argv, (char *)0, help);
 19:   PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);
 20:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 21:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 22:   n = m;
 23:   PetscOptionsHasName(NULL, NULL, "-rect1", &flg);
 24:   if (flg) {
 25:     n += 2;
 26:     rect = 1;
 27:   }
 28:   PetscOptionsHasName(NULL, NULL, "-rect2", &flg);
 29:   if (flg) {
 30:     n -= 2;
 31:     rect = 1;
 32:   }

 34:   /* Create and assemble matrix */
 35:   MatCreate(PETSC_COMM_WORLD, &mat);
 36:   MatSetSizes(mat, PETSC_DECIDE, PETSC_DECIDE, m, n);
 37:   MatSetFromOptions(mat);
 38:   MatSetUp(mat);
 39:   MatGetOwnershipRange(mat, &rstart, &rend);
 40:   for (i = rstart; i < rend; i++) {
 41:     for (j = 0; j < n; j++) {
 42:       v = 10 * i + j;
 43:       MatSetValues(mat, 1, &i, 1, &j, &v, INSERT_VALUES);
 44:     }
 45:   }
 46:   MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY);
 47:   MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY);

 49:   /* Print info about original matrix */
 50:   MatGetInfo(mat, MAT_GLOBAL_SUM, &info);
 51:   PetscPrintf(PETSC_COMM_WORLD, "original matrix nonzeros = %" PetscInt_FMT ", allocated nonzeros = %" PetscInt_FMT "\n", (PetscInt)info.nz_used, (PetscInt)info.nz_allocated);
 52:   MatNorm(mat, NORM_FROBENIUS, &normf);
 53:   MatNorm(mat, NORM_1, &norm1);
 54:   MatNorm(mat, NORM_INFINITY, &normi);
 55:   PetscPrintf(PETSC_COMM_WORLD, "original: Frobenious norm = %g, one norm = %g, infinity norm = %g\n", (double)normf, (double)norm1, (double)normi);
 56:   MatView(mat, PETSC_VIEWER_STDOUT_WORLD);

 58:   /* Form matrix transpose */
 59:   PetscOptionsHasName(NULL, NULL, "-in_place", &flg);
 60:   if (flg) {
 61:     MatTranspose(mat, MAT_INPLACE_MATRIX, &mat); /* in-place transpose */
 62:     tmat = mat;
 63:     mat  = 0;
 64:   } else { /* out-of-place transpose */
 65:     MatTranspose(mat, MAT_INITIAL_MATRIX, &tmat);
 66:   }

 68:   /* Print info about transpose matrix */
 69:   MatGetInfo(tmat, MAT_GLOBAL_SUM, &info);
 70:   PetscPrintf(PETSC_COMM_WORLD, "transpose matrix nonzeros = %" PetscInt_FMT ", allocated nonzeros = %" PetscInt_FMT "\n", (PetscInt)info.nz_used, (PetscInt)info.nz_allocated);
 71:   MatNorm(tmat, NORM_FROBENIUS, &normf);
 72:   MatNorm(tmat, NORM_1, &norm1);
 73:   MatNorm(tmat, NORM_INFINITY, &normi);
 74:   PetscPrintf(PETSC_COMM_WORLD, "transpose: Frobenious norm = %g, one norm = %g, infinity norm = %g\n", (double)normf, (double)norm1, (double)normi);
 75:   MatView(tmat, PETSC_VIEWER_STDOUT_WORLD);

 77:   /* Test MatAXPY */
 78:   if (mat && !rect) {
 79:     PetscScalar alpha = 1.0;
 80:     PetscOptionsGetScalar(NULL, NULL, "-alpha", &alpha, NULL);
 81:     PetscPrintf(PETSC_COMM_WORLD, "matrix addition:  B = B + alpha * A\n");
 82:     MatAXPY(tmat, alpha, mat, DIFFERENT_NONZERO_PATTERN);
 83:     MatView(tmat, PETSC_VIEWER_STDOUT_WORLD);
 84:   }

 86:   /* Free data structures */
 87:   MatDestroy(&tmat);
 88:   if (mat) MatDestroy(&mat);

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

 94: /*TEST

 96:    test:

 98:    testset:
 99:      args: -rect1
100:      test:
101:        suffix: r1
102:        output_file: output/ex49_r1.out
103:      test:
104:        suffix: r1_inplace
105:        args: -in_place
106:        output_file: output/ex49_r1.out
107:      test:
108:        suffix: r1_par
109:        nsize: 2
110:        output_file: output/ex49_r1_par.out
111:      test:
112:        suffix: r1_par_inplace
113:        args: -in_place
114:        nsize: 2
115:        output_file: output/ex49_r1_par.out

117: TEST*/