Actual source code: ex1.c


  2: static char help[] = "Used to benchmark changes to the PETSc VecScatter routines\n\n";
  3: #include <petscksp.h>
  4: extern PetscErrorCode  PetscLogView_VecScatter(PetscViewer);

  6: int main(int argc,char **args)
  7: {
  8:   KSP            ksp;
  9:   Mat            A;
 10:   Vec            x,b;
 11:   PetscViewer    fd;
 12:   char           file[PETSC_MAX_PATH_LEN];
 13:   PetscBool      flg,preload = PETSC_TRUE;

 15:   PetscInitialize(&argc,&args,(char*)0,help);
 16:   PetscLogDefaultBegin();
 17:   PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg);

 19:   PetscPreLoadBegin(preload,"Load system");

 21:   /*
 22:      Load the matrix and vector; then destroy the viewer.
 23:   */
 24:   MatCreate(PETSC_COMM_WORLD,&A);
 25:   MatSetFromOptions(A);
 26:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
 27:   MatLoad(A,fd);
 28:   PetscViewerDestroy(&fd);

 30:   MatCreateVecs(A,&x,&b);
 31:   VecSetFromOptions(b);
 32:   VecSet(b,1.0);

 34:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 35:   KSPSetFromOptions(ksp);
 36:   KSPSetOperators(ksp,A,A);
 37:   KSPSetUp(ksp);
 38:   KSPSetUpOnBlocks(ksp);

 40:   PetscPreLoadStage("KSPSolve");
 41:   KSPSolve(ksp,b,x);

 43:   MatDestroy(&A);
 44:   VecDestroy(&b);
 45:   VecDestroy(&x);
 46:   KSPDestroy(&ksp);
 47:   PetscPreLoadEnd();
 48:   PetscLogView_VecScatter(PETSC_VIEWER_STDOUT_WORLD);

 50:   PetscFinalize();
 51:   return 0;
 52: }

 54: #include <petsctime.h>
 55: #include <petsc/private/petscimpl.h>
 56: #include <petsc/private/vecimpl.h>
 57: #include <petsc/private/kspimpl.h>
 58: #include <petscmachineinfo.h>
 59: #include <petscconfiginfo.h>
 60: /*
 61:    This is a special log viewer that prints out detailed information only for the VecScatter routines
 62: */
 63: typedef enum { COUNT,TIME,NUMMESS,MESSLEN,REDUCT,FLOPS} Stats;
 64: PetscErrorCode  PetscLogView_VecScatter(PetscViewer viewer)
 65: {
 66:   MPI_Comm           comm       = PetscObjectComm((PetscObject) viewer);
 67:   PetscEventPerfInfo *eventInfo = NULL;
 68:   PetscLogDouble     locTotalTime,stats[6],maxstats[6],minstats[6],sumstats[6],avetime,ksptime;
 69:   PetscStageLog      stageLog;
 70:   const int          stage = 2;
 71:   int                event,events[] = {VEC_ScatterBegin,VEC_ScatterEnd};
 72:   PetscMPIInt        rank,size;
 73:   PetscErrorCode     ierr;
 74:   PetscInt           i;
 75:   char               arch[128],hostname[128],username[128],pname[PETSC_MAX_PATH_LEN],date[128],version[256];

 77:   PetscTime(&locTotalTime);  locTotalTime -= petsc_BaseTime;
 78:   MPI_Comm_size(comm, &size);
 79:   MPI_Comm_rank(comm, &rank);
 80:   PetscLogGetStageLog(&stageLog);
 81:   PetscViewerASCIIPrintf(viewer,"numProcs   = %d\n",size);

 83:   PetscGetArchType(arch,sizeof(arch));
 84:   PetscGetHostName(hostname,sizeof(hostname));
 85:   PetscGetUserName(username,sizeof(username));
 86:   PetscGetProgramName(pname,sizeof(pname));
 87:   PetscGetDate(date,sizeof(date));
 88:   PetscGetVersion(version,sizeof(version));
 89:   PetscViewerASCIIPrintf(viewer,"%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date);
 90:   PetscViewerASCIIPrintf(viewer, "Using %s\n", version);
 91:   PetscViewerASCIIPrintf(viewer, "Configure options: %s",petscconfigureoptions);
 92:   PetscViewerASCIIPrintf(viewer, "%s", petscmachineinfo);
 93:   PetscViewerASCIIPrintf(viewer, "%s", petsccompilerinfo);
 94:   PetscViewerASCIIPrintf(viewer, "%s", petsccompilerflagsinfo);
 95:   PetscViewerASCIIPrintf(viewer, "%s", petsclinkerinfo);
 96:   PetscViewerASCIIPrintf(viewer, "%s\n", PETSC_MPICC_SHOW);
 97:   PetscOptionsView(NULL,viewer);
 98: #if defined(PETSC_HAVE_HWLOC)
 99:   PetscProcessPlacementView(viewer);
100: #endif
101:   PetscViewerASCIIPrintf(viewer, "----------------------------------------------------\n");

103:   PetscViewerASCIIPrintf(viewer,"                Time     Min to Max Range   Proportion of KSP\n");

105:   eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo;
106:   MPI_Allreduce(&eventInfo[KSP_Solve].time,&ksptime,1,MPIU_PETSCLOGDOUBLE,MPI_SUM,PETSC_COMM_WORLD);
107:   ksptime = ksptime/size;

109:   for (i=0; i<(int)(sizeof(events)/sizeof(int)); i++) {
110:     event = events[i];
111:     stats[COUNT]   = eventInfo[event].count;
112:     stats[TIME]    = eventInfo[event].time;
113:     stats[NUMMESS] = eventInfo[event].numMessages;
114:     stats[MESSLEN] = eventInfo[event].messageLength;
115:     stats[REDUCT]  = eventInfo[event].numReductions;
116:     stats[FLOPS]   = eventInfo[event].flops;
117:     MPI_Allreduce(stats,maxstats,6,MPIU_PETSCLOGDOUBLE,MPI_MAX,PETSC_COMM_WORLD);
118:     MPI_Allreduce(stats,minstats,6,MPIU_PETSCLOGDOUBLE,MPI_MIN,PETSC_COMM_WORLD);
119:     MPI_Allreduce(stats,sumstats,6,MPIU_PETSCLOGDOUBLE,MPI_SUM,PETSC_COMM_WORLD);

121:     avetime  = sumstats[1]/size;
122:     PetscViewerASCIIPrintf(viewer,"%s %4.2e   -%5.1f %% %5.1f %%   %4.2e %%\n",stageLog->eventLog->eventInfo[event].name,
123:                                   avetime,100.*(avetime-minstats[1])/avetime,100.*(maxstats[1]-avetime)/avetime,100.*avetime/ksptime);
124:   }
125:   PetscViewerFlush(viewer);
126:   return 0;
127: }

129: /*TEST

131:    build:
132:      requires: defined(PETSC_USE_LOG)

134:    test:
135:      TODO: need to implement

137: TEST*/