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