1 static char help[] = "Test MatCreateSubMatrix with -mat_type nest and block sizes.\n"; 2 3 #include <petscmat.h> 4 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; 13 PetscErrorCode ierr; 14 15 ierr = PetscInitialize(&argc,&argv,NULL,help); if (ierr) return ierr; 16 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 17 PetscCheckFalse(size > 1,PETSC_COMM_WORLD,PETSC_ERR_ARG_WRONG,"Only coded for one process"); 18 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&A)); 19 CHKERRQ(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,2,1)); 20 CHKERRQ(MatSetBlockSizes(A,2,1)); 21 CHKERRQ(MatSetType(A,MATAIJ)); 22 CHKERRQ(ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD,2,1,indices,PETSC_COPY_VALUES,&rmap)); 23 CHKERRQ(ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD,1,1,indices,PETSC_COPY_VALUES,&cmap)); 24 CHKERRQ(MatSetLocalToGlobalMapping(A,rmap,cmap)); 25 CHKERRQ(ISLocalToGlobalMappingDestroy(&rmap)); 26 CHKERRQ(ISLocalToGlobalMappingDestroy(&cmap)); 27 CHKERRQ(ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD,1,2,indices,PETSC_COPY_VALUES,&rmap)); 28 CHKERRQ(ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD,1,3,indices,PETSC_COPY_VALUES,&cmap)); 29 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&B)); 30 CHKERRQ(MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,2,3)); 31 CHKERRQ(MatSetBlockSizes(B,1,1)); 32 CHKERRQ(MatSetType(B,MATAIJ)); 33 CHKERRQ(MatSetLocalToGlobalMapping(B,rmap,cmap)); 34 CHKERRQ(ISLocalToGlobalMappingDestroy(&rmap)); 35 CHKERRQ(ISLocalToGlobalMappingDestroy(&cmap)); 36 CHKERRQ(MatSetUp(A)); 37 CHKERRQ(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 CHKERRQ(MatCreateNest(PETSC_COMM_WORLD,2,NULL,3,NULL,mats,&C)); 45 CHKERRQ(MatSetUp(C)); 46 CHKERRQ(MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY)); 47 CHKERRQ(MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY)); 48 CHKERRQ(MatView(C,NULL)); 49 CHKERRQ(MatNestGetISs(C,rows,NULL)); 50 for (i=0; i<2; i++) { 51 Mat submat; 52 IS cols[3]; 53 PetscInt j; 54 CHKERRQ(ISView(rows[i],NULL)); 55 CHKERRQ(MatCreateSubMatrix(C,rows[i],NULL,MAT_INITIAL_MATRIX,&submat)); 56 CHKERRQ(MatView(submat,NULL)); 57 CHKERRQ(MatNestGetISs(submat,NULL,cols)); 58 for (j=0; j<3; j++) { 59 CHKERRQ(ISView(cols[j],NULL)); 60 } 61 CHKERRQ(MatDestroy(&submat)); 62 } 63 CHKERRQ(MatDestroy(&C)); 64 CHKERRQ(MatDestroy(&B)); 65 CHKERRQ(MatDestroy(&A)); 66 ierr = PetscFinalize(); 67 return ierr; 68 } 69 70 /*TEST 71 72 test: 73 74 TEST*/ 75