Actual source code: ex3.c
1: /*$Id: ex3.c,v 1.29 2001/08/10 03:34:29 bsmith Exp $*/
3: static char help[] = "Tests AOData.nn";
5: #include petscao.h
6: #include petscbt.h
8: int main(int argc,char **argv)
9: {
10: int n = 2,nglobal,bs = 2,*keys,*data,ierr,rank,size,i,start;
11: PetscReal *gd;
12: AOData aodata;
13: PetscViewer binary;
14: PetscBT ld;
16: PetscInitialize(&argc,&argv,(char*)0,help);
17: PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
19: MPI_Comm_rank(PETSC_COMM_WORLD,&rank); n = n + rank;
20: MPI_Allreduce(&n,&nglobal,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
21: MPI_Comm_size(PETSC_COMM_WORLD,&size);
23: /*
24: Create a database with two sets of keys
25: */
26: AODataCreateBasic(PETSC_COMM_WORLD,&aodata);
27: AODataKeyAdd(aodata,"key1",PETSC_DECIDE,nglobal);
28: AODataKeyAdd(aodata,"key2",PETSC_DECIDE,nglobal);
30: /* allocate space for the keys each processor will provide */
31: PetscMalloc(n*sizeof(int),&keys);
33: /*
34: We assign the first set of keys (0 to 2) to processor 0, etc.
35: This computes the first local key on each processor
36: */
37: MPI_Scan(&n,&start,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
38: start -= n;
40: for (i=0; i<n; i++) {
41: keys[i] = start + i;
42: }
44: /*
45: Allocate data for the first key and first segment
46: */
47: PetscMalloc(bs*n*sizeof(int),&data);
48: for (i=0; i<n; i++) {
49: data[2*i] = -(start + i);
50: data[2*i+1] = -(start + i) - 10000;
51: }
52: AODataSegmentAdd(aodata,"key1","seg1",bs,n,keys,data,PETSC_INT);
53: PetscFree(data);
55: /*
56: Allocate data for first key and second segment
57: */
58: bs = 3;
59: PetscMalloc(bs*n*sizeof(PetscReal),&gd);
60: for (i=0; i<n; i++) {
61: gd[3*i] = -(start + i);
62: gd[3*i+1] = -(start + i) - 10000;
63: gd[3*i+2] = -(start + i) - 100000;
64: }
65: AODataSegmentAdd(aodata,"key1","seg2",bs,n,keys,gd,PETSC_REAL);
67: /*
68: Allocate data for first key and third segment
69: */
70: bs = 1;
71: PetscBTCreate(n,ld);
72: for (i=0; i<n; i++) {
73: if (i % 2) PetscBTSet(ld,i);
74: }
75: AODataSegmentAdd(aodata,"key1","seg3",bs,n,keys,ld,PETSC_LOGICAL);
76: PetscBTDestroy(ld);
78: /*
79: Use same data for second key and first segment
80: */
81: bs = 3;
82: AODataSegmentAdd(aodata,"key2","seg1",bs,n,keys,gd,PETSC_REAL);
83: PetscFree(gd);
84: PetscFree(keys);
86: AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);
88: /*
89: Save the database to a file
90: */
91: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"dataoutput",PETSC_BINARY_CREATE,&binary);
92: AODataView(aodata,binary);
93: PetscViewerDestroy(binary);
94:
95: AODataDestroy(aodata);
97: PetscFinalize();
98: return 0;
99: }
100: