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 {
64: COUNT,
65: TIME,
66: NUMMESS,
67: MESSLEN,
68: REDUCT,
69: FLOPS
70: } Stats;
71: PetscErrorCode PetscLogView_VecScatter(PetscViewer viewer)
72: {
73: MPI_Comm comm = PetscObjectComm((PetscObject)viewer);
74: PetscEventPerfInfo *eventInfo = NULL;
75: PetscLogDouble locTotalTime, stats[6], maxstats[6], minstats[6], sumstats[6], avetime, ksptime;
76: PetscStageLog stageLog;
77: const int stage = 2;
78: int event, events[] = {VEC_ScatterBegin, VEC_ScatterEnd};
79: PetscMPIInt rank, size;
80: PetscInt i;
81: char arch[128], hostname[128], username[128], pname[PETSC_MAX_PATH_LEN], date[128], version[256];
83: PetscTime(&locTotalTime);
84: locTotalTime -= petsc_BaseTime;
85: MPI_Comm_size(comm, &size);
86: MPI_Comm_rank(comm, &rank);
87: PetscLogGetStageLog(&stageLog);
88: PetscViewerASCIIPrintf(viewer, "numProcs = %d\n", size);
90: PetscGetArchType(arch, sizeof(arch));
91: PetscGetHostName(hostname, sizeof(hostname));
92: PetscGetUserName(username, sizeof(username));
93: PetscGetProgramName(pname, sizeof(pname));
94: PetscGetDate(date, sizeof(date));
95: PetscGetVersion(version, sizeof(version));
96: PetscViewerASCIIPrintf(viewer, "%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date);
97: PetscViewerASCIIPrintf(viewer, "Using %s\n", version);
98: PetscViewerASCIIPrintf(viewer, "Configure options: %s", petscconfigureoptions);
99: PetscViewerASCIIPrintf(viewer, "%s", petscmachineinfo);
100: PetscViewerASCIIPrintf(viewer, "%s", petsccompilerinfo);
101: PetscViewerASCIIPrintf(viewer, "%s", petsccompilerflagsinfo);
102: PetscViewerASCIIPrintf(viewer, "%s", petsclinkerinfo);
103: PetscViewerASCIIPrintf(viewer, "%s\n", PETSC_MPICC_SHOW);
104: PetscOptionsView(NULL, viewer);
105: #if defined(PETSC_HAVE_HWLOC)
106: PetscProcessPlacementView(viewer);
107: #endif
108: PetscViewerASCIIPrintf(viewer, "----------------------------------------------------\n");
110: PetscViewerASCIIPrintf(viewer, " Time Min to Max Range Proportion of KSP\n");
112: eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo;
113: MPI_Allreduce(&eventInfo[KSP_Solve].time, &ksptime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_WORLD);
114: ksptime = ksptime / size;
116: for (i = 0; i < (int)(sizeof(events) / sizeof(int)); i++) {
117: event = events[i];
118: stats[COUNT] = eventInfo[event].count;
119: stats[TIME] = eventInfo[event].time;
120: stats[NUMMESS] = eventInfo[event].numMessages;
121: stats[MESSLEN] = eventInfo[event].messageLength;
122: stats[REDUCT] = eventInfo[event].numReductions;
123: stats[FLOPS] = eventInfo[event].flops;
124: MPI_Allreduce(stats, maxstats, 6, MPIU_PETSCLOGDOUBLE, MPI_MAX, PETSC_COMM_WORLD);
125: MPI_Allreduce(stats, minstats, 6, MPIU_PETSCLOGDOUBLE, MPI_MIN, PETSC_COMM_WORLD);
126: MPI_Allreduce(stats, sumstats, 6, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_WORLD);
128: avetime = sumstats[1] / size;
129: PetscViewerASCIIPrintf(viewer, "%s %4.2e -%5.1f %% %5.1f %% %4.2e %%\n", stageLog->eventLog->eventInfo[event].name, avetime, 100. * (avetime - minstats[1]) / avetime, 100. * (maxstats[1] - avetime) / avetime, 100. * avetime / ksptime);
130: }
131: PetscViewerFlush(viewer);
132: return 0;
133: }
135: /*TEST
137: build:
138: requires: defined(PETSC_USE_LOG)
140: test:
141: TODO: need to implement
143: TEST*/