1c4762a1bSJed Brown 2c4762a1bSJed Brown /* tests MatSeqSBAIJSetPreallocationCSR() and MatMPISBAIJSetPreallocationCSR() */ 3c4762a1bSJed Brown 4c4762a1bSJed Brown #include <petsc.h> 5c4762a1bSJed Brown 6c4762a1bSJed Brown int main(int argc, char **args) 7c4762a1bSJed Brown { 8c4762a1bSJed Brown PetscInt ia[3] = { 0, 2, 4}; 9c4762a1bSJed Brown PetscInt ja[4] = { 0, 1, 0, 1}; 10c4762a1bSJed Brown PetscScalar c[16] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15 }; 11c4762a1bSJed Brown Mat ssbaij; 12c4762a1bSJed Brown Mat msbaij; 13c4762a1bSJed Brown 14*9566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&args,(char*)0,(char*)0)); 15c4762a1bSJed Brown 16*9566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_SELF, &ssbaij)); 17*9566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_SELF, &msbaij)); 18*9566063dSJacob Faibussowitsch PetscCall(MatSetType(ssbaij, MATSEQSBAIJ)); 19*9566063dSJacob Faibussowitsch PetscCall(MatSetType(msbaij, MATMPISBAIJ)); 20*9566063dSJacob Faibussowitsch PetscCall(MatSetBlockSize(ssbaij, 2)); 21*9566063dSJacob Faibussowitsch PetscCall(MatSetSizes(ssbaij, 4, 4, 4, 4)); 22*9566063dSJacob Faibussowitsch PetscCall(MatSetBlockSize(msbaij, 2)); 23*9566063dSJacob Faibussowitsch PetscCall(MatSetSizes(msbaij, 4, 4, 4, 4)); 24*9566063dSJacob Faibussowitsch PetscCall(MatSeqSBAIJSetPreallocationCSR(ssbaij, 2, ia, ja, c)); 25*9566063dSJacob Faibussowitsch PetscCall(MatMPISBAIJSetPreallocationCSR(msbaij, 2, ia, ja, c)); 26*9566063dSJacob Faibussowitsch PetscCall(MatView(ssbaij, PETSC_VIEWER_STDOUT_(PETSC_COMM_SELF))); 27*9566063dSJacob Faibussowitsch PetscCall(MatView(msbaij, PETSC_VIEWER_STDOUT_(PETSC_COMM_SELF))); 28*9566063dSJacob Faibussowitsch PetscCall(MatDestroy(&ssbaij)); 29*9566063dSJacob Faibussowitsch PetscCall(MatDestroy(&msbaij)); 30*9566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 31b122ec5aSJacob Faibussowitsch return 0; 32c4762a1bSJed Brown } 33c4762a1bSJed Brown 34c4762a1bSJed Brown /*TEST 35c4762a1bSJed Brown 36c4762a1bSJed Brown test: 37c4762a1bSJed Brown filter: sed "s?\.??g" 38c4762a1bSJed Brown 39c4762a1bSJed Brown TEST*/ 40