1*c4762a1bSJed Brown 2*c4762a1bSJed Brown static char help[] = "Tests application ordering.\n\n"; 3*c4762a1bSJed Brown 4*c4762a1bSJed Brown #include <petscsys.h> 5*c4762a1bSJed Brown #include <petscao.h> 6*c4762a1bSJed Brown #include <petscviewer.h> 7*c4762a1bSJed Brown 8*c4762a1bSJed Brown int main(int argc,char **argv) 9*c4762a1bSJed Brown { 10*c4762a1bSJed Brown PetscErrorCode ierr; 11*c4762a1bSJed Brown PetscMPIInt rank,size; 12*c4762a1bSJed Brown PetscInt n,*ispetsc,*isapp,start,N,i; 13*c4762a1bSJed Brown AO ao; 14*c4762a1bSJed Brown 15*c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 16*c4762a1bSJed Brown ierr = PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);CHKERRQ(ierr); 17*c4762a1bSJed Brown ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); n = rank + 2; 18*c4762a1bSJed Brown ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 19*c4762a1bSJed Brown 20*c4762a1bSJed Brown /* create the orderings */ 21*c4762a1bSJed Brown ierr = PetscMalloc2(n,&ispetsc,n,&isapp);CHKERRQ(ierr); 22*c4762a1bSJed Brown 23*c4762a1bSJed Brown ierr = MPI_Scan(&n,&start,1,MPIU_INT,MPI_SUM,PETSC_COMM_WORLD);CHKERRQ(ierr); 24*c4762a1bSJed Brown ierr = MPI_Allreduce(&n,&N,1,MPIU_INT,MPI_SUM,PETSC_COMM_WORLD);CHKERRQ(ierr); 25*c4762a1bSJed Brown start -= n; 26*c4762a1bSJed Brown 27*c4762a1bSJed Brown for (i=0; i<n; i++) { 28*c4762a1bSJed Brown ispetsc[i] = start + i; 29*c4762a1bSJed Brown isapp[i] = N - start - i - 1; 30*c4762a1bSJed Brown } 31*c4762a1bSJed Brown 32*c4762a1bSJed Brown /* create the application ordering */ 33*c4762a1bSJed Brown ierr = AOCreateBasic(PETSC_COMM_WORLD,n,isapp,ispetsc,&ao);CHKERRQ(ierr); 34*c4762a1bSJed Brown ierr = AOView(ao,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 35*c4762a1bSJed Brown 36*c4762a1bSJed Brown /* check the mapping */ 37*c4762a1bSJed Brown ierr = AOPetscToApplication(ao,n,ispetsc);CHKERRQ(ierr); 38*c4762a1bSJed Brown for (i=0; i<n; i++) { 39*c4762a1bSJed Brown if (ispetsc[i] != isapp[i]) { 40*c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_WORLD,"[%d] Problem with mapping %D to %D\n",rank,i,ispetsc[i]);CHKERRQ(ierr); 41*c4762a1bSJed Brown } 42*c4762a1bSJed Brown } 43*c4762a1bSJed Brown ierr = PetscFree2(ispetsc,isapp);CHKERRQ(ierr); 44*c4762a1bSJed Brown 45*c4762a1bSJed Brown ierr = AODestroy(&ao);CHKERRQ(ierr); 46*c4762a1bSJed Brown ierr = PetscFinalize(); 47*c4762a1bSJed Brown return ierr; 48*c4762a1bSJed Brown } 49*c4762a1bSJed Brown 50*c4762a1bSJed Brown 51*c4762a1bSJed Brown 52*c4762a1bSJed Brown 53*c4762a1bSJed Brown 54*c4762a1bSJed Brown /*TEST 55*c4762a1bSJed Brown 56*c4762a1bSJed Brown test: 57*c4762a1bSJed Brown 58*c4762a1bSJed Brown test: 59*c4762a1bSJed Brown suffix: 2 60*c4762a1bSJed Brown nsize: 2 61*c4762a1bSJed Brown 62*c4762a1bSJed Brown test: 63*c4762a1bSJed Brown suffix: 3 64*c4762a1bSJed Brown nsize: 3 65*c4762a1bSJed Brown 66*c4762a1bSJed Brown TEST*/ 67