1 /* 2 Created by Huy Vo on 12/3/18. 3 */ 4 static char help[] = "Test memory scalable AO.\n\n"; 5 6 #include<petsc.h> 7 #include<petscvec.h> 8 #include<petscmat.h> 9 #include<petscao.h> 10 11 int main(int argc, char *argv[]) 12 { 13 PetscInt ierr; 14 PetscInt n_global = 16; 15 MPI_Comm comm; 16 PetscLayout layout; 17 PetscInt local_size; 18 PetscInt start, end; 19 PetscMPIInt rank; 20 PetscInt *app_indices,*petsc_indices,*ia,*ia0; 21 PetscInt i; 22 AO app2petsc; 23 IS app_is, petsc_is; 24 const PetscInt n_loc = 8; 25 26 ierr = PetscInitialize(&argc, &argv, (char *) 0, help); if (ierr) return ierr; 27 comm = PETSC_COMM_WORLD; 28 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 29 30 ierr = PetscLayoutCreate(comm, &layout);CHKERRQ(ierr); 31 ierr = PetscLayoutSetSize(layout, n_global);CHKERRQ(ierr); 32 ierr = PetscLayoutSetLocalSize(layout, PETSC_DECIDE);CHKERRQ(ierr); 33 ierr = PetscLayoutSetUp(layout);CHKERRQ(ierr); 34 ierr = PetscLayoutGetLocalSize(layout, &local_size);CHKERRQ(ierr); 35 ierr = PetscLayoutGetRange(layout, &start, &end);CHKERRQ(ierr); 36 37 ierr = PetscMalloc1(local_size,&app_indices);CHKERRQ(ierr); 38 ierr = PetscMalloc1(local_size,&petsc_indices);CHKERRQ(ierr); 39 /* Add values for local indices for usual states */ 40 for (i = 0; i < local_size; ++i) { 41 app_indices[i] = start + i; 42 petsc_indices[i] = end -1 - i; 43 } 44 45 /* Create the AO object that maps from lexicographic ordering to Petsc Vec ordering */ 46 ierr = ISCreateGeneral(comm, local_size, &app_indices[0], PETSC_COPY_VALUES, &app_is);CHKERRQ(ierr); 47 ierr = ISCreateGeneral(comm, local_size, &petsc_indices[0], PETSC_COPY_VALUES, &petsc_is);CHKERRQ(ierr); 48 ierr = AOCreate(comm, &app2petsc);CHKERRQ(ierr); 49 ierr = AOSetIS(app2petsc, app_is, petsc_is);CHKERRQ(ierr); 50 ierr = AOSetType(app2petsc, AOMEMORYSCALABLE);CHKERRQ(ierr); 51 ierr = AOSetFromOptions(app2petsc);CHKERRQ(ierr); 52 ierr = ISDestroy(&app_is);CHKERRQ(ierr); 53 ierr = ISDestroy(&petsc_is);CHKERRQ(ierr); 54 ierr = AOView(app2petsc, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 55 56 /* Test AOApplicationToPetsc */ 57 ierr = PetscMalloc1(n_loc,&ia);CHKERRQ(ierr); 58 ierr = PetscMalloc1(n_loc,&ia0);CHKERRQ(ierr); 59 if (!rank) { 60 ia[0] = 0; 61 ia[1] = -1; 62 ia[2] = 1; 63 ia[3] = 2; 64 ia[4] = -1; 65 ia[5] = 4; 66 ia[6] = 5; 67 ia[7] = 6; 68 } else { 69 ia[0] = -1; 70 ia[1] = 8; 71 ia[2] = 9; 72 ia[3] = 10; 73 ia[4] = -1; 74 ia[5] = 12; 75 ia[6] = 13; 76 ia[7] = 14; 77 } 78 ierr = PetscArraycpy(ia0,ia,n_loc);CHKERRQ(ierr); 79 80 ierr = AOApplicationToPetsc(app2petsc, n_loc, ia);CHKERRQ(ierr); 81 82 for (i=0; i<n_loc; ++i) { 83 ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD,"proc = %d : %D -> %D \n", rank, ia0[i], ia[i]);CHKERRQ(ierr); 84 } 85 ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);CHKERRQ(ierr); 86 ierr = AODestroy(&app2petsc);CHKERRQ(ierr); 87 ierr = PetscLayoutDestroy(&layout);CHKERRQ(ierr); 88 ierr = PetscFree(app_indices);CHKERRQ(ierr); 89 ierr = PetscFree(petsc_indices);CHKERRQ(ierr); 90 ierr = PetscFree(ia);CHKERRQ(ierr); 91 ierr = PetscFree(ia0);CHKERRQ(ierr); 92 ierr = PetscFinalize(); 93 return ierr; 94 } 95 96 /*TEST 97 98 test: 99 nsize: 2 100 101 TEST*/ 102