Actual source code: ex46.c


  2: static char help[] = "Tests late MatSetBlockSizes.\n\n";

  4: #include <petscmat.h>

  6: int main(int argc,char **args)
  7: {
  8:   Mat                    A;
  9:   Vec                    x[4];
 10:   IS                     is;
 11:   ISLocalToGlobalMapping rmap,cmap;
 12:   PetscInt               bs[4],l2gbs[4],rbs,cbs,l2grbs,l2gcbs,i;

 14:   PetscInitialize(&argc,&args,(char*)0,help);
 15:   MatCreate(PETSC_COMM_WORLD,&A);
 16:   MatSetSizes(A,12,12,PETSC_DECIDE,PETSC_DECIDE);
 17:   MatSetType(A,MATAIJ);
 18:   ISCreateStride(PETSC_COMM_WORLD,12,0,1,&is);
 19:   ISLocalToGlobalMappingCreateIS(is,&rmap);
 20:   ISLocalToGlobalMappingSetBlockSize(rmap,2);
 21:   ISLocalToGlobalMappingCreateIS(is,&cmap);
 22:   ISLocalToGlobalMappingSetBlockSize(cmap,2);

 24:   MatSetLocalToGlobalMapping(A,rmap,cmap);
 25:   ISLocalToGlobalMappingDestroy(&rmap);
 26:   ISLocalToGlobalMappingDestroy(&cmap);
 27:   ISDestroy(&is);
 28:   MatSetUp(A);

 30:   MatCreateVecs(A,&x[1],&x[0]);
 31:   MatSetBlockSizes(A,6,3);
 32:   MatCreateVecs(A,&x[3],&x[2]);
 33:   for (i=0;i<4;i++) {
 34:     ISLocalToGlobalMapping l2g;

 36:     VecGetBlockSize(x[i],&bs[i]);
 37:     VecGetLocalToGlobalMapping(x[i],&l2g);
 38:     ISLocalToGlobalMappingGetBlockSize(l2g,&l2gbs[i]);
 39:     VecDestroy(&x[i]);
 40:   }
 41:   MatGetBlockSizes(A,&rbs,&cbs);
 42:   MatGetLocalToGlobalMapping(A,&rmap,&cmap);
 43:   ISLocalToGlobalMappingGetBlockSize(rmap,&l2grbs);
 44:   ISLocalToGlobalMappingGetBlockSize(cmap,&l2gcbs);
 45:   MatDestroy(&A);
 46:   PetscPrintf(PETSC_COMM_WORLD,"Mat Block sizes: %" PetscInt_FMT " %" PetscInt_FMT " (l2g %" PetscInt_FMT " %" PetscInt_FMT ")\n",rbs,cbs,l2grbs,l2gcbs);
 47:   PetscPrintf(PETSC_COMM_WORLD,"Vec Block sizes: %" PetscInt_FMT " %" PetscInt_FMT " (l2g %" PetscInt_FMT " %" PetscInt_FMT ")\n",bs[0],bs[1],l2gbs[0],l2gbs[1]);
 48:   PetscPrintf(PETSC_COMM_WORLD,"Vec Block sizes: %" PetscInt_FMT " %" PetscInt_FMT " (l2g %" PetscInt_FMT " %" PetscInt_FMT ")\n",bs[2],bs[3],l2gbs[2],l2gbs[3]);
 49:   PetscFinalize();
 50:   return 0;
 51: }

 53: /*TEST

 55:    test:
 56:       nsize: 2

 58: TEST*/