Actual source code: ex5.c
1: /*$Id: ex5.c,v 1.19 2001/04/10 19:37:18 bsmith Exp $*/
3: static char help[] = "Tests AODataRemap(). nn";
5: #include petscao.h
7: int main(int argc,char **argv)
8: {
9: int n,nglobal,bs = 1,*keys,*data,ierr,rank,size,i,start,*news;
10: AOData aodata;
11: AO ao;
13: PetscInitialize(&argc,&argv,(char*)0,help);
14: PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
16: MPI_Comm_rank(PETSC_COMM_WORLD,&rank); n = rank + 2;
17: MPI_Allreduce(&n,&nglobal,1,MPI_INT,MPI_SUM,PETSC_COMM_WORLD);
18: MPI_Comm_size(PETSC_COMM_WORLD,&size);
20: /*
21: Create a database with one key and one segment
22: */
23: AODataCreateBasic(PETSC_COMM_WORLD,&aodata);
25: /*
26: Put one segment in the key
27: */
28: AODataKeyAdd(aodata,"key1",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[i] = start + i + 1; /* the data is the neighbor to the right */
50: }
51: data[n-1] = 0; /* make it periodic */
52: AODataSegmentAdd(aodata,"key1","key1",bs,n,keys,data,PETSC_INT);
53: PetscFree(data);
54: PetscFree(keys);
56: /*
57: View the database
58: */
59: AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);
60:
61: /*
62: Remap the database so that i -> nglobal - i - 1
63: */
64: PetscMalloc(n*sizeof(int),&news);
65: for (i=0; i<n; i++) {
66: news[i] = nglobal - i - start - 1;
67: }
68: AOCreateBasic(PETSC_COMM_WORLD,n,news,PETSC_NULL,&ao);
69: PetscFree(news);
70: AODataKeyRemap(aodata,"key1",ao);
71: AODestroy(ao);
72: AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);
73: AODataDestroy(aodata);
75: PetscFinalize();
76: return 0;
77: }
78: