Actual source code: ex1.c
1: static char help[] = "Used to benchmark changes to the PETSc VecScatter routines\n\n";
2: #include <petscksp.h>
3: extern PetscErrorCode PetscLogView_VecScatter(PetscViewer);
5: int main(int argc, char **args)
6: {
7: KSP ksp;
8: Mat A;
9: Vec x, b;
10: PetscViewer fd;
11: char file[PETSC_MAX_PATH_LEN];
12: PetscBool flg, preload = PETSC_TRUE;
14: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
15: PetscCall(PetscLogDefaultBegin());
16: PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
18: PetscPreLoadBegin(preload, "Load system");
20: /*
21: Load the matrix and vector; then destroy the viewer.
22: */
23: PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
24: PetscCall(MatSetFromOptions(A));
25: PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd));
26: PetscCall(MatLoad(A, fd));
27: PetscCall(PetscViewerDestroy(&fd));
29: PetscCall(MatCreateVecs(A, &x, &b));
30: PetscCall(VecSetFromOptions(b));
31: PetscCall(VecSet(b, 1.0));
33: PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
34: PetscCall(KSPSetFromOptions(ksp));
35: PetscCall(KSPSetOperators(ksp, A, A));
36: PetscCall(KSPSetUp(ksp));
37: PetscCall(KSPSetUpOnBlocks(ksp));
39: PetscPreLoadStage("KSPSolve");
40: PetscCall(KSPSolve(ksp, b, x));
42: PetscCall(MatDestroy(&A));
43: PetscCall(VecDestroy(&b));
44: PetscCall(VecDestroy(&x));
45: PetscCall(KSPDestroy(&ksp));
46: PetscPreLoadEnd();
47: PetscCall(PetscLogView_VecScatter(PETSC_VIEWER_STDOUT_WORLD));
49: PetscCall(PetscFinalize());
50: return 0;
51: }
53: #include <petsctime.h>
54: #include <petsc/private/petscimpl.h>
55: #include <petsc/private/vecimpl.h>
56: #include <petsc/private/kspimpl.h>
57: #include <petscmachineinfo.h>
58: #include <petscconfiginfo.h>
59: /*
60: This is a special log viewer that prints out detailed information only for the VecScatter routines
61: */
62: typedef enum {
63: COUNT,
64: TIME,
65: NUMMESS,
66: MESSLEN,
67: REDUCT,
68: FLOPS
69: } Stats;
71: PetscErrorCode PetscLogView_VecScatter(PetscViewer viewer)
72: {
73: MPI_Comm comm = PetscObjectComm((PetscObject)viewer);
74: PetscEventPerfInfo eventInfo;
75: PetscLogDouble locTotalTime, stats[6], maxstats[6], minstats[6], sumstats[6], avetime, ksptime;
76: const int stage = 2;
77: int events[] = {VEC_ScatterBegin, VEC_ScatterEnd};
78: PetscMPIInt rank, size;
79: char arch[128], hostname[128], username[128], pname[PETSC_MAX_PATH_LEN], date[128], version[256];
81: PetscFunctionBegin;
82: PetscCall(PetscTime(&locTotalTime));
83: locTotalTime -= petsc_BaseTime;
84: PetscCallMPI(MPI_Comm_size(comm, &size));
85: PetscCallMPI(MPI_Comm_rank(comm, &rank));
86: PetscCall(PetscViewerASCIIPrintf(viewer, "numProcs = %d\n", size));
88: PetscCall(PetscGetArchType(arch, sizeof(arch)));
89: PetscCall(PetscGetHostName(hostname, sizeof(hostname)));
90: PetscCall(PetscGetUserName(username, sizeof(username)));
91: PetscCall(PetscGetProgramName(pname, sizeof(pname)));
92: PetscCall(PetscGetDate(date, sizeof(date)));
93: PetscCall(PetscGetVersion(version, sizeof(version)));
94: PetscCall(PetscViewerASCIIPrintf(viewer, "%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date));
95: PetscCall(PetscViewerASCIIPrintf(viewer, "Using %s\n", version));
96: PetscCall(PetscViewerASCIIPrintf(viewer, "Configure options: %s", petscconfigureoptions));
97: PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petscmachineinfo));
98: PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petsccompilerinfo));
99: PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petsccompilerflagsinfo));
100: PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petsclinkerinfo));
101: PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", PETSC_MPICC_SHOW));
102: PetscCall(PetscOptionsView(NULL, viewer));
103: #if defined(PETSC_HAVE_HWLOC)
104: PetscCall(PetscProcessPlacementView(viewer));
105: #endif
106: PetscCall(PetscViewerASCIIPrintf(viewer, "----------------------------------------------------\n"));
108: PetscCall(PetscViewerASCIIPrintf(viewer, " Time Min to Max Range Proportion of KSP\n"));
110: PetscCall(PetscLogEventGetPerfInfo(stage, KSP_Solve, &eventInfo));
111: PetscCall(MPIU_Allreduce(&eventInfo.time, &ksptime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_WORLD));
112: ksptime = ksptime / size;
114: for (size_t i = 0; i < PETSC_STATIC_ARRAY_LENGTH(events); i++) {
115: PetscEventPerfInfo eventInfo;
116: const char *name;
118: PetscCall(PetscLogEventGetPerfInfo(stage, events[i], &eventInfo));
119: PetscCall(PetscLogEventGetName(events[i], &name));
120: stats[COUNT] = eventInfo.count;
121: stats[TIME] = eventInfo.time;
122: stats[NUMMESS] = eventInfo.numMessages;
123: stats[MESSLEN] = eventInfo.messageLength;
124: stats[REDUCT] = eventInfo.numReductions;
125: stats[FLOPS] = eventInfo.flops;
126: PetscCall(MPIU_Allreduce(stats, maxstats, 6, MPIU_PETSCLOGDOUBLE, MPI_MAX, PETSC_COMM_WORLD));
127: PetscCall(MPIU_Allreduce(stats, minstats, 6, MPIU_PETSCLOGDOUBLE, MPI_MIN, PETSC_COMM_WORLD));
128: PetscCall(MPIU_Allreduce(stats, sumstats, 6, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_WORLD));
130: avetime = sumstats[1] / size;
131: PetscCall(PetscViewerASCIIPrintf(viewer, "%s %4.2e -%5.1f %% %5.1f %% %4.2e %%\n", name, avetime, 100. * (avetime - minstats[1]) / avetime, 100. * (maxstats[1] - avetime) / avetime, 100. * avetime / ksptime));
132: }
133: PetscCall(PetscViewerFlush(viewer));
134: PetscFunctionReturn(PETSC_SUCCESS);
135: }
137: /*TEST
139: build:
140: requires: defined(PETSC_USE_LOG)
142: test:
143: TODO: need to implement
145: TEST*/