1*c4762a1bSJed Brown 2*c4762a1bSJed Brown static char help[] = "Tests late MatSetBlockSizes.\n\n"; 3*c4762a1bSJed Brown 4*c4762a1bSJed Brown #include <petscmat.h> 5*c4762a1bSJed Brown 6*c4762a1bSJed Brown int main(int argc,char **args) 7*c4762a1bSJed Brown { 8*c4762a1bSJed Brown Mat A; 9*c4762a1bSJed Brown Vec x[4]; 10*c4762a1bSJed Brown IS is; 11*c4762a1bSJed Brown ISLocalToGlobalMapping rmap,cmap; 12*c4762a1bSJed Brown PetscInt bs[4],l2gbs[4],rbs,cbs,l2grbs,l2gcbs,i; 13*c4762a1bSJed Brown PetscErrorCode ierr; 14*c4762a1bSJed Brown 15*c4762a1bSJed Brown ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 16*c4762a1bSJed Brown ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 17*c4762a1bSJed Brown ierr = MatSetSizes(A,12,12,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr); 18*c4762a1bSJed Brown ierr = MatSetType(A,MATAIJ);CHKERRQ(ierr); 19*c4762a1bSJed Brown ierr = ISCreateStride(PETSC_COMM_WORLD,12,0,1,&is);CHKERRQ(ierr); 20*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingCreateIS(is,&rmap);CHKERRQ(ierr); 21*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingSetBlockSize(rmap,2);CHKERRQ(ierr); 22*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingCreateIS(is,&cmap);CHKERRQ(ierr); 23*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingSetBlockSize(cmap,2);CHKERRQ(ierr); 24*c4762a1bSJed Brown 25*c4762a1bSJed Brown ierr = MatSetLocalToGlobalMapping(A,rmap,cmap);CHKERRQ(ierr); 26*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingDestroy(&rmap);CHKERRQ(ierr); 27*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingDestroy(&cmap);CHKERRQ(ierr); 28*c4762a1bSJed Brown ierr = ISDestroy(&is);CHKERRQ(ierr); 29*c4762a1bSJed Brown ierr = MatSetUp(A);CHKERRQ(ierr); 30*c4762a1bSJed Brown 31*c4762a1bSJed Brown ierr = MatCreateVecs(A,&x[1],&x[0]);CHKERRQ(ierr); 32*c4762a1bSJed Brown ierr = MatSetBlockSizes(A,6,3);CHKERRQ(ierr); 33*c4762a1bSJed Brown ierr = MatCreateVecs(A,&x[3],&x[2]);CHKERRQ(ierr); 34*c4762a1bSJed Brown for (i=0;i<4;i++) { 35*c4762a1bSJed Brown ISLocalToGlobalMapping l2g; 36*c4762a1bSJed Brown 37*c4762a1bSJed Brown ierr = VecGetBlockSize(x[i],&bs[i]);CHKERRQ(ierr); 38*c4762a1bSJed Brown ierr = VecGetLocalToGlobalMapping(x[i],&l2g);CHKERRQ(ierr); 39*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingGetBlockSize(l2g,&l2gbs[i]);CHKERRQ(ierr); 40*c4762a1bSJed Brown ierr = VecDestroy(&x[i]);CHKERRQ(ierr); 41*c4762a1bSJed Brown } 42*c4762a1bSJed Brown ierr = MatGetBlockSizes(A,&rbs,&cbs);CHKERRQ(ierr); 43*c4762a1bSJed Brown ierr = MatGetLocalToGlobalMapping(A,&rmap,&cmap);CHKERRQ(ierr); 44*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingGetBlockSize(rmap,&l2grbs);CHKERRQ(ierr); 45*c4762a1bSJed Brown ierr = ISLocalToGlobalMappingGetBlockSize(cmap,&l2gcbs);CHKERRQ(ierr); 46*c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 47*c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_WORLD,"Mat Block sizes: %D %D (l2g %D %D)\n",rbs,cbs,l2grbs,l2gcbs);CHKERRQ(ierr); 48*c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_WORLD,"Vec Block sizes: %D %D (l2g %D %D)\n",bs[0],bs[1],l2gbs[0],l2gbs[1]);CHKERRQ(ierr); 49*c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_WORLD,"Vec Block sizes: %D %D (l2g %D %D)\n",bs[2],bs[3],l2gbs[2],l2gbs[3]);CHKERRQ(ierr); 50*c4762a1bSJed Brown ierr = PetscFinalize(); 51*c4762a1bSJed Brown return ierr; 52*c4762a1bSJed Brown } 53*c4762a1bSJed Brown 54*c4762a1bSJed Brown 55*c4762a1bSJed Brown 56*c4762a1bSJed Brown /*TEST 57*c4762a1bSJed Brown 58*c4762a1bSJed Brown test: 59*c4762a1bSJed Brown nsize: 2 60*c4762a1bSJed Brown 61*c4762a1bSJed Brown TEST*/ 62