Actual source code: ex1.c
2: static char help[] = "Demonstrates constructing an application ordering.\n\n";
4: #include <petscsys.h>
5: #include <petscao.h>
6: #include <petscviewer.h>
8: int main(int argc, char **argv)
9: {
10: PetscInt i, n = 5;
11: PetscInt getpetsc[] = {0, 3, 4}, getapp[] = {2, 1, 9, 7};
12: PetscInt getpetsc1[] = {0, 3, 4}, getapp1[] = {2, 1, 9, 7};
13: PetscInt getpetsc2[] = {0, 3, 4}, getapp2[] = {2, 1, 9, 7};
14: PetscInt getpetsc3[] = {0, 3, 4}, getapp3[] = {2, 1, 9, 7};
15: PetscInt getpetsc4[] = {0, 3, 4}, getapp4[] = {2, 1, 9, 7};
16: PetscMPIInt rank, size;
17: IS ispetsc, isapp;
18: AO ao;
19: const PetscInt *app;
22: PetscInitialize(&argc, &argv, (char *)0, help);
23: PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
24: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
25: MPI_Comm_size(PETSC_COMM_WORLD, &size);
27: /* create the index sets */
28: ISCreateStride(PETSC_COMM_WORLD, n, rank, size, &isapp);
29: ISCreateStride(PETSC_COMM_WORLD, n, n * rank, 1, &ispetsc); /* natural numbering */
31: /* create the application ordering */
32: AOCreateBasicIS(isapp, ispetsc, &ao);
33: AOView(ao, PETSC_VIEWER_STDOUT_WORLD);
35: AOPetscToApplication(ao, 4, getapp);
36: AOApplicationToPetsc(ao, 3, getpetsc);
38: PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] 2,1,9,7 PetscToApplication %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, getapp[0], getapp[1], getapp[2], getapp[3]);
39: PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] 0,3,4 ApplicationToPetsc %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, getpetsc[0], getpetsc[1], getpetsc[2]);
40: PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT);
41: AODestroy(&ao);
43: /* test MemoryScalable ao */
44: /*-------------------------*/
45: PetscPrintf(PETSC_COMM_WORLD, "\nTest AOCreateMemoryScalable: \n");
46: AOCreateMemoryScalableIS(isapp, ispetsc, &ao);
47: AOView(ao, PETSC_VIEWER_STDOUT_WORLD);
49: AOPetscToApplication(ao, 4, getapp1);
50: AOApplicationToPetsc(ao, 3, getpetsc1);
52: /* Check accuracy */;
56: AODestroy(&ao);
58: /* test MemoryScalable ao: ispetsc = NULL */
59: /*-----------------------------------------------*/
60: PetscPrintf(PETSC_COMM_WORLD, "\nTest AOCreateMemoryScalable with ispetsc=NULL:\n");
61: AOCreateMemoryScalableIS(isapp, NULL, &ao);
63: AOView(ao, PETSC_VIEWER_STDOUT_WORLD);
65: AOPetscToApplication(ao, 4, getapp2);
66: AOApplicationToPetsc(ao, 3, getpetsc2);
68: /* Check accuracy */;
71: AODestroy(&ao);
73: /* test AOCreateMemoryScalable() ao: */
74: ISGetIndices(isapp, &app);
75: AOCreateMemoryScalable(PETSC_COMM_WORLD, n, app, NULL, &ao);
76: ISRestoreIndices(isapp, &app);
78: AOPetscToApplication(ao, 4, getapp4);
79: AOApplicationToPetsc(ao, 3, getpetsc4);
81: /* Check accuracy */;
84: AODestroy(&ao);
86: /* test general API */
87: /*------------------*/
88: PetscPrintf(PETSC_COMM_WORLD, "\nTest general API: \n");
89: AOCreate(PETSC_COMM_WORLD, &ao);
90: AOSetIS(ao, isapp, ispetsc);
91: AOSetType(ao, AOMEMORYSCALABLE);
92: AOSetFromOptions(ao);
94: /* ispetsc and isapp are nolonger used. */
95: ISDestroy(&ispetsc);
96: ISDestroy(&isapp);
98: AOPetscToApplication(ao, 4, getapp3);
99: AOApplicationToPetsc(ao, 3, getpetsc3);
101: PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] 2,1,9,7 PetscToApplication %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, getapp3[0], getapp3[1], getapp3[2], getapp3[3]);
102: PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] 0,3,4 ApplicationToPetsc %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n", rank, getpetsc3[0], getpetsc3[1], getpetsc3[2]);
103: PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT);
105: /* Check accuracy */;
109: AODestroy(&ao);
110: PetscFinalize();
111: return 0;
112: }
114: /*TEST
116: test:
118: test:
119: suffix: 2
120: nsize: 2
122: test:
123: suffix: 3
124: nsize: 3
126: test:
127: suffix: 4
128: nsize: 3
129: args: -ao_type basic
130: output_file: output/ex1_3.out
132: TEST*/