1 2 static char help[] = "Demonstrates constructing an application ordering.\n\n"; 3 4 #include <petscsys.h> 5 #include <petscao.h> 6 #include <petscviewer.h> 7 8 int main(int argc,char **argv) 9 { 10 PetscErrorCode ierr; 11 PetscInt i,n = 5; 12 PetscInt getpetsc[] = {0,3,4},getapp[] = {2,1,9,7}; 13 PetscInt getpetsc1[] = {0,3,4},getapp1[] = {2,1,9,7}; 14 PetscInt getpetsc2[] = {0,3,4},getapp2[] = {2,1,9,7}; 15 PetscInt getpetsc3[] = {0,3,4},getapp3[] = {2,1,9,7}; 16 PetscInt getpetsc4[] = {0,3,4},getapp4[] = {2,1,9,7}; 17 PetscMPIInt rank,size; 18 IS ispetsc,isapp; 19 AO ao; 20 const PetscInt *app; 21 22 ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 23 CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 24 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 25 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 26 27 /* create the index sets */ 28 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,n,rank,size,&isapp)); 29 CHKERRQ(ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&ispetsc)); /* natural numbering */ 30 31 /* create the application ordering */ 32 CHKERRQ(AOCreateBasicIS(isapp,ispetsc,&ao)); 33 CHKERRQ(AOView(ao,PETSC_VIEWER_STDOUT_WORLD)); 34 35 CHKERRQ(AOPetscToApplication(ao,4,getapp)); 36 CHKERRQ(AOApplicationToPetsc(ao,3,getpetsc)); 37 38 CHKERRQ(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 CHKERRQ(PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] 0,3,4 ApplicationToPetsc %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",rank,getpetsc[0],getpetsc[1],getpetsc[2])); 40 CHKERRQ(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT)); 41 CHKERRQ(AODestroy(&ao)); 42 43 /* test MemoryScalable ao */ 44 /*-------------------------*/ 45 CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"\nTest AOCreateMemoryScalable: \n")); 46 CHKERRQ(AOCreateMemoryScalableIS(isapp,ispetsc,&ao)); 47 CHKERRQ(AOView(ao,PETSC_VIEWER_STDOUT_WORLD)); 48 49 CHKERRQ(AOPetscToApplication(ao,4,getapp1)); 50 CHKERRQ(AOApplicationToPetsc(ao,3,getpetsc1)); 51 52 /* Check accuracy */; 53 for (i=0; i<4; i++) { 54 PetscCheckFalse(getapp1[i] != getapp[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getapp1 %" PetscInt_FMT " != getapp %" PetscInt_FMT,getapp1[i],getapp[i]); 55 } 56 for (i=0; i<3; i++) { 57 PetscCheckFalse(getpetsc1[i] != getpetsc[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getpetsc1 %" PetscInt_FMT " != getpetsc %" PetscInt_FMT,getpetsc1[i],getpetsc[i]); 58 } 59 60 CHKERRQ(AODestroy(&ao)); 61 62 /* test MemoryScalable ao: ispetsc = NULL */ 63 /*-----------------------------------------------*/ 64 CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"\nTest AOCreateMemoryScalable with ispetsc=NULL:\n")); 65 CHKERRQ(AOCreateMemoryScalableIS(isapp,NULL,&ao)); 66 67 CHKERRQ(AOView(ao,PETSC_VIEWER_STDOUT_WORLD)); 68 69 CHKERRQ(AOPetscToApplication(ao,4,getapp2)); 70 CHKERRQ(AOApplicationToPetsc(ao,3,getpetsc2)); 71 72 /* Check accuracy */; 73 for (i=0; i<4; i++) { 74 PetscCheckFalse(getapp2[i] != getapp[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getapp2 %" PetscInt_FMT " != getapp %" PetscInt_FMT,getapp2[i],getapp[i]); 75 } 76 for (i=0; i<3; i++) { 77 PetscCheckFalse(getpetsc2[i] != getpetsc[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getpetsc2 %" PetscInt_FMT " != getpetsc %" PetscInt_FMT,getpetsc2[i],getpetsc[i]); 78 } 79 CHKERRQ(AODestroy(&ao)); 80 81 /* test AOCreateMemoryScalable() ao: */ 82 CHKERRQ(ISGetIndices(isapp,&app)); 83 CHKERRQ(AOCreateMemoryScalable(PETSC_COMM_WORLD,n,app,NULL,&ao)); 84 CHKERRQ(ISRestoreIndices(isapp,&app)); 85 86 CHKERRQ(AOPetscToApplication(ao,4,getapp4)); 87 CHKERRQ(AOApplicationToPetsc(ao,3,getpetsc4)); 88 89 /* Check accuracy */; 90 for (i=0; i<4; i++) { 91 PetscCheckFalse(getapp4[i] != getapp[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getapp4 %" PetscInt_FMT " != getapp %" PetscInt_FMT,getapp4[i],getapp[i]); 92 } 93 for (i=0; i<3; i++) { 94 PetscCheckFalse(getpetsc4[i] != getpetsc[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getpetsc4 %" PetscInt_FMT " != getpetsc %" PetscInt_FMT,getpetsc4[i],getpetsc[i]); 95 } 96 CHKERRQ(AODestroy(&ao)); 97 98 /* test general API */ 99 /*------------------*/ 100 CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"\nTest general API: \n")); 101 CHKERRQ(AOCreate(PETSC_COMM_WORLD,&ao)); 102 CHKERRQ(AOSetIS(ao,isapp,ispetsc)); 103 CHKERRQ(AOSetType(ao,AOMEMORYSCALABLE)); 104 CHKERRQ(AOSetFromOptions(ao)); 105 106 /* ispetsc and isapp are nolonger used. */ 107 CHKERRQ(ISDestroy(&ispetsc)); 108 CHKERRQ(ISDestroy(&isapp)); 109 110 CHKERRQ(AOPetscToApplication(ao,4,getapp3)); 111 CHKERRQ(AOApplicationToPetsc(ao,3,getpetsc3)); 112 113 CHKERRQ(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])); 114 CHKERRQ(PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] 0,3,4 ApplicationToPetsc %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",rank,getpetsc3[0],getpetsc3[1],getpetsc3[2])); 115 CHKERRQ(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT)); 116 117 /* Check accuracy */; 118 for (i=0; i<4; i++) { 119 PetscCheckFalse(getapp3[i] != getapp[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getapp3 %" PetscInt_FMT " != getapp %" PetscInt_FMT,getapp3[i],getapp[i]); 120 } 121 for (i=0; i<3; i++) { 122 PetscCheckFalse(getpetsc3[i] != getpetsc[i],PETSC_COMM_SELF,PETSC_ERR_USER,"getpetsc3 %" PetscInt_FMT " != getpetsc %" PetscInt_FMT,getpetsc3[i],getpetsc[i]); 123 } 124 125 CHKERRQ(AODestroy(&ao)); 126 ierr = PetscFinalize(); 127 return ierr; 128 } 129 130 /*TEST 131 132 test: 133 134 test: 135 suffix: 2 136 nsize: 2 137 138 test: 139 suffix: 3 140 nsize: 3 141 142 test: 143 suffix: 4 144 nsize: 3 145 args: -ao_type basic 146 output_file: output/ex1_3.out 147 148 TEST*/ 149