1c4762a1bSJed Brown static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown 5c4762a1bSJed Brown int main(int argc,char **argv) 6c4762a1bSJed Brown { 7c4762a1bSJed Brown Mat A,B; 8c4762a1bSJed Brown PetscInt i,j,column,*ooj; 9c4762a1bSJed Brown PetscInt *di,*dj,*oi,*oj,nd; 10c4762a1bSJed Brown const PetscInt *garray; 11c4762a1bSJed Brown PetscScalar *oa,*da; 12c4762a1bSJed Brown PetscScalar value; 13c4762a1bSJed Brown PetscRandom rctx; 14c4762a1bSJed Brown PetscErrorCode ierr; 15c4762a1bSJed Brown PetscBool equal,done; 16c4762a1bSJed Brown Mat AA,AB; 17c4762a1bSJed Brown PetscMPIInt size,rank; 18c4762a1bSJed Brown 19c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 20ffc4695bSBarry Smith ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRMPI(ierr); 21*2c71b3e2SJacob Faibussowitsch PetscCheckFalse(size == 1,PETSC_COMM_WORLD,PETSC_ERR_ARG_OUTOFRANGE,"Must run with 2 or more processes"); 22ffc4695bSBarry Smith ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr); 23c4762a1bSJed Brown 24c4762a1bSJed Brown /* Create a mpiaij matrix for checking */ 25c4762a1bSJed Brown ierr = MatCreateAIJ(PETSC_COMM_WORLD,5,5,PETSC_DECIDE,PETSC_DECIDE,0,NULL,0,NULL,&A);CHKERRQ(ierr); 26c4762a1bSJed Brown ierr = MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);CHKERRQ(ierr); 27c4762a1bSJed Brown ierr = MatSetUp(A);CHKERRQ(ierr); 28c4762a1bSJed Brown ierr = PetscRandomCreate(PETSC_COMM_WORLD,&rctx);CHKERRQ(ierr); 29c4762a1bSJed Brown ierr = PetscRandomSetFromOptions(rctx);CHKERRQ(ierr); 30c4762a1bSJed Brown 31c4762a1bSJed Brown for (i=5*rank; i<5*rank+5; i++) { 32c4762a1bSJed Brown for (j=0; j<5*size; j++) { 33c4762a1bSJed Brown ierr = PetscRandomGetValue(rctx,&value);CHKERRQ(ierr); 34c4762a1bSJed Brown column = (PetscInt) (5*size*PetscRealPart(value)); 35c4762a1bSJed Brown ierr = PetscRandomGetValue(rctx,&value);CHKERRQ(ierr); 36c4762a1bSJed Brown ierr = MatSetValues(A,1,&i,1,&column,&value,INSERT_VALUES);CHKERRQ(ierr); 37c4762a1bSJed Brown } 38c4762a1bSJed Brown } 39c4762a1bSJed Brown ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 40c4762a1bSJed Brown ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 41c4762a1bSJed Brown 42c4762a1bSJed Brown ierr = MatMPIAIJGetSeqAIJ(A,&AA,&AB,&garray);CHKERRQ(ierr); 43c4762a1bSJed Brown ierr = MatGetRowIJ(AA,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&di,(const PetscInt**)&dj,&done);CHKERRQ(ierr); 44c4762a1bSJed Brown ierr = MatSeqAIJGetArray(AA,&da);CHKERRQ(ierr); 45c4762a1bSJed Brown ierr = MatGetRowIJ(AB,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&oi,(const PetscInt**)&oj,&done);CHKERRQ(ierr); 46c4762a1bSJed Brown ierr = MatSeqAIJGetArray(AB,&oa);CHKERRQ(ierr); 47c4762a1bSJed Brown 48c4762a1bSJed Brown ierr = PetscMalloc1(oi[5],&ooj);CHKERRQ(ierr); 49c4762a1bSJed Brown ierr = PetscArraycpy(ooj,oj,oi[5]);CHKERRQ(ierr); 50c4762a1bSJed Brown /* modify the column entries in the non-diagonal portion back to global numbering */ 51c4762a1bSJed Brown for (i=0; i<oi[5]; i++) { 52c4762a1bSJed Brown ooj[i] = garray[ooj[i]]; 53c4762a1bSJed Brown } 54c4762a1bSJed Brown 55c4762a1bSJed Brown ierr = MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD,5,5,PETSC_DETERMINE,PETSC_DETERMINE,di,dj,da,oi,ooj,oa,&B);CHKERRQ(ierr); 56c4762a1bSJed Brown ierr = MatSetUp(B);CHKERRQ(ierr); 57c4762a1bSJed Brown ierr = MatEqual(A,B,&equal);CHKERRQ(ierr); 58c4762a1bSJed Brown 59c4762a1bSJed Brown ierr = MatRestoreRowIJ(AA,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&di,(const PetscInt**)&dj,&done);CHKERRQ(ierr); 60c4762a1bSJed Brown ierr = MatSeqAIJRestoreArray(AA,&da);CHKERRQ(ierr); 61c4762a1bSJed Brown ierr = MatRestoreRowIJ(AB,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&oi,(const PetscInt**)&oj,&done);CHKERRQ(ierr); 62c4762a1bSJed Brown ierr = MatSeqAIJRestoreArray(AB,&oa);CHKERRQ(ierr); 63c4762a1bSJed Brown 64*2c71b3e2SJacob Faibussowitsch PetscCheckFalse(!equal,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Likely a bug in MatCreateMPIAIJWithSplitArrays()"); 65c4762a1bSJed Brown 66c4762a1bSJed Brown /* Free spaces */ 67c4762a1bSJed Brown ierr = PetscFree(ooj);CHKERRQ(ierr); 68c4762a1bSJed Brown ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr); 69c4762a1bSJed Brown ierr = MatDestroy(&A);CHKERRQ(ierr); 70c4762a1bSJed Brown ierr = MatDestroy(&B);CHKERRQ(ierr); 71c4762a1bSJed Brown ierr = PetscFree(oj);CHKERRQ(ierr); 72c4762a1bSJed Brown ierr = PetscFinalize(); 73c4762a1bSJed Brown return ierr; 74c4762a1bSJed Brown } 75c4762a1bSJed Brown 76c4762a1bSJed Brown /*TEST 77c4762a1bSJed Brown 78c4762a1bSJed Brown test: 79c4762a1bSJed Brown nsize: 3 80c4762a1bSJed Brown 81c4762a1bSJed Brown TEST*/ 82