Actual source code: ex177.c


  2: static char help[] = "Tests various routines in MatKAIJ format.\n";

  4: #include <petscmat.h>
  5: #define IMAX 15

  7: int main(int argc, char **args)
  8: {
  9:   Mat          A, B, TA;
 10:   PetscScalar *S, *T;
 11:   PetscViewer  fd;
 12:   char         file[PETSC_MAX_PATH_LEN];
 13:   PetscInt     m, n, M, N, p = 1, q = 1, i, j;
 14:   PetscMPIInt  rank, size;
 15:   PetscBool    flg;

 18:   PetscInitialize(&argc, &args, (char *)0, help);
 19:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 20:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 22:   /* Load AIJ matrix A */
 23:   PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg);
 25:   PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd);
 26:   MatCreate(PETSC_COMM_WORLD, &A);
 27:   MatLoad(A, fd);
 28:   PetscViewerDestroy(&fd);

 30:   /* Get dof, then create S and T */
 31:   PetscOptionsGetInt(NULL, NULL, "-p", &p, NULL);
 32:   PetscOptionsGetInt(NULL, NULL, "-q", &q, NULL);
 33:   PetscMalloc2(p * q, &S, p * q, &T);
 34:   for (i = 0; i < p * q; i++) S[i] = 0;

 36:   for (i = 0; i < p; i++) {
 37:     for (j = 0; j < q; j++) {
 38:       /* Set some random non-zero values */
 39:       S[i + p * j] = ((PetscReal)((i + 1) * (j + 1))) / ((PetscReal)(p + q));
 40:       T[i + p * j] = ((PetscReal)((p - i) + j)) / ((PetscReal)(p * q));
 41:     }
 42:   }

 44:   /* Test KAIJ when both S & T are not NULL */

 46:   /* Create KAIJ matrix TA */
 47:   MatCreateKAIJ(A, p, q, S, T, &TA);
 48:   MatGetLocalSize(A, &m, &n);
 49:   MatGetSize(A, &M, &N);

 51:   MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B);

 53:   /* Test MatKAIJGetScaledIdentity() */
 54:   MatKAIJGetScaledIdentity(TA, &flg);
 56:   /* Test MatMult() */
 57:   MatMultEqual(TA, B, 10, &flg);
 59:   /* Test MatMultAdd() */
 60:   MatMultAddEqual(TA, B, 10, &flg);

 63:   MatDestroy(&TA);
 64:   MatDestroy(&B);

 66:   /* Test KAIJ when S is NULL */

 68:   /* Create KAIJ matrix TA */
 69:   MatCreateKAIJ(A, p, q, NULL, T, &TA);
 70:   MatGetLocalSize(A, &m, &n);
 71:   MatGetSize(A, &M, &N);

 73:   MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B);

 75:   /* Test MatKAIJGetScaledIdentity() */
 76:   MatKAIJGetScaledIdentity(TA, &flg);
 78:   /* Test MatMult() */
 79:   MatMultEqual(TA, B, 10, &flg);
 81:   /* Test MatMultAdd() */
 82:   MatMultAddEqual(TA, B, 10, &flg);

 85:   MatDestroy(&TA);
 86:   MatDestroy(&B);

 88:   /* Test KAIJ when T is NULL */

 90:   /* Create KAIJ matrix TA */
 91:   MatCreateKAIJ(A, p, q, S, NULL, &TA);
 92:   MatGetLocalSize(A, &m, &n);
 93:   MatGetSize(A, &M, &N);

 95:   MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B);

 97:   /* Test MatKAIJGetScaledIdentity() */
 98:   MatKAIJGetScaledIdentity(TA, &flg);
100:   /* Test MatMult() */
101:   MatMultEqual(TA, B, 10, &flg);
103:   /* Test MatMultAdd() */
104:   MatMultAddEqual(TA, B, 10, &flg);

107:   MatDestroy(&TA);
108:   MatDestroy(&B);

110:   /* Test KAIJ when T is is an identity matrix */

112:   if (p == q) {
113:     for (i = 0; i < p; i++) {
114:       for (j = 0; j < q; j++) {
115:         if (i == j) T[i + j * p] = 1.0;
116:         else T[i + j * p] = 0.0;
117:       }
118:     }

120:     /* Create KAIJ matrix TA */
121:     MatCreateKAIJ(A, p, q, S, T, &TA);
122:     MatGetLocalSize(A, &m, &n);
123:     MatGetSize(A, &M, &N);

125:     MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B);

127:     /* Test MatKAIJGetScaledIdentity() */
128:     MatKAIJGetScaledIdentity(TA, &flg);
130:     /* Test MatMult() */
131:     MatMultEqual(TA, B, 10, &flg);
133:     /* Test MatMultAdd() */
134:     MatMultAddEqual(TA, B, 10, &flg);

137:     MatDestroy(&TA);
138:     MatDestroy(&B);

140:     MatCreateKAIJ(A, p, q, NULL, T, &TA);
141:     /* Test MatKAIJGetScaledIdentity() */
142:     MatKAIJGetScaledIdentity(TA, &flg);
144:     MatDestroy(&TA);

146:     for (i = 0; i < p; i++) {
147:       for (j = 0; j < q; j++) {
148:         if (i == j) S[i + j * p] = T[i + j * p] = 2.0;
149:         else S[i + j * p] = T[i + j * p] = 0.0;
150:       }
151:     }
152:     MatCreateKAIJ(A, p, q, S, T, &TA);
153:     /* Test MatKAIJGetScaledIdentity() */
154:     MatKAIJGetScaledIdentity(TA, &flg);
156:     MatDestroy(&TA);
157:   }

159:   /* Done with all tests */

161:   PetscFree2(S, T);
162:   MatDestroy(&A);
163:   PetscFinalize();
164:   return 0;
165: }

167: /*TEST

169:   build:
170:     requires: !complex

172:   test:
173:     requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
174:     output_file: output/ex177.out
175:     nsize: {{1 2 3 4}}
176:     args: -f ${DATAFILESPATH}/matrices/small -p {{2 3 7}} -q {{3 7}} -viewer_binary_skip_info

178: TEST*/