Actual source code: ex5.c

  1: /*
  2:  Created by Huy Vo on 12/3/18.
  3: */
  4: static char help[] = "Test memory scalable AO.\n\n";

  6: #include<petsc.h>
  7: #include<petscvec.h>
  8: #include<petscmat.h>
  9: #include<petscao.h>

 11: int main(int argc, char *argv[])
 12: {
 13:   PetscInt              n_global = 16;
 14:   MPI_Comm              comm;
 15:   PetscLayout           layout;
 16:   PetscInt              local_size;
 17:   PetscInt              start, end;
 18:   PetscMPIInt           rank;
 19:   PetscInt              *app_indices,*petsc_indices,*ia,*ia0;
 20:   PetscInt              i;
 21:   AO                    app2petsc;
 22:   IS                    app_is, petsc_is;
 23:   const PetscInt        n_loc = 8;

 25:   PetscInitialize(&argc, &argv, (char *) 0, help);
 26:   comm = PETSC_COMM_WORLD;
 27:   MPI_Comm_rank(comm, &rank);

 29:   PetscLayoutCreate(comm, &layout);
 30:   PetscLayoutSetSize(layout, n_global);
 31:   PetscLayoutSetLocalSize(layout, PETSC_DECIDE);
 32:   PetscLayoutSetUp(layout);
 33:   PetscLayoutGetLocalSize(layout, &local_size);
 34:   PetscLayoutGetRange(layout, &start, &end);

 36:   PetscMalloc1(local_size,&app_indices);
 37:   PetscMalloc1(local_size,&petsc_indices);
 38:   /*  Add values for local indices for usual states */
 39:   for (i = 0; i < local_size; ++i) {
 40:     app_indices[i] = start + i;
 41:     petsc_indices[i] = end -1 - i;
 42:   }

 44:   /* Create the AO object that maps from lexicographic ordering to Petsc Vec ordering */
 45:   ISCreateGeneral(comm, local_size, &app_indices[0], PETSC_COPY_VALUES, &app_is);
 46:   ISCreateGeneral(comm, local_size, &petsc_indices[0], PETSC_COPY_VALUES, &petsc_is);
 47:   AOCreate(comm, &app2petsc);
 48:   AOSetIS(app2petsc, app_is, petsc_is);
 49:   AOSetType(app2petsc, AOMEMORYSCALABLE);
 50:   AOSetFromOptions(app2petsc);
 51:   ISDestroy(&app_is);
 52:   ISDestroy(&petsc_is);
 53:   AOView(app2petsc, PETSC_VIEWER_STDOUT_WORLD);

 55:   /* Test AOApplicationToPetsc */
 56:   PetscMalloc1(n_loc,&ia);
 57:   PetscMalloc1(n_loc,&ia0);
 58:   if (rank == 0) {
 59:     ia[0] = 0;
 60:     ia[1] = -1;
 61:     ia[2] = 1;
 62:     ia[3] = 2;
 63:     ia[4] = -1;
 64:     ia[5] = 4;
 65:     ia[6] = 5;
 66:     ia[7] = 6;
 67:   } else {
 68:     ia[0] = -1;
 69:     ia[1] = 8;
 70:     ia[2] = 9;
 71:     ia[3] = 10;
 72:     ia[4] = -1;
 73:     ia[5] = 12;
 74:     ia[6] = 13;
 75:     ia[7] = 14;
 76:   }
 77:   PetscArraycpy(ia0,ia,n_loc);

 79:   AOApplicationToPetsc(app2petsc, n_loc, ia);

 81:   for (i=0; i<n_loc; ++i) {
 82:     PetscSynchronizedPrintf(PETSC_COMM_WORLD,"proc = %d : %" PetscInt_FMT " -> %" PetscInt_FMT " \n", rank, ia0[i], ia[i]);
 83:   }
 84:   PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
 85:   AODestroy(&app2petsc);
 86:   PetscLayoutDestroy(&layout);
 87:   PetscFree(app_indices);
 88:   PetscFree(petsc_indices);
 89:   PetscFree(ia);
 90:   PetscFree(ia0);
 91:   PetscFinalize();
 92:   return 0;
 93: }

 95: /*TEST

 97:    test:
 98:      nsize: 2

100: TEST*/