Actual source code: ex1.c
petsc-3.12.2 2019-11-22
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];
14: PetscBool flg,preload = PETSC_TRUE;
16: PetscInitialize(&argc,&args,(char*)0,help);
17: PetscLogDefaultBegin();
18: PetscOptionsGetString(NULL,NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);
20: PetscPreLoadBegin(preload,"Load system");
22: /*
23: Load the matrix and vector; then destroy the viewer.
24: */
25: MatCreate(PETSC_COMM_WORLD,&A);
26: MatSetFromOptions(A);
27: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
28: MatLoad(A,fd);
29: PetscViewerDestroy(&fd);
31: MatCreateVecs(A,&x,&b);
32: VecSetFromOptions(b);
33: VecSet(b,1.0);
35: KSPCreate(PETSC_COMM_WORLD,&ksp);
36: KSPSetFromOptions(ksp);
37: KSPSetOperators(ksp,A,A);
38: KSPSetUp(ksp);
39: KSPSetUpOnBlocks(ksp);
41: PetscPreLoadStage("KSPSolve");
42: KSPSolve(ksp,b,x);
44: MatDestroy(&A);
45: VecDestroy(&b);
46: VecDestroy(&x);
47: KSPDestroy(&ksp);
48: PetscPreLoadEnd();
49: PetscLogView_VecScatter(PETSC_VIEWER_STDOUT_WORLD);
51: PetscFinalize();
52: return ierr;
53: }
55: #include <petsctime.h>
56: #include <petsc/private/petscimpl.h>
57: #include <petsc/private/vecimpl.h>
58: #include <petsc/private/kspimpl.h>
59: #include <petscmachineinfo.h>
60: #include <petscconfiginfo.h>
61: /*
62: This is a special log viewer that prints out detailed information only for the VecScatter routines
63: */
64: typedef enum { COUNT,TIME,NUMMESS,MESSLEN,REDUCT,FLOPS} Stats;
65: PetscErrorCode PetscLogView_VecScatter(PetscViewer viewer)
66: {
67: MPI_Comm comm = PetscObjectComm((PetscObject) viewer);
68: PetscEventPerfInfo *eventInfo = NULL;
69: PetscLogDouble locTotalTime,stats[6],maxstats[6],minstats[6],sumstats[6],avetime,ksptime;
70: PetscStageLog stageLog;
71: const int stage = 2;
72: int event,events[] = {VEC_ScatterBegin,VEC_ScatterEnd};
73: PetscMPIInt rank,size;
74: PetscErrorCode ierr;
75: PetscInt i;
76: char arch[128],hostname[128],username[128],pname[PETSC_MAX_PATH_LEN],date[128],version[256];
79: PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime;
80: MPI_Comm_size(comm, &size);
81: MPI_Comm_rank(comm, &rank);
82: PetscLogGetStageLog(&stageLog);
83: PetscViewerASCIIPrintf(viewer,"numProcs = %d\n",size);
85: PetscGetArchType(arch,sizeof(arch));
86: PetscGetHostName(hostname,sizeof(hostname));
87: PetscGetUserName(username,sizeof(username));
88: PetscGetProgramName(pname,sizeof(pname));
89: PetscGetDate(date,sizeof(date));
90: PetscGetVersion(version,sizeof(version));
91: PetscViewerASCIIPrintf(viewer,"%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date);
92: PetscViewerASCIIPrintf(viewer, "Using %s\n", version);
93: PetscViewerASCIIPrintf(viewer, "Configure options: %s",petscconfigureoptions);
94: PetscViewerASCIIPrintf(viewer, "%s", petscmachineinfo);
95: PetscViewerASCIIPrintf(viewer, "%s", petsccompilerinfo);
96: PetscViewerASCIIPrintf(viewer, "%s", petsccompilerflagsinfo);
97: PetscViewerASCIIPrintf(viewer, "%s", petsclinkerinfo);
98: PetscViewerASCIIPrintf(viewer, "%s\n", PETSC_MPICC_SHOW);
99: PetscOptionsView(NULL,viewer);
100: #if defined(PETSC_HAVE_HWLOC)
101: PetscProcessPlacementView(viewer);
102: #endif
103: PetscViewerASCIIPrintf(viewer, "----------------------------------------------------\n");
105: PetscViewerASCIIPrintf(viewer," Time Min to Max Range Proportion of KSP\n");
107: eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo;
108: MPI_Allreduce(&eventInfo[KSP_Solve].time,&ksptime,1,MPIU_PETSCLOGDOUBLE,MPI_SUM,PETSC_COMM_WORLD);
109: ksptime = ksptime/size;
111: for (i=0; i<(int)(sizeof(events)/sizeof(int)); i++) {
112: event = events[i];
113: stats[COUNT] = eventInfo[event].count;
114: stats[TIME] = eventInfo[event].time;
115: stats[NUMMESS] = eventInfo[event].numMessages;
116: stats[MESSLEN] = eventInfo[event].messageLength;
117: stats[REDUCT] = eventInfo[event].numReductions;
118: stats[FLOPS] = eventInfo[event].flops;
119: MPI_Allreduce(stats,maxstats,6,MPIU_PETSCLOGDOUBLE,MPI_MAX,PETSC_COMM_WORLD);
120: MPI_Allreduce(stats,minstats,6,MPIU_PETSCLOGDOUBLE,MPI_MIN,PETSC_COMM_WORLD);
121: MPI_Allreduce(stats,sumstats,6,MPIU_PETSCLOGDOUBLE,MPI_SUM,PETSC_COMM_WORLD);
123: avetime = sumstats[1]/size;
124: PetscViewerASCIIPrintf(viewer,"%s %4.2e -%5.1f %% %5.1f %% %4.2e %%\n",stageLog->eventLog->eventInfo[event].name,
125: avetime,100.*(avetime-minstats[1])/avetime,100.*(maxstats[1]-avetime)/avetime,100.*avetime/ksptime);
126: }
127: PetscViewerFlush(viewer);
128: return(0);
129: }
131: /*TEST
133: build:
134: requires: define(PETSC_USE_LOG) c99
136: test:
137: TODO: need to implement
139: TEST*/