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