1 static const char help[] = "Test parallel assembly of SBAIJ matrices\n\n"; 2 3 #include <petscmat.h> 4 5 PetscErrorCode Assemble(MPI_Comm comm,PetscInt n,MatType mtype) 6 { 7 Mat A; 8 PetscInt first,last,i; 9 PetscMPIInt rank,size; 10 11 PetscFunctionBegin; 12 PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 13 PetscCall(MatSetSizes(A, PETSC_DECIDE,PETSC_DECIDE,n,n)); 14 PetscCall(MatSetType(A,MATMPISBAIJ)); 15 PetscCall(MatSetFromOptions(A)); 16 PetscCallMPI(MPI_Comm_size(comm,&size)); 17 PetscCallMPI(MPI_Comm_rank(comm,&rank)); 18 if (rank < size-1) { 19 PetscCall(MatMPISBAIJSetPreallocation(A,1,1,NULL,1,NULL)); 20 } else { 21 PetscCall(MatMPISBAIJSetPreallocation(A,1,2,NULL,0,NULL)); 22 } 23 PetscCall(MatGetOwnershipRange(A,&first,&last)); 24 PetscCall(MatSetOption(A,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE)); 25 last--; 26 for (i=first; i<=last; i++) { 27 PetscCall(MatSetValue(A,i,i,2.,INSERT_VALUES)); 28 if (i != n-1) PetscCall(MatSetValue(A,i,n-1,-1.,INSERT_VALUES)); 29 } 30 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 31 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 32 PetscCall(MatDestroy(&A)); 33 PetscFunctionReturn(0); 34 } 35 36 int main(int argc,char *argv[]) 37 { 38 MPI_Comm comm; 39 PetscInt n = 6; 40 41 PetscCall(PetscInitialize(&argc,&argv,NULL,help)); 42 comm = PETSC_COMM_WORLD; 43 PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 44 PetscCall(Assemble(comm,n,MATMPISBAIJ)); 45 PetscCall(PetscFinalize()); 46 return 0; 47 } 48 49 /*TEST 50 51 test: 52 nsize: 4 53 args: -n 1000 -mat_view ascii::ascii_info_detail 54 requires: double !complex !defined(PETSC_USE_64BIT_INDICES) 55 56 TEST*/ 57