1: static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n";

  3: #include <petscmat.h>

  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;

 19:   PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
 20:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 21:   if (size == 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_ARG_OUTOFRANGE,"Must run with 2 or more processes");
 22:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);

 24:   /* Create a mpiaij matrix for checking */
 25:   MatCreateAIJ(PETSC_COMM_WORLD,5,5,PETSC_DECIDE,PETSC_DECIDE,0,NULL,0,NULL,&A);
 26:   MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
 27:   MatSetUp(A);
 28:   PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
 29:   PetscRandomSetFromOptions(rctx);

 31:   for (i=5*rank; i<5*rank+5; i++) {
 32:     for (j=0; j<5*size; j++) {
 33:       PetscRandomGetValue(rctx,&value);
 34:       column = (PetscInt) (5*size*PetscRealPart(value));
 35:       PetscRandomGetValue(rctx,&value);
 36:       MatSetValues(A,1,&i,1,&column,&value,INSERT_VALUES);
 37:     }
 38:   }
 39:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 40:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

 42:   MatMPIAIJGetSeqAIJ(A,&AA,&AB,&garray);
 43:   MatGetRowIJ(AA,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&di,(const PetscInt**)&dj,&done);
 44:   MatSeqAIJGetArray(AA,&da);
 45:   MatGetRowIJ(AB,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&oi,(const PetscInt**)&oj,&done);
 46:   MatSeqAIJGetArray(AB,&oa);

 48:   PetscMalloc1(oi[5],&ooj);
 49:   PetscArraycpy(ooj,oj,oi[5]);
 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:   }

 55:   MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD,5,5,PETSC_DETERMINE,PETSC_DETERMINE,di,dj,da,oi,ooj,oa,&B);
 56:   MatSetUp(B);
 57:   MatEqual(A,B,&equal);

 59:   MatRestoreRowIJ(AA,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&di,(const PetscInt**)&dj,&done);
 60:   MatSeqAIJRestoreArray(AA,&da);
 61:   MatRestoreRowIJ(AB,0,PETSC_FALSE,PETSC_FALSE,&nd,(const PetscInt**)&oi,(const PetscInt**)&oj,&done);
 62:   MatSeqAIJRestoreArray(AB,&oa);

 64:   if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Likely a bug in MatCreateMPIAIJWithSplitArrays()");

 66:   /* Free spaces */
 67:   PetscFree(ooj);
 68:   PetscRandomDestroy(&rctx);
 69:   MatDestroy(&A);
 70:   MatDestroy(&B);
 71:   PetscFree(oj);
 72:   PetscFinalize();
 73:   return ierr;
 74: }


 77: /*TEST

 79:    test:
 80:       nsize: 3

 82: TEST*/