Actual source code: ex3.c

  2: static char help[] = "Tests AOData.\n\n";

 4:  #include petscao.h
 5:  #include petscbt.h

  9: int main(int argc,char **argv)
 10: {
 11:   PetscInt       n = 2,nglobal,bs = 2,*keys,*data,i,start;
 13:   PetscMPIInt    rank,size;
 14:   PetscReal      *gd;
 15:   AOData         aodata;
 16:   PetscViewer    binary;
 17:   PetscBT        ld;

 19:   PetscInitialize(&argc,&argv,(char*)0,help);
 20:   PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);

 22:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank); n = n + rank;
 23:   MPI_Allreduce(&n,&nglobal,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
 24:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 26:   /*
 27:        Create a database with two sets of keys 
 28:   */
 29:   AODataCreateBasic(PETSC_COMM_WORLD,&aodata);
 30:   AODataKeyAdd(aodata,"key1",PETSC_DECIDE,nglobal);
 31:   AODataKeyAdd(aodata,"key2",PETSC_DECIDE,nglobal);

 33:   /* allocate space for the keys each processor will provide */
 34:   PetscMalloc(n*sizeof(PetscInt),&keys);

 36:   /*
 37:      We assign the first set of keys (0 to 2) to processor 0, etc.
 38:      This computes the first local key on each processor
 39:   */
 40:   MPI_Scan(&n,&start,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
 41:   start -= n;

 43:   for (i=0; i<n; i++) {
 44:     keys[i]     = start + i;
 45:   }

 47:   /* 
 48:       Allocate data for the first key and first segment 
 49:   */
 50:   PetscMalloc(bs*n*sizeof(PetscInt),&data);
 51:   for (i=0; i<n; i++) {
 52:     data[2*i]   = -(start + i);
 53:     data[2*i+1] = -(start + i) - 10000;
 54:   }
 55:   AODataSegmentAdd(aodata,"key1","seg1",bs,n,keys,data,PETSC_INT);
 56:   PetscFree(data);

 58:   /*
 59:       Allocate data for first key and second segment 
 60:   */
 61:   bs   = 3;
 62:   PetscMalloc(bs*n*sizeof(PetscReal),&gd);
 63:   for (i=0; i<n; i++) {
 64:     gd[3*i]   = -(start + i);
 65:     gd[3*i+1] = -(start + i) - 10000;
 66:     gd[3*i+2] = -(start + i) - 100000;
 67:   }
 68:   AODataSegmentAdd(aodata,"key1","seg2",bs,n,keys,gd,PETSC_REAL);

 70:   /*
 71:       Allocate data for first key and third segment 
 72:   */
 73:   bs   = 1;
 74:   PetscBTCreate(n,ld);
 75:   for (i=0; i<n; i++) {
 76:     if (i % 2) {PetscBTSet(ld,i);}
 77:   }
 78:   AODataSegmentAdd(aodata,"key1","seg3",bs,n,keys,ld,PETSC_LOGICAL);
 79:   PetscBTDestroy(ld);

 81:   /*
 82:        Use same data for second key and first segment 
 83:   */
 84:   bs   = 3;
 85:   AODataSegmentAdd(aodata,"key2","seg1",bs,n,keys,gd,PETSC_REAL);
 86:   PetscFree(gd);
 87:   PetscFree(keys);

 89:   AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);

 91:   /*
 92:         Save the database to a file
 93:   */
 94:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,"dataoutput",PETSC_FILE_CREATE,&binary);
 95:   AODataView(aodata,binary);
 96:   PetscViewerDestroy(binary);
 97: 
 98:   AODataDestroy(aodata);

100:   PetscFinalize();
101:   return 0;
102: }
103: