Actual source code: ex110.c

petsc-3.3-p6 2013-02-11
  1: static char help[] = "Testing MatCreateMPIAIJWithSplitArrays().\n\n";

  3: #include <petscmat.h>
  4: #include <../src/mat/impls/aij/mpi/mpiaij.h>

  8: int main(int argc,char **argv) {
  9:   Mat            A,B;
 10:   PetscInt       i,j,column;
 11:   PetscInt       *di,*dj,*oi,*oj;
 12:   PetscScalar    *oa,*da,value;
 13:   PetscRandom    rctx;
 15:   PetscBool      equal;
 16:   Mat_SeqAIJ     *daij,*oaij;
 17:   Mat_MPIAIJ     *Ampiaij;
 18:   PetscMPIInt    size,rank;

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

 25:   /* Create a mpiaij matrix for checking */
 26:   MatCreateAIJ(PETSC_COMM_WORLD,5,5,PETSC_DECIDE,PETSC_DECIDE,0,PETSC_NULL,0,PETSC_NULL,&A);
 27:   MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
 28:   MatSetUp(A);
 29:   PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
 30:   PetscRandomSetFromOptions(rctx);
 31: 
 32:   for (i=5*rank; i<5*rank+5; i++) {
 33:     for (j=0; j<5*size; j++){
 34:       PetscRandomGetValue(rctx,&value);
 35:       column = (PetscInt) (5*size*PetscRealPart(value));
 36:       PetscRandomGetValue(rctx,&value);
 37:       MatSetValues(A,1,&i,1,&column,&value,INSERT_VALUES);
 38:     }
 39:   }
 40:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 41:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);

 43:   Ampiaij = (Mat_MPIAIJ*) A->data;
 44:   daij    = (Mat_SeqAIJ*) Ampiaij->A->data;
 45:   oaij    = (Mat_SeqAIJ*) Ampiaij->B->data;
 46: 
 47:   di      = daij->i;
 48:   dj      = daij->j;
 49:   da      = daij->a;

 51:   oi      = oaij->i;
 52:   oa      = oaij->a;
 53:   PetscMalloc(oi[5]*sizeof(PetscInt),&oj);
 54:   PetscMemcpy(oj,oaij->j,oi[5]*sizeof(PetscInt));
 55:   /* modify the column entries in the non-diagonal portion back to global numbering */
 56:   for (i=0; i<oi[5]; i++) {
 57:     oj[i] = Ampiaij->garray[oj[i]];
 58:   }

 60:   MatCreateMPIAIJWithSplitArrays(PETSC_COMM_WORLD,5,5,PETSC_DETERMINE,PETSC_DETERMINE,di,dj,da,oi,oj,oa,&B);
 61:   MatSetUp(B);
 62:   MatEqual(A,B,&equal);

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

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