1 static char help[] = "Tests MatMPI{AIJ,BAIJ,SBAIJ}SetPreallocationCSR\n\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc,char **argv) 6 { 7 PetscErrorCode ierr; 8 PetscInt ia[3]={0,2,4}; 9 PetscInt ja[4]={0,1,0,1}; 10 PetscScalar c[16]={0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15}; 11 PetscInt ia2[5]={0,4,8,12,16}; 12 PetscInt ja2[16]={0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3}; 13 PetscScalar c2[16]={0,1,4,5,2,3,6,7,8,9,12,13,10,11,14,15}; 14 PetscMPIInt size,rank; 15 Mat ssbaij; 16 PetscBool rect = PETSC_FALSE; 17 18 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 19 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 20 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 21 PetscCheckFalse(size < 2,PETSC_COMM_WORLD,PETSC_ERR_SUP,"This is an example with more then one processors"); 22 if (rank) { 23 PetscInt i; 24 for (i = 0; i < 3; i++) ia[i] = 0; 25 for (i = 0; i < 5; i++) ia2[i] = 0; 26 } 27 CHKERRQ(PetscOptionsGetBool(NULL,NULL,"-rect",&rect,NULL)); 28 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&ssbaij)); 29 CHKERRQ(MatSetBlockSize(ssbaij,2)); 30 if (rect) { 31 CHKERRQ(MatSetType(ssbaij,MATMPIBAIJ)); 32 CHKERRQ(MatSetSizes(ssbaij,4,6,PETSC_DECIDE,PETSC_DECIDE)); 33 } else { 34 CHKERRQ(MatSetType(ssbaij,MATMPISBAIJ)); 35 CHKERRQ(MatSetSizes(ssbaij,4,4,PETSC_DECIDE,PETSC_DECIDE)); 36 } 37 CHKERRQ(MatSetFromOptions(ssbaij)); 38 CHKERRQ(MatMPIAIJSetPreallocationCSR(ssbaij,ia2,ja2,c2)); 39 CHKERRQ(MatMPIBAIJSetPreallocationCSR(ssbaij,2,ia,ja,c)); 40 CHKERRQ(MatMPISBAIJSetPreallocationCSR(ssbaij,2,ia,ja,c)); 41 CHKERRQ(MatViewFromOptions(ssbaij,NULL,"-view")); 42 CHKERRQ(MatDestroy(&ssbaij)); 43 ierr = PetscFinalize(); 44 return ierr; 45 } 46 47 /*TEST 48 49 test: 50 filter: grep -v type | sed -e "s/\.//g" 51 suffix: aijbaij_csr 52 nsize: 2 53 args: -mat_type {{aij baij}} -view -rect {{0 1}} 54 55 test: 56 filter: sed -e "s/\.//g" 57 suffix: sbaij_csr 58 nsize: 2 59 args: -mat_type sbaij -view -rect {{0 1}} 60 61 TEST*/ 62