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