Actual source code: ex58.c


  2: static char help[] = "Tests MatTranspose() and MatEqual() for MPIAIJ matrices.\n\n";

  4: #include <petscmat.h>

  6: int main(int argc, char **argv)
  7: {
  8:   Mat         A, B;
  9:   PetscInt    m = 7, n, i, rstart, rend, cols[3];
 10:   PetscScalar v[3];
 11:   PetscBool   equal;
 12:   const char *eq[2];

 15:   PetscInitialize(&argc, &argv, (char *)0, help);
 16:   PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);
 17:   n = m;

 19:   /* ------- Assemble matrix, --------- */

 21:   MatCreateAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, m, n, 0, 0, 0, 0, &A);
 22:   MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE);
 23:   MatGetOwnershipRange(A, &rstart, &rend);
 24:   if (!rstart) {
 25:     cols[0] = 0;
 26:     cols[1] = 1;
 27:     v[0]    = 2.0;
 28:     v[1]    = -1.0;
 29:     MatSetValues(A, 1, &rstart, 2, cols, v, INSERT_VALUES);
 30:     rstart++;
 31:   }
 32:   if (rend == m) {
 33:     rend--;
 34:     cols[0] = rend - 1;
 35:     cols[1] = rend;
 36:     v[0]    = -1.0;
 37:     v[1]    = 2.0;
 38:     MatSetValues(A, 1, &rend, 2, cols, v, INSERT_VALUES);
 39:   }
 40:   v[0] = -1.0;
 41:   v[1] = 2.0;
 42:   v[2] = -1.0;
 43:   for (i = rstart; i < rend; i++) {
 44:     cols[0] = i - 1;
 45:     cols[1] = i;
 46:     cols[2] = i + 1;
 47:     MatSetValues(A, 1, &i, 3, cols, v, INSERT_VALUES);
 48:   }
 49:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 50:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);

 52:   MatTranspose(A, MAT_INITIAL_MATRIX, &B);

 54:   MatEqual(A, B, &equal);

 56:   eq[0] = "not equal";
 57:   eq[1] = "equal";
 58:   PetscPrintf(PETSC_COMM_WORLD, "Matrices are %s\n", eq[equal]);

 60:   MatTranspose(A, MAT_REUSE_MATRIX, &B);
 61:   MatEqual(A, B, &equal);
 62:   if (!equal) PetscPrintf(PETSC_COMM_WORLD, "MatTranspose with MAT_REUSE_MATRIX failed");

 64:   /* Free data structures */
 65:   MatDestroy(&A);
 66:   MatDestroy(&B);

 68:   PetscFinalize();
 69:   return 0;
 70: }

 72: /*TEST

 74:     test:

 76:     test:
 77:       suffix: 2
 78:       nsize: 2
 79:       output_file: output/ex58_1.out

 81: TEST*/