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;

 21:   PetscInitialize(&argc,&argv,(char*)0,help);
 22:   PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
 23:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 24:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 26:   /* create the index sets */
 27:   ISCreateStride(PETSC_COMM_WORLD,n,rank,size,&isapp);
 28:   ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&ispetsc); /* natural numbering */

 30:   /* create the application ordering */
 31:   AOCreateBasicIS(isapp,ispetsc,&ao);
 32:   AOView(ao,PETSC_VIEWER_STDOUT_WORLD);

 34:   AOPetscToApplication(ao,4,getapp);
 35:   AOApplicationToPetsc(ao,3,getpetsc);

 37:   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]);
 38:   PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] 0,3,4 ApplicationToPetsc %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",rank,getpetsc[0],getpetsc[1],getpetsc[2]);
 39:   PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
 40:   AODestroy(&ao);

 42:   /* test MemoryScalable ao */
 43:   /*-------------------------*/
 44:   PetscPrintf(PETSC_COMM_WORLD,"\nTest AOCreateMemoryScalable: \n");
 45:   AOCreateMemoryScalableIS(isapp,ispetsc,&ao);
 46:   AOView(ao,PETSC_VIEWER_STDOUT_WORLD);

 48:   AOPetscToApplication(ao,4,getapp1);
 49:   AOApplicationToPetsc(ao,3,getpetsc1);

 51:   /* Check accuracy */;
 52:   for (i=0; i<4; i++) {
 54:   }
 55:   for (i=0; i<3; i++) {
 57:   }

 59:   AODestroy(&ao);

 61:   /* test MemoryScalable ao: ispetsc = NULL */
 62:   /*-----------------------------------------------*/
 63:   PetscPrintf(PETSC_COMM_WORLD,"\nTest AOCreateMemoryScalable with ispetsc=NULL:\n");
 64:   AOCreateMemoryScalableIS(isapp,NULL,&ao);

 66:   AOView(ao,PETSC_VIEWER_STDOUT_WORLD);

 68:   AOPetscToApplication(ao,4,getapp2);
 69:   AOApplicationToPetsc(ao,3,getpetsc2);

 71:   /* Check accuracy */;
 72:   for (i=0; i<4; i++) {
 74:   }
 75:   for (i=0; i<3; i++) {
 77:   }
 78:   AODestroy(&ao);

 80:   /* test AOCreateMemoryScalable() ao: */
 81:   ISGetIndices(isapp,&app);
 82:   AOCreateMemoryScalable(PETSC_COMM_WORLD,n,app,NULL,&ao);
 83:   ISRestoreIndices(isapp,&app);

 85:   AOPetscToApplication(ao,4,getapp4);
 86:   AOApplicationToPetsc(ao,3,getpetsc4);

 88:   /* Check accuracy */;
 89:   for (i=0; i<4; i++) {
 91:   }
 92:   for (i=0; i<3; i++) {
 94:   }
 95:   AODestroy(&ao);

 97:   /* test general API */
 98:   /*------------------*/
 99:   PetscPrintf(PETSC_COMM_WORLD,"\nTest general API: \n");
100:   AOCreate(PETSC_COMM_WORLD,&ao);
101:   AOSetIS(ao,isapp,ispetsc);
102:   AOSetType(ao,AOMEMORYSCALABLE);
103:   AOSetFromOptions(ao);

105:   /* ispetsc and isapp are nolonger used. */
106:   ISDestroy(&ispetsc);
107:   ISDestroy(&isapp);

109:   AOPetscToApplication(ao,4,getapp3);
110:   AOApplicationToPetsc(ao,3,getpetsc3);

112:   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]);
113:   PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] 0,3,4 ApplicationToPetsc %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",rank,getpetsc3[0],getpetsc3[1],getpetsc3[2]);
114:   PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);

116:   /* Check accuracy */;
117:   for (i=0; i<4; i++) {
119:   }
120:   for (i=0; i<3; i++) {
122:   }

124:   AODestroy(&ao);
125:   PetscFinalize();
126:   return 0;
127: }

129: /*TEST

131:    test:

133:    test:
134:       suffix: 2
135:       nsize: 2

137:    test:
138:       suffix: 3
139:       nsize: 3

141:    test:
142:       suffix: 4
143:       nsize: 3
144:       args: -ao_type basic
145:       output_file: output/ex1_3.out

147: TEST*/