Actual source code: ex5.c
2: static char help[] = "Tests AODataRemap(). \n\n";
4: #include petscao.h
8: int main(int argc,char **argv)
9: {
10: PetscInt n,nglobal,bs = 1,*keys,*data,i,start,*news;
12: PetscMPIInt rank,size;
13: AOData aodata;
14: AO ao;
16: PetscInitialize(&argc,&argv,(char*)0,help);
17: PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
19: MPI_Comm_rank(PETSC_COMM_WORLD,&rank); n = rank + 2;
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 one key and one segment
25: */
26: AODataCreateBasic(PETSC_COMM_WORLD,&aodata);
28: /*
29: Put one segment in the key
30: */
31: AODataKeyAdd(aodata,"key1",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[i] = start + i + 1; /* the data is the neighbor to the right */
53: }
54: data[n-1] = 0; /* make it periodic */
55: AODataSegmentAdd(aodata,"key1","key1",bs,n,keys,data,PETSC_INT);
56: PetscFree(data);
57: PetscFree(keys);
59: /*
60: View the database
61: */
62: AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);
63:
64: /*
65: Remap the database so that i -> nglobal - i - 1
66: */
67: PetscMalloc(n*sizeof(PetscInt),&news);
68: for (i=0; i<n; i++) {
69: news[i] = nglobal - i - start - 1;
70: }
71: AOCreateBasic(PETSC_COMM_WORLD,n,news,PETSC_NULL,&ao);
72: PetscFree(news);
73: AODataKeyRemap(aodata,"key1",ao);
74: AODestroy(ao);
75: AODataView(aodata,PETSC_VIEWER_STDOUT_WORLD);
76: AODataDestroy(aodata);
78: PetscFinalize();
79: return 0;
80: }
81: