Actual source code: ex26.c


  2: static char help[] = "Tests MatGetRowIJ for SeqAIJ, SeqBAIJ and SeqSBAIJ\n\n";

  4: #include <petscmat.h>

  6: PetscErrorCode DumpCSR(Mat A, PetscInt shift, PetscBool symmetric, PetscBool compressed)
  7: {
  8:   MatType         type;
  9:   PetscInt        i, j, nr, bs = 1;
 10:   const PetscInt *ia, *ja;
 11:   PetscBool       done, isseqbaij, isseqsbaij;

 14:   PetscObjectBaseTypeCompare((PetscObject)A, MATSEQBAIJ, &isseqbaij);
 15:   PetscObjectBaseTypeCompare((PetscObject)A, MATSEQSBAIJ, &isseqsbaij);
 16:   if (isseqbaij || isseqsbaij) MatGetBlockSize(A, &bs);
 17:   MatGetType(A, &type);
 18:   MatGetRowIJ(A, shift, symmetric, compressed, &nr, &ia, &ja, &done);
 19:   PetscPrintf(PETSC_COMM_SELF, "===========================================================\n");
 20:   PetscPrintf(PETSC_COMM_SELF, "CSR for %s: shift %" PetscInt_FMT " symmetric %" PetscInt_FMT " compressed %" PetscInt_FMT "\n", type, shift, (PetscInt)symmetric, (PetscInt)compressed);
 21:   for (i = 0; i < nr; i++) {
 22:     PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ":", i + shift);
 23:     for (j = ia[i]; j < ia[i + 1]; j++) PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT, ja[j - shift]);
 24:     PetscPrintf(PETSC_COMM_SELF, "\n");
 25:   }
 26:   MatRestoreRowIJ(A, shift, symmetric, compressed, &nr, &ia, &ja, &done);
 27:   return 0;
 28: }

 30: int main(int argc, char **args)
 31: {
 32:   Mat         A, B, C;
 33:   PetscInt    i, j, k, m = 3, n = 3, bs = 1;
 34:   PetscMPIInt size;

 37:   PetscInitialize(&argc, &args, (char *)0, help);
 38:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 40:   PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);
 41:   PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
 42:   PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL);
 43:   /* adjust sizes by block size */
 44:   if (m % bs) m += bs - m % bs;
 45:   if (n % bs) n += bs - n % bs;

 47:   MatCreate(PETSC_COMM_SELF, &A);
 48:   MatSetSizes(A, m * n, m * n, PETSC_DECIDE, PETSC_DECIDE);
 49:   MatSetBlockSize(A, bs);
 50:   MatSetType(A, MATSEQAIJ);
 51:   MatSetUp(A);
 52:   MatCreate(PETSC_COMM_SELF, &B);
 53:   MatSetSizes(B, m * n, m * n, PETSC_DECIDE, PETSC_DECIDE);
 54:   MatSetBlockSize(B, bs);
 55:   MatSetType(B, MATSEQBAIJ);
 56:   MatSetUp(B);
 57:   MatCreate(PETSC_COMM_SELF, &C);
 58:   MatSetSizes(C, m * n, m * n, PETSC_DECIDE, PETSC_DECIDE);
 59:   MatSetBlockSize(C, bs);
 60:   MatSetType(C, MATSEQSBAIJ);
 61:   MatSetUp(C);
 62:   MatSetOption(C, MAT_IGNORE_LOWER_TRIANGULAR, PETSC_TRUE);

 64:   for (i = 0; i < m; i++) {
 65:     for (j = 0; j < n; j++) {
 66:       PetscScalar v  = -1.0;
 67:       PetscInt    Ii = j + n * i, J;
 68:       J              = Ii - n;
 69:       if (J >= 0) {
 70:         MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 71:         MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 72:         MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 73:       }
 74:       J = Ii + n;
 75:       if (J < m * n) {
 76:         MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 77:         MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 78:         MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 79:       }
 80:       J = Ii - 1;
 81:       if (J >= 0) {
 82:         MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 83:         MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 84:         MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 85:       }
 86:       J = Ii + 1;
 87:       if (J < m * n) {
 88:         MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 89:         MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 90:         MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES);
 91:       }
 92:       v = 4.0;
 93:       MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES);
 94:       MatSetValues(B, 1, &Ii, 1, &Ii, &v, INSERT_VALUES);
 95:       MatSetValues(C, 1, &Ii, 1, &Ii, &v, INSERT_VALUES);
 96:     }
 97:   }
 98:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 99:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
100:   MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
101:   MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
102:   MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY);
103:   MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY);

105:   /* test MatGetRowIJ for the three Mat types */
106:   MatView(A, NULL);
107:   MatView(B, NULL);
108:   MatView(C, NULL);
109:   for (i = 0; i < 2; i++) {
110:     PetscInt shift = i;
111:     for (j = 0; j < 2; j++) {
112:       PetscBool symmetric = ((j > 0) ? PETSC_FALSE : PETSC_TRUE);
113:       for (k = 0; k < 2; k++) {
114:         PetscBool compressed = ((k > 0) ? PETSC_FALSE : PETSC_TRUE);
115:         DumpCSR(A, shift, symmetric, compressed);
116:         DumpCSR(B, shift, symmetric, compressed);
117:         DumpCSR(C, shift, symmetric, compressed);
118:       }
119:     }
120:   }
121:   MatDestroy(&A);
122:   MatDestroy(&B);
123:   MatDestroy(&C);
124:   PetscFinalize();
125:   return 0;
126: }

128: /*TEST

130:    test:

132:    test:
133:       suffix: 2
134:       args: -bs 2

136: TEST*/