1#!/usr/bin/env python3 2from __future__ import print_function 3import os,sys 4sys.path.append(os.path.join(os.environ['PETSC_DIR'], 'config')) 5sys.path.append(os.getcwd()) 6 7class CSVLog(object): pass 8 9class PETSc(object): 10 def __init__(self): 11 return 12 13 def dir(self): 14 '''Return the root directory for the PETSc tree (usually $PETSC_DIR)''' 15 # This should search for a valid PETSc 16 return os.environ['PETSC_DIR'] 17 18 def arch(self): 19 '''Return the PETSc build label (usually $PETSC_ARCH)''' 20 # This should be configurable 21 return os.environ['PETSC_ARCH'] 22 23 def mpiexec(self): 24 '''Return the path for the mpi launch executable''' 25 mpiexec = os.path.join(self.dir(), self.arch(), 'bin', 'mpiexec') 26 if not os.path.isfile(mpiexec): 27 return None 28 return mpiexec 29 30 def example(self, library, num): 31 '''Return the path to the executable for a given example number''' 32 return os.path.join(self.dir(), self.arch(), 'lib', library.lower(), 'ex'+str(num)) 33 34 def source(self, library, num, filenametail): 35 '''Return the path to the sources for a given example number''' 36 d = os.path.join(self.dir(), 'src', library.lower(), 'examples', 'tutorials') 37 name = 'ex'+str(num) 38 sources = [] 39 for f in os.listdir(d): 40 if f == name+'.c': 41 sources.insert(0, f) 42 elif f.startswith(name) and f.endswith(filenametail): 43 sources.append(f) 44 return map(lambda f: os.path.join(d, f), sources) 45 46class PETScExample(object): 47 def __init__(self, library, num, **defaultOptions): 48 self.petsc = PETSc() 49 self.library = library 50 self.num = num 51 self.opts = defaultOptions 52 return 53 54 @staticmethod 55 def runShellCommand(command, cwd = None, log = True): 56 import subprocess 57 58 Popen = subprocess.Popen 59 PIPE = subprocess.PIPE 60 if log: print('Executing: %s\n' % (command,)) 61 pipe = Popen(command, cwd=cwd, stdin=None, stdout=PIPE, stderr=PIPE, bufsize=-1, shell=True, universal_newlines=True) 62 (out, err) = pipe.communicate() 63 ret = pipe.returncode 64 return (out, err, ret) 65 66 def optionsToString(self, **opts): 67 '''Convert a dictionary of options to a command line argument string''' 68 a = [] 69 for key,value in opts.items(): 70 if value is None: 71 a.append('-'+key) 72 else: 73 a.append('-'+key+' '+str(value)) 74 return ' '.join(a) 75 76 def build(self, log = True): 77 sdir = os.path.join(self.petsc.dir(), 'src', self.library.lower(), 'examples', 'tutorials') 78 odir = os.getcwd() 79 os.chdir(sdir) 80 cmd = 'make ex'+str(self.num) 81 out, err, ret = self.runShellCommand(cmd, cwd = sdir, log = log) 82 if err: 83 raise RuntimeError('Unable to build example:\n'+err+out) 84 os.chdir(odir) 85 bdir = os.path.dirname(self.petsc.example(self.library, self.num)) 86 try: 87 os.makedirs(bdir) 88 except OSError: 89 if not os.path.isdir(bdir): 90 raise 91 cmd = 'mv '+os.path.join(sdir, 'ex'+str(self.num))+' '+self.petsc.example(self.library, self.num) 92 out, err, ret = self.runShellCommand(cmd, log = log) 93 if ret: 94 print(err) 95 print(out) 96 return 97 98 def run(self, numProcs = 1, log = True, **opts): 99 cmd = '' 100 if self.petsc.mpiexec() is not None: 101 cmd += self.petsc.mpiexec() + ' ' 102 numProcs = os.environ.get('NUM_RANKS', numProcs) 103 cmd += ' -n ' + str(numProcs) + ' ' 104 if 'PE_HOSTFILE' in os.environ: 105 cmd += ' -hostfile hostfile ' 106 cmd += ' '.join([self.petsc.example(self.library, self.num), self.optionsToString(**self.opts), self.optionsToString(**opts)]) 107 if 'batch' in opts and opts['batch']: 108 del opts['batch'] 109 from benchmarkBatch import generateBatchScript 110 filename = generateBatchScript(self.num, numProcs, 120, ' '+self.optionsToString(**self.opts)+' '+self.optionsToString(**opts)) 111 # Submit job 112 out, err, ret = self.runShellCommand('qsub -q gpu '+filename, log = log) 113 if ret: 114 print(err) 115 print(out) 116 else: 117 out, err, ret = self.runShellCommand(cmd, log = log) 118 if ret: 119 print(err) 120 print(out) 121 return out 122 123def processSummaryCSV(filename, defaultStage, eventNames, sizes, times, errors, events): 124 '''Process the CSV log summary into plot data''' 125 import csv 126 m = CSVLog() 127 setattr(m, 'Stages', {}) 128 with open(filename) as csvfile: 129 reader = csv.DictReader(csvfile) 130 for row in reader: 131 stageName = row["Stage Name"] 132 eventName = row["Event Name"] 133 rank = int(row["Rank"]) 134 if not stageName in m.Stages: m.Stages[stageName] = {} 135 if not eventName in m.Stages[stageName]: m.Stages[stageName][eventName] = {} 136 m.Stages[stageName][eventName][rank] = {"time" : float(row["Time"]), "numMessages": float(row["Num Messages"]), "messageLength": float(row["Message Length"]), "numReductions" : float(row["Num Reductions"]), "flop" : float(row["FLOP"])} 137 for i in range(8): 138 dname = "dof"+str(i) 139 ename = "e"+str(i) 140 if row[dname]: 141 if not "dof" in m.Stages[stageName][eventName][rank]: 142 m.Stages[stageName][eventName][rank]["dof"] = [] 143 m.Stages[stageName][eventName][rank]["error"] = [] 144 m.Stages[stageName][eventName][rank]["dof"].append(int(float(row[dname]))) 145 m.Stages[stageName][eventName][rank]["error"].append(float(row[ename])) 146 return m 147 148def processSummary(moduleName, defaultStage, eventNames, sizes, times, errors, events): 149 '''Process the Python log summary into plot data''' 150 if os.path.isfile(moduleName+'.csv'): 151 m = processSummaryCSV(moduleName+'.csv', defaultStage, eventNames, sizes, times, errors, events) 152 else: 153 m = __import__(moduleName) 154 # Total Time 155 times.append(m.Stages[defaultStage]["summary"][0]["time"]) 156 # Particular events 157 for name in eventNames: 158 if name.find(':') >= 0: 159 stageName, name = name.split(':', 1) 160 else: 161 stageName = defaultStage 162 stage = m.Stages[stageName] 163 if name in stage: 164 if not name in events: 165 events[name] = [] 166 event = stage[name][0] 167 etimes = [stage[name][p]["time"] for p in stage[name]] 168 eflops = [stage[name][p]["flop"] for p in stage[name]] 169 if "dof" in event: 170 sizes.append(event["dof"][0]) 171 errors.append(event["error"][0]) 172 try: 173 events[name].append((max(etimes), sum(eflops)/(max(etimes) * 1e6))) 174 except ZeroDivisionError: 175 events[name].append((max(etimes), 0)) 176 return 177 178def plotTime(library, num, eventNames, sizes, times, events): 179 from pylab import legend, plot, show, title, xlabel, ylabel 180 import numpy as np 181 182 arches = sizes.keys() 183 data = [] 184 for arch in arches: 185 data.append(sizes[arch]) 186 data.append(times[arch]) 187 plot(*data) 188 title('Performance on '+library+' Example '+str(num)) 189 xlabel('Number of Dof') 190 ylabel('Time (s)') 191 legend(arches, 'upper left', shadow = True) 192 show() 193 return 194 195def plotEventTime(library, num, eventNames, sizes, times, events, filename = None): 196 from pylab import close, legend, plot, savefig, show, title, xlabel, ylabel 197 import numpy as np 198 199 close() 200 arches = sizes.keys() 201 bs = events[arches[0]].keys()[0] 202 data = [] 203 names = [] 204 for event, color in zip(eventNames, ['b', 'g', 'r', 'y']): 205 for arch, style in zip(arches, ['-', ':']): 206 if event in events[arch][bs]: 207 names.append(arch+'-'+str(bs)+' '+event) 208 data.append(sizes[arch][bs]) 209 data.append(np.array(events[arch][bs][event])[:,0]) 210 data.append(color+style) 211 else: 212 print('Could not find %s in %s-%d events' % (event, arch, bs)) 213 print(data) 214 plot(*data) 215 title('Performance on '+library+' Example '+str(num)) 216 xlabel('Number of Dof') 217 ylabel('Time (s)') 218 legend(names, 'upper left', shadow = True) 219 if filename is None: 220 show() 221 else: 222 savefig(filename) 223 return 224 225def plotEventFlop(library, num, eventNames, sizes, times, events, filename = None): 226 from pylab import legend, plot, savefig, semilogy, show, title, xlabel, ylabel 227 import numpy as np 228 229 arches = sizes.keys() 230 bs = events[arches[0]].keys()[0] 231 data = [] 232 names = [] 233 for event, color in zip(eventNames, ['b', 'g', 'r', 'y']): 234 for arch, style in zip(arches, ['-', ':']): 235 if event in events[arch][bs]: 236 names.append(arch+'-'+str(bs)+' '+event) 237 data.append(sizes[arch][bs]) 238 data.append(1e-3*np.array(events[arch][bs][event])[:,1]) 239 data.append(color+style) 240 else: 241 print('Could not find %s in %s-%d events' % (event, arch, bs)) 242 semilogy(*data) 243 title('Performance on '+library+' Example '+str(num)) 244 xlabel('Number of Dof') 245 ylabel('Computation Rate (GF/s)') 246 legend(names, 'upper left', shadow = True) 247 if filename is None: 248 show() 249 else: 250 savefig(filename) 251 return 252 253def plotEventScaling(library, num, eventNames, procs, events, filename = None): 254 from pylab import legend, plot, savefig, semilogy, show, title, xlabel, ylabel 255 import numpy as np 256 257 arches = procs.keys() 258 bs = events[arches[0]].keys()[0] 259 data = [] 260 names = [] 261 for arch, style in zip(arches, ['-', ':']): 262 for event, color in zip(eventNames, ['b', 'g', 'r', 'y']): 263 if event in events[arch][bs]: 264 names.append(arch+'-'+str(bs)+' '+event) 265 data.append(procs[arch][bs]) 266 data.append(1e-3*np.array(events[arch][bs][event])[:,1]) 267 data.append(color+style) 268 else: 269 print('Could not find %s in %s-%d events' % (event, arch, bs)) 270 plot(*data) 271 title('Performance on '+library+' Example '+str(num)) 272 xlabel('Number of Processors') 273 ylabel('Computation Rate (GF/s)') 274 legend(names, 'upper left', shadow = True) 275 if filename is None: 276 show() 277 else: 278 savefig(filename) 279 return 280 281def plotSummaryLine(library, num, eventNames, sizes, times, events): 282 from pylab import legend, plot, show, title, xlabel, ylabel 283 import numpy as np 284 showTime = False 285 showEventTime = True 286 showEventFlops = True 287 arches = sizes.keys() 288 # Time 289 if showTime: 290 data = [] 291 for arch in arches: 292 data.append(sizes[arch]) 293 data.append(times[arch]) 294 plot(*data) 295 title('Performance on '+library+' Example '+str(num)) 296 xlabel('Number of Dof') 297 ylabel('Time (s)') 298 legend(arches, 'upper left', shadow = True) 299 show() 300 # Common event time 301 # We could make a stacked plot like Rio uses here 302 if showEventTime: 303 data = [] 304 names = [] 305 for event, color in zip(eventNames, ['b', 'g', 'r', 'y']): 306 for arch, style in zip(arches, ['-', ':']): 307 if event in events[arch]: 308 names.append(arch+' '+event) 309 data.append(sizes[arch]) 310 data.append(np.array(events[arch][event])[:,0]) 311 data.append(color+style) 312 else: 313 print('Could not find %s in %s events' % (event, arch)) 314 print(data) 315 plot(*data) 316 title('Performance on '+library+' Example '+str(num)) 317 xlabel('Number of Dof') 318 ylabel('Time (s)') 319 legend(names, 'upper left', shadow = True) 320 show() 321 # Common event flops 322 # We could make a stacked plot like Rio uses here 323 if showEventFlops: 324 data = [] 325 names = [] 326 for event, color in zip(eventNames, ['b', 'g', 'r', 'y']): 327 for arch, style in zip(arches, ['-', ':']): 328 if event in events[arch]: 329 names.append(arch+' '+event) 330 data.append(sizes[arch]) 331 data.append(np.array(events[arch][event])[:,1]) 332 data.append(color+style) 333 else: 334 print('Could not find %s in %s events' % (event, arch)) 335 plot(*data) 336 title('Performance on '+library+' Example '+str(num)) 337 xlabel('Number of Dof') 338 ylabel('Computation Rate (MF/s)') 339 legend(names, 'upper left', shadow = True) 340 show() 341 return 342 343def plotMeshConvergence(library, num, eventNames, sizes, times, errors, events): 344 import numpy as np 345 import matplotlib.pyplot as plt 346 data = [] 347 legends = [] 348 print(sizes) 349 print(errors) 350 for run in sizes: 351 rsizes = np.array(sizes[run]) 352 data.extend([rsizes, errors[run], rsizes, (errors[run][0]*rsizes[0]*2)*rsizes**(meshExp[run]/-2.0)]) 353 legends.extend(['Experiment '+run, r'Synthetic '+run+r' $\alpha = '+str(meshExp[run])+'$']) 354 SizeError = plt.loglog(*data) 355 plt.title(library+' ex'+str(num)+' Mesh Convergence') 356 plt.xlabel('Size') 357 plt.ylabel(r'Error $\|x - x^*\|_2$') 358 plt.legend(legends) 359 plt.show() 360 return 361 362def plotWorkPrecision(library, num, eventNames, sizes, times, errors, events): 363 import numpy as np 364 import matplotlib.pyplot as plt 365 data = [] 366 legends = [] 367 for run in times: 368 rtimes = np.array(times[run]) 369 data.extend([rtimes, errors[run], rtimes, (errors[run][0]*rtimes[0]*2)*rtimes**(timeExp[run])]) 370 legends.extend(['Experiment '+run, 'Synthetic '+run+' exponent '+str(timeExp[run])]) 371 TimeError = plt.loglog(*data) 372 plt.title(library+' ex'+str(num)+' Work Precision') 373 plt.xlabel('Time (s)') 374 plt.ylabel(r'Error $\|x - x^*\|_2$') 375 plt.legend(legends) 376 plt.show() 377 return 378 379def plotWorkPrecisionPareto(library, num, eventNames, sizes, times, errors, events): 380 import numpy as np 381 import matplotlib.pyplot as plt 382 data = [] 383 legends = [] 384 for run in times: 385 rtimes = np.array(times[run]) 386 data.extend([rtimes, errors[run]]) 387 legends.append('Experiment '+run) 388 TimeErrorPareto = plt.semilogy(*data) 389 plt.title(library+' ex'+str(num)+' Work Precision: Pareto Front') 390 plt.xlabel('Time (s)') 391 plt.ylabel(r'Error $\|x - x^*\|_2$') 392 plt.legend(legends) 393 plt.show() 394 return 395 396def plotSummaryBar(library, num, eventNames, sizes, times, events): 397 import numpy as np 398 import matplotlib.pyplot as plt 399 400 eventColors = ['b', 'g', 'r', 'y'] 401 arches = sizes.keys() 402 names = [] 403 N = len(sizes[arches[0]]) 404 width = 0.2 405 ind = np.arange(N) - 0.25 406 bars = {} 407 for arch in arches: 408 bars[arch] = [] 409 bottom = np.zeros(N) 410 for event, color in zip(eventNames, eventColors): 411 names.append(arch+' '+event) 412 times = np.array(events[arch][event])[:,0] 413 bars[arch].append(plt.bar(ind, times, width, color=color, bottom=bottom)) 414 bottom += times 415 ind += 0.3 416 417 plt.xlabel('Number of Dof') 418 plt.ylabel('Time (s)') 419 plt.title('GPU vs. CPU Performance on '+library+' Example '+str(num)) 420 plt.xticks(np.arange(N), map(str, sizes[arches[0]])) 421 #plt.yticks(np.arange(0,81,10)) 422 #plt.legend( (p1[0], p2[0]), ('Men', 'Women') ) 423 plt.legend([bar[0] for bar in bars[arches[0]]], eventNames, 'upper right', shadow = True) 424 425 plt.show() 426 return 427 428def processOptions(opts, name, n): 429 newopts = {} 430 for key, val in opts.items(): 431 val = opts[key] 432 if val and val.find('%') >= 0: 433 newval = val % (name.replace('/', '-'), n) 434 newopts[key] = newval 435 else: 436 newopts[key] = val 437 return newopts 438 439def getLogName(opts): 440 if 'log_view' in opts: 441 val = opts['log_view'] 442 s = val.find(':') 443 e = val.find(':', s+1) 444 logName = os.path.splitext(val[s+1:e])[0] 445 return logName 446 return None 447 448def run_DMDA(ex, name, opts, args, sizes, times, events, log=True, execute=True): 449 for n in map(int, args.size): 450 newopts = processOptions(opts, name, n) 451 if execute: 452 ex.run(log=log, da_grid_x=n, da_grid_y=n, **newopts) 453 processSummary(getLogName(newopts), args.stage, args.events, sizes[name], times[name], errors[name], events[name]) 454 return 455 456def run_DMPlex(ex, name, opts, args, sizes, times, events, log=True, execute=True): 457 newopts = processOptions(opts, name, args.refine) 458 if execute: 459 ex.run(log=log, dim=args.dim, snes_convergence_estimate=None, convest_num_refine=args.refine, interpolate=1, **newopts) 460 for r in range(args.refine+1): 461 stage = args.stage 462 if stage.find('%') >= 0: stage = stage % (r) 463 processSummary(getLogName(newopts), stage, args.events, sizes[name], times[name], errors[name], events[name]) 464 return 465 466def outputData(sizes, times, events, name = 'output.py'): 467 if os.path.exists(name): 468 base, ext = os.path.splitext(name) 469 num = 1 470 while os.path.exists(base+str(num)+ext): 471 num += 1 472 name = base+str(num)+ext 473 with file(name, 'w') as f: 474 f.write('#PETSC_ARCH='+os.environ['PETSC_ARCH']+' '+' '.join(sys.argv)+'\n') 475 f.write('sizes = '+repr(sizes)+'\n') 476 f.write('times = '+repr(times)+'\n') 477 f.write('events = '+repr(events)+'\n') 478 return 479 480if __name__ == '__main__': 481 import argparse 482 import __main__ 483 484 parser = argparse.ArgumentParser(description = 'PETSc Benchmarking', 485 epilog = 'This script runs src/<library>/tutorials/ex<num>, For more information, visit https://petsc.org/', 486 formatter_class = argparse.ArgumentDefaultsHelpFormatter) 487 parser.add_argument('--library', default='SNES', help='The PETSc library used in this example') 488 parser.add_argument('--num', type = int, default='5', help='The example number') 489 parser.add_argument('--module', default='summary', help='The module for timing output') 490 parser.add_argument('--stage', default='Main Stage', help='The default logging stage') 491 parser.add_argument('--events', nargs='+', help='Events to process') 492 parser.add_argument('--plotOnly',action='store_true', default=False, help='Flag to only plot existing data') 493 parser.add_argument('--batch', action='store_true', default=False, help='Generate batch files for the runs instead') 494 parser.add_argument('--daemon', action='store_true', default=False, help='Run as a daemon') 495 parser.add_argument('--gpulang', default='OpenCL', help='GPU Language to use: Either CUDA or OpenCL (default)') 496 parser.add_argument('--plots', nargs='+', help='List of plots to show') 497 subparsers = parser.add_subparsers(help='DM types') 498 499 parser_dmda = subparsers.add_parser('DMDA', help='Use a DMDA for the problem geometry') 500 parser_dmda.add_argument('--size', nargs='+', default=['10'], help='Grid size (implementation dependent)') 501 parser_dmda.add_argument('--comp', type = int, default='1', help='Number of field components') 502 parser_dmda.add_argument('runs', nargs='*', help='Run descriptions: <name>=<args>') 503 504 parser_dmmesh = subparsers.add_parser('DMPlex', help='Use a DMPlex for the problem geometry') 505 parser_dmmesh.add_argument('--dim', type = int, default='2', help='Spatial dimension') 506 parser_dmmesh.add_argument('--refine', type = int, default='0', help='Number of refinements') 507 parser_dmmesh.add_argument('runs', nargs='*', help='Run descriptions: <name>=<args>') 508 509 args = parser.parse_args() 510 if hasattr(args, 'comp'): 511 args.dmType = 'DMDA' 512 else: 513 args.dmType = 'DMPlex' 514 515 ex = PETScExample(args.library, args.num, preload='off') 516 if args.gpulang == 'CUDA': 517 source = ex.petsc.source(args.library, args.num, '.cu') 518 else: 519 source = ex.petsc.source(args.library, args.num, 'OpenCL.c') # Using the convention of OpenCL code residing in source files ending in 'OpenCL.c' (at least for snes/ex52) 520 sizes = {} 521 times = {} 522 errors = {} 523 meshExp = {} 524 timeExp = {} 525 events = {} 526 log = not args.daemon 527 528 if args.daemon: 529 import daemon 530 print('Starting daemon') 531 daemon.createDaemon('.') 532 533 for run in args.runs: 534 name, stropts = run.split('=', 1) 535 opts = dict([t if len(t) == 2 else (t[0], None) for t in [arg.split('=', 1) for arg in stropts.split(' ')]]) 536 #opts['log_view'] = 'summary.dat' if args.batch else ':'+args.module+'%s%d.py:ascii_info_detail' 537 opts['log_view'] = 'summary.dat' if args.batch else ':'+args.module+'%s%d.csv:ascii_csv' 538 meshExp[name] = float(opts['meshExp']) 539 timeExp[name] = float(opts['timeExp']) 540 sizes[name] = [] 541 times[name] = [] 542 errors[name] = [] 543 events[name] = {} 544 getattr(__main__, 'run_'+args.dmType)(ex, name, opts, args, sizes, times, events, log=log, execute=(not args.plotOnly)) 545 outputData(sizes, times, events) 546 if not args.batch and log: 547 for plot in args.plots: 548 print('Plotting ',plot) 549 getattr(__main__, 'plot'+plot)(args.library, args.num, args.events, sizes, times, errors, events) 550 551# ./src/benchmarks/benchmarkExample.py --events SNESSolve --plots MeshConvergence WorkPrecision WorkPrecisionPareto --num 5 DMDA --size 32 64 128 256 512 1024 --comp 1 GMRES/ILU0="snes_monitor par=0.0 ksp_rtol=1.0e-9 mms=2 ksp_type=gmres pc_type=ilu meshExp=2.0 timeExp=-0.5" GMRES/LU="snes_monitor par=0.0 ksp_rtol=1.0e-9 mms=2 ksp_type=gmres pc_type=lu meshExp=2.0 timeExp=-0.75" GMRES/GAMG="snes_monitor par=0.0 ksp_rtol=1.0e-9 mms=2 ksp_type=gmres pc_type=gamg meshExp=2.0 timeExp=-1.0" 552# ./src/benchmarks/benchmarkExample.py --stage "ConvEst Refinement Level %d" --events SNESSolve "ConvEst Error" --plots MeshConvergence WorkPrecision WorkPrecisionPareto --num 13 DMPlex --refine 5 --dim 2 GMRES/ILU0="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=ilu meshExp=2.0 timeExp=-0.75 dm_refine=4 potential_petscspace_order=1" GMRES/LU="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=lu meshExp=2.0 timeExp=-0.75 dm_refine=4 potential_petscspace_order=1" GMRES/GAMG="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=gamg meshExp=2.0 timeExp=-1.0 dm_refine=4 potential_petscspace_order=1" 553# ./src/benchmarks/benchmarkExample.py --stage "ConvEst Refinement Level %d" --events SNESSolve "ConvEst Error" --plots MeshConvergence WorkPrecision WorkPrecisionPareto --num 13 DMPlex --refine 5 --dim 2 GMRES/ILU0="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=ilu meshExp=3.0 timeExp=-0.75 dm_refine=3 potential_petscspace_order=2" GMRES/LU="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=lu meshExp=3.0 timeExp=-0.75 dm_refine=3 potential_petscspace_order=2" GMRES/GAMG="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=gamg meshExp=3.0 timeExp=-1.0 dm_refine=3 potential_petscspace_order=2" 554# ./src/benchmarks/benchmarkExample.py --stage "ConvEst Refinement Level %d" --events SNESSolve "ConvEst Error" --plots MeshConvergence WorkPrecision WorkPrecisionPareto --num 13 DMPlex --refine 5 --dim 2 GMRES/ILU0="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=ilu meshExp=3.0 timeExp=-0.75 dm_refine=3 potential_petscspace_order=2" GMRES/LU="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=lu meshExp=3.0 timeExp=-0.75 dm_refine=3 potential_petscspace_order=2" GMRES/GAMG="snes_monitor ksp_rtol=1.0e-9 ksp_type=gmres pc_type=gamg meshExp=3.0 timeExp=-1.0 dm_refine=3 potential_petscspace_order=2" 555 556# Old GPU benchmarks 557# Benchmark for ex50 558# ./src/benchmarks/benchmarkExample.py --events VecMDot VecMAXPY KSPGMRESOrthog MatMult VecCUSPCopyTo VecCUSPCopyFrom MatCUSPCopyTo --num 50 DMDA --size 10 20 50 100 --comp 4 CPU='pc_type=none mat_no_inode dm_vec_type=seq dm_mat_type=seqaij' GPU='pc_type=none mat_no_inode dm_vec_type=seqcusp dm_mat_type=seqaijcusp cusp_synchronize' 559# Benchmark for ex52 560# ./src/benchmarks/benchmarkExample.py --events IntegBatchCPU IntegBatchGPU IntegGPUOnly --num 52 DMPlex --refine 0.0625 0.00625 0.000625 0.0000625 --blockExp 4 --order=1 CPU='dm_view show_residual=0 compute_function batch' GPU='dm_view show_residual=0 compute_function batch gpu gpu_batches=8' 561# ./src/benchmarks/benchmarkExample.py --events IntegBatchCPU IntegBatchGPU IntegGPUOnly --num 52 DMPlex --refine 0.0625 0.00625 0.000625 0.0000625 --blockExp 4 --order=1 --operator=elasticity CPU='dm_view op_type=elasticity show_residual=0 compute_function batch' GPU='dm_view op_type=elasticity show_residual=0 compute_function batch gpu gpu_batches=8' 562# ./src/benchmarks/benchmarkExample.py --events IntegBatchCPU IntegBatchGPU IntegGPUOnly --num 52 DMPlex --dim=3 --refine 0.0625 0.00625 0.000625 0.0000625 --blockExp 4 --order=1 CPU='dim=3 dm_view show_residual=0 compute_function batch' GPU='dim=3 dm_view show_residual=0 compute_function batch gpu gpu_batches=8' 563