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