Actual source code: ex212.c

  1: static char help[] = "Test MatCreateSubMatrix with -mat_type nest and block sizes.\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **argv)
  6: {
  7:   Mat                    A, B, C, mats[6];
  8:   IS                     rows[2];
  9:   ISLocalToGlobalMapping cmap, rmap;
 10:   const PetscInt         indices[3] = {0, 1, 2};
 11:   PetscInt               i;
 12:   PetscMPIInt            size;

 15:   PetscInitialize(&argc, &argv, NULL, help);
 16:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 18:   MatCreate(PETSC_COMM_WORLD, &A);
 19:   MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, 2, 1);
 20:   MatSetBlockSizes(A, 2, 1);
 21:   MatSetType(A, MATAIJ);
 22:   ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 2, 1, indices, PETSC_COPY_VALUES, &rmap);
 23:   ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, 1, indices, PETSC_COPY_VALUES, &cmap);
 24:   MatSetLocalToGlobalMapping(A, rmap, cmap);
 25:   ISLocalToGlobalMappingDestroy(&rmap);
 26:   ISLocalToGlobalMappingDestroy(&cmap);
 27:   ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, 2, indices, PETSC_COPY_VALUES, &rmap);
 28:   ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, 3, indices, PETSC_COPY_VALUES, &cmap);
 29:   MatCreate(PETSC_COMM_WORLD, &B);
 30:   MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, 2, 3);
 31:   MatSetBlockSizes(B, 1, 1);
 32:   MatSetType(B, MATAIJ);
 33:   MatSetLocalToGlobalMapping(B, rmap, cmap);
 34:   ISLocalToGlobalMappingDestroy(&rmap);
 35:   ISLocalToGlobalMappingDestroy(&cmap);
 36:   MatSetUp(A);
 37:   MatSetUp(B);
 38:   mats[0] = A;
 39:   mats[1] = B;
 40:   mats[2] = A;
 41:   mats[3] = NULL;
 42:   mats[4] = B;
 43:   mats[5] = A;
 44:   MatCreateNest(PETSC_COMM_WORLD, 2, NULL, 3, NULL, mats, &C);
 45:   MatSetUp(C);
 46:   MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY);
 47:   MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY);
 48:   MatView(C, NULL);
 49:   MatNestGetISs(C, rows, NULL);
 50:   for (i = 0; i < 2; i++) {
 51:     Mat      submat;
 52:     IS       cols[3];
 53:     PetscInt j;
 54:     ISView(rows[i], NULL);
 55:     MatCreateSubMatrix(C, rows[i], NULL, MAT_INITIAL_MATRIX, &submat);
 56:     MatView(submat, NULL);
 57:     MatNestGetISs(submat, NULL, cols);
 58:     for (j = 0; j < 3; j++) ISView(cols[j], NULL);
 59:     MatDestroy(&submat);
 60:   }
 61:   MatDestroy(&C);
 62:   MatDestroy(&B);
 63:   MatDestroy(&A);
 64:   PetscFinalize();
 65:   return 0;
 66: }

 68: /*TEST

 70:     test:

 72: TEST*/