1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Demonstrates constructing an application ordering.\n\n"; 3c4762a1bSJed Brown 4c4762a1bSJed Brown #include <petscao.h> 5c4762a1bSJed Brown #include <petscviewer.h> 6c4762a1bSJed Brown 7c4762a1bSJed Brown int main(int argc,char **argv) 8c4762a1bSJed Brown { 9c4762a1bSJed Brown PetscInt n = 5; 10c4762a1bSJed Brown PetscMPIInt rank,size; 11c4762a1bSJed Brown IS ispetsc,isapp; 12c4762a1bSJed Brown AO ao; 13c4762a1bSJed Brown 14*327415f7SBarry Smith PetscFunctionBeginUser; 159566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 169566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 179566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 189566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 19c4762a1bSJed Brown 20c4762a1bSJed Brown /* create the index sets */ 219566063dSJacob Faibussowitsch PetscCall(ISCreateStride(PETSC_COMM_WORLD,n,rank,size,&ispetsc)); 229566063dSJacob Faibussowitsch PetscCall(ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&isapp)); 23c4762a1bSJed Brown 24c4762a1bSJed Brown /* create the application ordering */ 259566063dSJacob Faibussowitsch PetscCall(AOCreateBasicIS(isapp,ispetsc,&ao)); 26c4762a1bSJed Brown 279566063dSJacob Faibussowitsch PetscCall(AOView(ao,PETSC_VIEWER_STDOUT_WORLD)); 28c4762a1bSJed Brown 299566063dSJacob Faibussowitsch PetscCall(ISView(ispetsc,PETSC_VIEWER_STDOUT_WORLD)); 309566063dSJacob Faibussowitsch PetscCall(ISView(isapp,PETSC_VIEWER_STDOUT_WORLD)); 319566063dSJacob Faibussowitsch PetscCall(AOPetscToApplicationIS(ao,ispetsc)); 329566063dSJacob Faibussowitsch PetscCall(ISView(isapp,PETSC_VIEWER_STDOUT_WORLD)); 339566063dSJacob Faibussowitsch PetscCall(ISView(ispetsc,PETSC_VIEWER_STDOUT_WORLD)); 34c4762a1bSJed Brown 359566063dSJacob Faibussowitsch PetscCall(ISDestroy(&ispetsc)); 369566063dSJacob Faibussowitsch PetscCall(ISDestroy(&isapp)); 37c4762a1bSJed Brown 389566063dSJacob Faibussowitsch PetscCall(AODestroy(&ao)); 399566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 40b122ec5aSJacob Faibussowitsch return 0; 41c4762a1bSJed Brown } 42c4762a1bSJed Brown 43c4762a1bSJed Brown /*TEST 44c4762a1bSJed Brown 45c4762a1bSJed Brown test: 46c4762a1bSJed Brown nsize: 2 47c4762a1bSJed Brown 48c4762a1bSJed Brown TEST*/ 49