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