1import numpy as np 2import time 3 4import benchmark1 5import use_adolc 6import use_scientific 7import use_uncertainties 8import use_algopy 9import use_numdifftools 10import use_funcdesigner 11import use_theano 12 13 14method = {'algopy_reverse_utps':0, 'algopy_reverse_utpm':1, 'pyadolc':2, 'scientific':3, 'uncertainties':4, 'numdifftools':5, 'algopy_forward_utpm':6, 'python':7, 'algopy_forward_utps':8, 'funcdesigner':9, 'theano':10} 15 16 17# FUNCTION COMPUTATION 18# -------------------- 19 20function_N_list = [1,2,4,8,16,32,64,128,256] 21# function_N_list = [2] 22 23 24results_function_list = [] 25for N in function_N_list: 26 print('N=',N) 27 results_function = np.zeros((11,3)) 28 29 # pure python 30 f = benchmark1.F(N) 31 t = time.time(); pass ; preproc_time = time.time() - t 32 t = time.time(); ref_f = f(3*np.ones(N)); run_time = time.time() - t 33 results_function[method['python']] = run_time, abs(ref_f - ref_f), preproc_time 34 print('ref_f=',ref_f) 35 36 # pyadolc 37 f = benchmark1.F(N) 38 t = time.time(); a = use_adolc.EVAL(f, np.ones(N), test='f'); preproc_time = time.time() - t 39 t = time.time(); f = a.function(3*np.ones(N)); run_time = time.time() - t 40 results_function[method['pyadolc']] = run_time, abs(ref_f - f), preproc_time 41 42 # algopy utps version 43 f = benchmark1.F(N) 44 t = time.time(); a = use_algopy.EVAL(f, np.ones(N), test='f'); preproc_time = time.time() - t 45 t = time.time(); f = a.function(3*np.ones(N)); run_time = time.time() - t 46 results_function[method['algopy_forward_utps']] = run_time, abs(ref_f - f), preproc_time 47 48 # algopy utpm version 49 f = benchmark1.G(N) 50 t = time.time(); a = use_algopy.EVAL2(f, np.ones(N), test='f'); preproc_time = time.time() - t 51 t = time.time(); f = a.function(3*np.ones(N)); run_time = time.time() - t 52 results_function[method['algopy_forward_utpm']] = run_time, abs(ref_f - f), preproc_time 53 54 # funcdesigner 55 f = benchmark1.F(N) 56 t = time.time(); a = use_funcdesigner.EVAL(f, np.ones(N), test='f'); preproc_time = time.time() - t 57 t = time.time(); f = a.function(3*np.ones(N)); run_time = time.time() - t 58 results_function[method['funcdesigner']] = run_time, abs(ref_f - f), preproc_time 59 60 # theano 61 f = benchmark1.F(N) 62 t = time.time(); a = use_theano.EVAL(f, np.ones(N), test='f'); preproc_time = time.time() - t 63 t = time.time(); f = a.function(3*np.ones(N)); run_time = time.time() - t 64 results_function[method['theano']] = run_time, np.abs(ref_f - f), preproc_time 65 results_function_list.append(results_function) 66 # print f 67 68results_functions = np.array(results_function_list) 69print('results_functions=\n',results_functions) 70 71# GRADIENT COMPUTATION 72# -------------------- 73 74gradient_N_list = [2,4,8,16,32,64,96] 75# gradient_N_list = [20] 76 77results_gradient_list = [] 78for N in gradient_N_list: 79 print('N=',N) 80 results_gradient = np.zeros((11,3)) 81 82 # pyadolc 83 f = benchmark1.F(N) 84 t = time.time(); a = use_adolc.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 85 t = time.time(); ref_g = a.gradient(3*np.ones(N)); run_time = time.time() - t 86 results_gradient[method['pyadolc']] = run_time, 0, preproc_time 87 88 # scientifc 89 f = benchmark1.F(N) 90 t = time.time(); a = use_scientific.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 91 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 92 results_gradient[method['scientific']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 93 94 # algopy, UTPS variant 95 f = benchmark1.F(N) 96 t = time.time(); a = use_algopy.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 97 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 98 results_gradient[method['algopy_reverse_utps']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 99 100 # algopy, UTPM variant 101 f = benchmark1.G(N) 102 t = time.time(); a = use_algopy.EVAL2(f, np.ones(N), test='g'); preproc_time = time.time() - t 103 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 104 results_gradient[method['algopy_reverse_utpm']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 105 106 # algopy, forward UTPM variant 107 f = benchmark1.G(N) 108 t = time.time(); a = use_algopy.EVAL2(f, np.ones(N), test='fg'); preproc_time = time.time() - t 109 t = time.time(); g = a.forwardgradient(3*np.ones(N)); run_time = time.time() - t 110 results_gradient[method['algopy_forward_utpm']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 111 112 # numdifftools 113 f = benchmark1.F(N) 114 t = time.time(); a = use_numdifftools.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 115 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 116 results_gradient[method['numdifftools']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 117 118 # uncertainties 119 f = benchmark1.F(N) 120 t = time.time(); a = use_uncertainties.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 121 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 122 results_gradient[method['uncertainties']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 123 124 # funcdesigner 125 f = benchmark1.F(N) 126 t = time.time(); a = use_funcdesigner.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 127 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 128 results_gradient[method['funcdesigner']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 129 130 # theano 131 f = benchmark1.F(N) 132 t = time.time(); a = use_theano.EVAL(f, np.ones(N), test='g'); preproc_time = time.time() - t 133 t = time.time(); g = a.gradient(3*np.ones(N)); run_time = time.time() - t 134 results_gradient[method['theano']] = run_time, np.linalg.norm(g - ref_g)/np.linalg.norm(ref_g), preproc_time 135 136 results_gradient_list.append(results_gradient) 137 138results_gradients = np.array(results_gradient_list) 139print('results_gradients=\n',results_gradients) 140 141# HESSIAN COMPUTATION 142# ------------------- 143print('starting hessian computation ') 144results_hessian_list = [] 145hessian_N_list = [1,2,4,8,16,32,64] 146# hessian_N_list = [2] 147 148for N in hessian_N_list: 149 print('N=',N) 150 results_hessian = np.zeros((11,3)) 151 152 # pyadolc 153 f = benchmark1.F(N) 154 t = time.time(); a = use_adolc.EVAL(f, np.ones(N), test='h'); preproc_time = time.time() - t 155 t = time.time(); ref_H = a.hessian(3*np.ones(N)); run_time = time.time() - t 156 results_hessian[method['pyadolc']] = run_time, 0, preproc_time 157 158 # numdifftools 159 f = benchmark1.F(N) 160 t = time.time(); a = use_numdifftools.EVAL(f, np.ones(N), test='h'); preproc_time = time.time() - t 161 t = time.time(); H = a.hessian(3*np.ones(N)); run_time = time.time() - t 162 results_hessian[method['numdifftools']] = run_time, np.linalg.norm( (H-ref_H).ravel())/ np.linalg.norm( (ref_H).ravel()), preproc_time 163 164 # algopy forward utpm variant 165 f = benchmark1.G(N) 166 t = time.time(); a = use_algopy.EVAL2(f, np.ones(N), test='fh'); preproc_time = time.time() - t 167 t = time.time(); H = a.forwardhessian(3*np.ones(N)); run_time = time.time() - t 168 results_hessian[method['algopy_forward_utpm']] = run_time, np.linalg.norm( (H-ref_H).ravel())/ np.linalg.norm( (ref_H).ravel()), preproc_time 169 170 # algopy forward/reverse utpm variant 171 f = benchmark1.G(N) 172 t = time.time(); a = use_algopy.EVAL2(f, np.ones(N), test='h'); preproc_time = time.time() - t 173 t = time.time(); H = a.hessian(3*np.ones(N)); run_time = time.time() - t 174 results_hessian[method['algopy_reverse_utpm']] = run_time, np.linalg.norm( (H-ref_H).ravel())/ np.linalg.norm( (ref_H).ravel()), preproc_time 175 176 # theano 177 f = benchmark1.F(N) 178 t = time.time(); a = use_theano.EVAL(f, np.ones(N), test='h'); preproc_time = time.time() - t 179 t = time.time(); H = a.hessian(3*np.ones(N)); run_time = time.time() - t 180 results_hessian[method['theano']] = run_time, np.linalg.norm( (H-ref_H).ravel())/ np.linalg.norm( (ref_H).ravel()), preproc_time 181 182 results_hessian_list.append(results_hessian) 183 184 185results_hessians = np.array(results_hessian_list) 186 187print(hessian_N_list) 188print('results_hessians=\n',results_hessians) 189 190 191# PLOT RESULTS 192 193print(results_gradients.shape) 194 195import matplotlib.pyplot as pyplot 196import prettyplotting 197 198# plot function run times 199pyplot.figure() 200pyplot.title('Function run times') 201pyplot.loglog(function_N_list, results_functions[:,method['pyadolc'],0], '-ko', markerfacecolor='None', label = 'pyadolc') 202pyplot.loglog(function_N_list, results_functions[:,method['theano'],0], '-.k1', markerfacecolor='None', label = 'theano') 203 204pyplot.loglog(function_N_list, results_functions[:,method['python'],0], '-ks', markerfacecolor='None', label = 'python') 205pyplot.loglog(function_N_list, results_functions[:,method['algopy_forward_utps'],0], '-.k<', markerfacecolor='None', label = 'algopy scalar') 206pyplot.loglog(function_N_list, results_functions[:,method['algopy_forward_utpm'],0], '-.k>', markerfacecolor='None', label = 'algopy matrix') 207pyplot.loglog(function_N_list, results_functions[:,method['funcdesigner'],0], '-.ks', markerfacecolor='None', label = 'funcdesigner') 208pyplot.ylabel('time $t$ [seconds]') 209pyplot.xlabel('problem size $N$') 210pyplot.grid() 211leg = pyplot.legend(loc=2) 212frame= leg.get_frame() 213frame.set_alpha(0.4) 214pyplot.savefig('function_runtimes.png') 215 216# plot gradient run times 217pyplot.figure() 218pyplot.title('Gradient run times') 219pyplot.loglog(gradient_N_list, results_gradients[:,method['pyadolc'],0], '-ko', markerfacecolor='None', label = 'pyadolc') 220pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_reverse_utps'],0], '-.k<', markerfacecolor='None', label = 'algopy reverse utps') 221pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_reverse_utpm'],0], '-.k>', markerfacecolor='None', label = 'algopy reverse utpm') 222pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_forward_utpm'],0], '-.kv', markerfacecolor='None', label = 'algopy forward utpm') 223pyplot.loglog(gradient_N_list, results_gradients[:,method['uncertainties'],0], '--kh', markerfacecolor='None', label = 'uncertainties') 224pyplot.loglog(gradient_N_list, results_gradients[:,method['numdifftools'],0], '--ks', markerfacecolor='None', label = 'numdifftools') 225pyplot.loglog(gradient_N_list, results_gradients[:,method['funcdesigner'],0], '--kd', markerfacecolor='None', label = 'funcdesigner') 226pyplot.loglog(gradient_N_list, results_gradients[:,method['theano'],0], '-.k1', markerfacecolor='None', label = 'theano') 227pyplot.ylabel('time $t$ [seconds]') 228pyplot.xlabel('problem size $N$') 229pyplot.grid() 230leg = pyplot.legend(loc=2) 231frame= leg.get_frame() 232frame.set_alpha(0.4) 233pyplot.savefig('gradient_runtimes.png') 234 235# plot hessian run times 236pyplot.figure() 237pyplot.title('Hessian run times') 238pyplot.loglog(hessian_N_list, results_hessians[:,method['pyadolc'],0], '-ko', markerfacecolor='None', label = 'pyadolc (fo/rev)') 239pyplot.loglog(hessian_N_list, results_hessians[:,method['algopy_forward_utpm'],0], '-.k>', markerfacecolor='None', label = 'algopy (fo)') 240pyplot.loglog(hessian_N_list, results_hessians[:,method['algopy_reverse_utpm'],0], '-.k<', markerfacecolor='None', label = 'algopy (fo/rev)') 241pyplot.loglog(hessian_N_list, results_hessians[:,method['numdifftools'],0], '--ks', markerfacecolor='None', label = 'numdifftools') 242pyplot.loglog(hessian_N_list, results_hessians[:,method['theano'],0], '-.k1', markerfacecolor='None', label = 'theano') 243pyplot.ylabel('time $t$ [seconds]') 244pyplot.xlabel('problem size $N$') 245pyplot.grid() 246leg = pyplot.legend(loc=2) 247frame= leg.get_frame() 248frame.set_alpha(0.4) 249pyplot.savefig('hessian_runtimes.png') 250 251 252# plot gradient preprocessing times 253pyplot.figure() 254pyplot.title('Gradient preprocessing times') 255pyplot.loglog(gradient_N_list, results_gradients[:,method['pyadolc'],2], '-ko', markerfacecolor='None', label = 'pyadolc') 256pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_reverse_utps'],2], '-.k<', markerfacecolor='None', label = 'algopy reverse utps') 257pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_reverse_utpm'],2], '-.k>', markerfacecolor='None', label = 'algopy reverse utpm') 258pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_forward_utpm'],2], '-.kv', markerfacecolor='None', label = 'algopy forward utpm') 259pyplot.loglog(gradient_N_list, results_gradients[:,method['uncertainties'],2], '--kh', markerfacecolor='None', label = 'uncertainties') 260pyplot.loglog(gradient_N_list, results_gradients[:,method['numdifftools'],2], '--ks', markerfacecolor='None', label = 'numdifftools') 261pyplot.loglog(gradient_N_list, results_gradients[:,method['funcdesigner'],2], '--kd', markerfacecolor='None', label = 'funcdesigner') 262pyplot.loglog(gradient_N_list, results_gradients[:,method['theano'],2], '-.k1', markerfacecolor='None', label = 'theano') 263pyplot.ylabel('time $t$ [seconds]') 264pyplot.xlabel('problem size $N$') 265pyplot.grid() 266leg = pyplot.legend(loc=2) 267frame= leg.get_frame() 268frame.set_alpha(0.4) 269pyplot.savefig('gradient_preprocessingtimes.png') 270 271# plot hessian preprocessing times 272pyplot.figure() 273pyplot.title('Hessian preprocessing times') 274pyplot.loglog(hessian_N_list, results_hessians[:,method['pyadolc'],2], '-ko', markerfacecolor='None', label = 'pyadolc (fo/rev)') 275pyplot.loglog(hessian_N_list, results_hessians[:,method['algopy_forward_utpm'],2], '-.k>', markerfacecolor='None', label = 'algopy (fo)') 276pyplot.loglog(hessian_N_list, results_hessians[:,method['algopy_reverse_utpm'],2], '-.k<', markerfacecolor='None', label = 'algopy (fo/rev)') 277pyplot.loglog(hessian_N_list, results_hessians[:,method['numdifftools'],2], '--ks', markerfacecolor='None', label = 'numdifftools') 278pyplot.loglog(hessian_N_list, results_hessians[:,method['theano'],2], '-.k1', markerfacecolor='None', label = 'theano') 279pyplot.ylabel('time $t$ [seconds]') 280pyplot.xlabel('problem size $N$') 281pyplot.grid() 282leg = pyplot.legend(loc=2) 283frame= leg.get_frame() 284frame.set_alpha(0.4) 285pyplot.savefig('hessian_preprocessingtimes.png') 286 287# plot gradient errors 288pyplot.figure() 289pyplot.title('Gradient Correctness') 290pyplot.loglog(gradient_N_list, results_gradients[:,method['numdifftools'],1], '--ks', markerfacecolor='None', label = 'numdifftools') 291pyplot.loglog(gradient_N_list, results_gradients[:,method['pyadolc'],1], '-ko', markerfacecolor='None', label = 'pyadolc') 292pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_reverse_utps'],1], '-.k<', markerfacecolor='None', label = 'algopy reverse utps') 293pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_reverse_utpm'],1], '-.k>', markerfacecolor='None', label = 'algopy reverse utpm') 294pyplot.loglog(gradient_N_list, results_gradients[:,method['algopy_forward_utpm'],1], '-.kv', markerfacecolor='None', label = 'algopy forward utpm') 295pyplot.loglog(gradient_N_list, results_gradients[:,method['uncertainties'],1], '--kh', markerfacecolor='None', label = 'uncertainties') 296pyplot.loglog(gradient_N_list, results_gradients[:,method['funcdesigner'],1], '--kd', markerfacecolor='None', label = 'funcdesigner') 297pyplot.loglog(gradient_N_list, results_gradients[:,method['theano'],1], '-.k1', markerfacecolor='None', label = 'theano') 298pyplot.ylabel(r'relative error $\|g_{ref} - g\|/\|g_{ref}\}$') 299pyplot.xlabel('problem size $N$') 300pyplot.grid() 301leg = pyplot.legend(loc=2) 302frame= leg.get_frame() 303frame.set_alpha(0.4) 304pyplot.savefig('gradient_errors.png') 305 306# plot hessian errors 307pyplot.figure() 308pyplot.title('Hessian Correctness') 309pyplot.loglog(hessian_N_list, results_hessians[:,method['numdifftools'],1], '--ks', markerfacecolor='None', label = 'numdifftools') 310pyplot.loglog(hessian_N_list, results_hessians[:,method['pyadolc'],1], '-ko', markerfacecolor='None', label = 'pyadolc (fo/rev)') 311pyplot.loglog(hessian_N_list, results_hessians[:,method['algopy_forward_utpm'],1], '-.k>', markerfacecolor='None', label = 'algopy (fo)') 312pyplot.loglog(hessian_N_list, results_hessians[:,method['algopy_reverse_utpm'],1], '-.k<', markerfacecolor='None', label = 'algopy (fo/rev)') 313pyplot.loglog(hessian_N_list, results_hessians[:,method['theano'],1], '-.k1', markerfacecolor='None', label = 'theano') 314pyplot.ylabel(r'relative error $\|H_{ref} - H\|/\|H_{ref}\|$') 315pyplot.xlabel('problem size $N$') 316pyplot.grid() 317leg = pyplot.legend(loc=2) 318frame= leg.get_frame() 319frame.set_alpha(0.4) 320pyplot.savefig('hessian_preprocessingtimes.png') 321 322 323 324# pyplot.show() 325 326 327