1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68
|
"""
MathWorld: "Hundred-Dollar, Hundred-Digit Challenge Problems", Challenge #3.
http://mathworld.wolfram.com/Hundred-DollarHundred-DigitChallengeProblems.html
The Computer Language Benchmarks Game
http://benchmarksgame.alioth.debian.org/u64q/spectralnorm-description.html#spectralnorm
Contributed by Sebastien Loisel
Fixed by Isaac Gouy
Sped up by Josh Goldfoot
Dirtily sped up by Simon Descarpentries
Concurrency by Jason Stitt
"""
DEFAULT_N = 130
def eval_A(i, j):
return 1.0 / ((i + j) * (i + j + 1) // 2 + i + 1)
def eval_times_u(func, u):
return [func((i, u)) for i in range(len(list(u)))]
def eval_AtA_times_u(u):
return eval_times_u(part_At_times_u, eval_times_u(part_A_times_u, u))
def part_A_times_u(i_u):
i, u = i_u
partial_sum = 0
for j, u_j in enumerate(u):
partial_sum += eval_A(i, j) * u_j
return partial_sum
def part_At_times_u(i_u):
i, u = i_u
partial_sum = 0
for j, u_j in enumerate(u):
partial_sum += eval_A(j, i) * u_j
return partial_sum
def bench_spectral_norm(loops):
range_it = range(loops)
for _ in range_it:
u = [1] * DEFAULT_N
for dummy in range(10):
v = eval_AtA_times_u(u)
u = eval_AtA_times_u(v)
vBv = vv = 0
for ue, ve in zip(u, v):
vBv += ue * ve
vv += ve * ve
def run_benchmark():
bench_spectral_norm(1)
if __name__ == "__main__":
run_benchmark()
|