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*9566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 15*9566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 16*9566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 17*9566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 18c4762a1bSJed Brown 19c4762a1bSJed Brown /* create the index sets */ 20*9566063dSJacob Faibussowitsch PetscCall(ISCreateStride(PETSC_COMM_WORLD,n,rank,size,&ispetsc)); 21*9566063dSJacob Faibussowitsch PetscCall(ISCreateStride(PETSC_COMM_WORLD,n,n*rank,1,&isapp)); 22c4762a1bSJed Brown 23c4762a1bSJed Brown /* create the application ordering */ 24*9566063dSJacob Faibussowitsch PetscCall(AOCreateBasicIS(isapp,ispetsc,&ao)); 25c4762a1bSJed Brown 26*9566063dSJacob Faibussowitsch PetscCall(AOView(ao,PETSC_VIEWER_STDOUT_WORLD)); 27c4762a1bSJed Brown 28*9566063dSJacob Faibussowitsch PetscCall(ISView(ispetsc,PETSC_VIEWER_STDOUT_WORLD)); 29*9566063dSJacob Faibussowitsch PetscCall(ISView(isapp,PETSC_VIEWER_STDOUT_WORLD)); 30*9566063dSJacob Faibussowitsch PetscCall(AOPetscToApplicationIS(ao,ispetsc)); 31*9566063dSJacob Faibussowitsch PetscCall(ISView(isapp,PETSC_VIEWER_STDOUT_WORLD)); 32*9566063dSJacob Faibussowitsch PetscCall(ISView(ispetsc,PETSC_VIEWER_STDOUT_WORLD)); 33c4762a1bSJed Brown 34*9566063dSJacob Faibussowitsch PetscCall(ISDestroy(&ispetsc)); 35*9566063dSJacob Faibussowitsch PetscCall(ISDestroy(&isapp)); 36c4762a1bSJed Brown 37*9566063dSJacob Faibussowitsch PetscCall(AODestroy(&ao)); 38*9566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 39b122ec5aSJacob Faibussowitsch return 0; 40c4762a1bSJed Brown } 41c4762a1bSJed Brown 42c4762a1bSJed Brown /*TEST 43c4762a1bSJed Brown 44c4762a1bSJed Brown test: 45c4762a1bSJed Brown nsize: 2 46c4762a1bSJed Brown 47c4762a1bSJed Brown TEST*/ 48