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.\n\n";
5: #include petscao.h
6: #include petscbt.h
10: int main(int argc,char **argv)
11: {
12: int n = 2,nglobal,bs = 2,*keys,*data,ierr,rank,size,i,start;
13: PetscReal *gd;
14: AOData aodata;
15: PetscViewer binary;
16: PetscBT ld;
18: PetscInitialize(&argc,&argv,(char*)0,help);
19: PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank); n = n + rank;
22: MPI_Allreduce(&n,&nglobal,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
23: MPI_Comm_size(PETSC_COMM_WORLD,&size);
25: /*
26: Create a database with two sets of keys
27: */
28: AODataCreateBasic(PETSC_COMM_WORLD,&aodata);
29: AODataKeyAdd(aodata,"key1",PETSC_DECIDE,nglobal);
30: AODataKeyAdd(aodata,"key2",PETSC_DECIDE,nglobal);
32: /* allocate space for the keys each processor will provide */
33: PetscMalloc(n*sizeof(int),&keys);
35: /*
36: We assign the first set of keys (0 to 2) to processor 0, etc.
37: This computes the first local key on each processor
38: */
39: MPI_Scan(&n,&start,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
40: start -= n;
42: for (i=0; i<n; i++) {
43: keys[i] = start + i;
44: }
46: /*
47: Allocate data for the first key and first segment
48: */
49: PetscMalloc(bs*n*sizeof(int),&data);
50: for (i=0; i<n; i++) {
51: data[2*i] = -(start + i);
52: data[2*i+1] = -(start + i) - 10000;
53: }
54: AODataSegmentAdd(aodata,"key1","seg1",bs,n,keys,data,PETSC_INT);
55: PetscFree(data);
57: /*
58: Allocate data for first key and second segment
59: */
60: bs = 3;
61: PetscMalloc(bs*n*sizeof(PetscReal),&gd);
62: for (i=0; i<n; i++) {
63: gd[3*i] = -(start + i);
64: gd[3*i+1] = -(start + i) - 10000;
65: gd[3*i+2] = -(start + i) - 100000;
66: }
67: AODataSegmentAdd(aodata,"key1","seg2",bs,n,keys,gd,PETSC_REAL);
69: /*
70: Allocate data for first key and third segment
71: */
72: bs = 1;
73: PetscBTCreate(n,ld);
74: for (i=0; i<n; i++) {
75: if (i % 2) PetscBTSet(ld,i);
76: }
77: AODataSegmentAdd(aodata,"key1","seg3",bs,n,keys,ld,PETSC_LOGICAL);
78: PetscBTDestroy(ld);
80: /*
81: Use same data for second key and first segment
82: */
83: bs = 3;
84: AODataSegmentAdd(aodata,"key2","seg1",bs,n,keys,gd,PETSC_REAL);
85: PetscFree(gd);
86: PetscFree(keys);
88: AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);
90: /*
91: Save the database to a file
92: */
93: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"dataoutput",PETSC_BINARY_CREATE,&binary);
94: AODataView(aodata,binary);
95: PetscViewerDestroy(binary);
96:
97: AODataDestroy(aodata);
99: PetscFinalize();
100: return 0;
101: }
102: