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