Addition-chain exponentiation: Difference between revisions
Content added Content deleted
m (typo) |
m (Python example) |
||
Line 1,751: | Line 1,751: | ||
--1.00003 ^ 27182 (20) = 1.9999999999727968238298855 |
--1.00003 ^ 27182 (20) = 1.9999999999727968238298855 |
||
--1.00003 ^ 27182 (19) = 1.9999999999727968238298527</pre> |
--1.00003 ^ 27182 (19) = 1.9999999999727968238298527</pre> |
||
=={{header|Python}}== |
|||
<syntaxhighlight lang="python">''' Rosetta Code task Addition-chain_exponentiation ''' |
|||
from math import sqrt |
|||
from numpy import array |
|||
from mpmath import mpf |
|||
class AdditionChains: |
|||
''' two methods of calculating addition chains ''' |
|||
def __init__(self): |
|||
''' memoization for knuth_path ''' |
|||
self.chains, self.idx, self.pos = [[1]], 0, 0 |
|||
self.pat, self.lvl = {1: 0}, [[1]] |
|||
def add_chain(self): |
|||
''' method 1: add chains depth then breadth first until done ''' |
|||
newchain = self.chains[self.idx].copy() |
|||
newchain.append(self.chains[self.idx][-1] + |
|||
self.chains[self.idx][self.pos]) |
|||
self.chains.append(newchain) |
|||
if self.pos == len(self.chains[self.idx])-1: |
|||
self.idx += 1 |
|||
self.pos = 0 |
|||
else: |
|||
self.pos += 1 |
|||
return newchain |
|||
def find_chain(self, nexp): |
|||
''' method 1 interface: search for chain ending with n, adding more as needed ''' |
|||
assert nexp > 0 |
|||
if nexp == 1: |
|||
return [1] |
|||
chn = next((a for a in self.chains if a[-1] == nexp), None) |
|||
if chn is None: |
|||
while True: |
|||
chn = self.add_chain() |
|||
if chn[-1] == nexp: |
|||
break |
|||
return chn |
|||
def knuth_path(self, ngoal): |
|||
''' method 2: knuth method, uses memoization to search for a shorter chain ''' |
|||
if ngoal < 1: |
|||
return [] |
|||
while not ngoal in self.pat: |
|||
new_lvl = [] |
|||
for i in self.lvl[0]: |
|||
for j in self.knuth_path(i): |
|||
if not i + j in self.pat: |
|||
self.pat[i + j] = i |
|||
new_lvl.append(i + j) |
|||
self.lvl[0] = new_lvl |
|||
returnpath = self.knuth_path(self.pat[ngoal]) |
|||
returnpath.append(ngoal) |
|||
return returnpath |
|||
def cpow(xbase, chain): |
|||
''' raise xbase by an addition exponentiation chain for what becomes x**chain[-1] ''' |
|||
pows, products = 0, {0: 1, 1: xbase} |
|||
for i in chain: |
|||
products[i] = products[pows] * products[i - pows] |
|||
pows = i |
|||
return products[chain[-1]] |
|||
if __name__ == '__main__': |
|||
# test both addition chain methods |
|||
acs = AdditionChains() |
|||
print('First one hundred addition chain lengths:') |
|||
for k in range(1, 101): |
|||
print(f'{len(acs.find_chain(k))-1:3}', end='\n'if k % 10 == 0 else '') |
|||
print('\nKnuth chains for addition chains of 31415 and 27182:') |
|||
chns = {m: acs.knuth_path(m) for m in [31415, 27182]} |
|||
for (num, cha) in chns.items(): |
|||
print(f'Exponent: {num:10}\n Addition Chain: {cha[:-1]}') |
|||
print('\n1.00002206445416^31415 =', cpow(1.00002206445416, chns[31415])) |
|||
print('1.00002550055251^27182 =', cpow(1.00002550055251, chns[27182])) |
|||
print('1.000025 + 0.000058i)^27182 =', |
|||
cpow(complex(1.000025, 0.000058), chns[27182])) |
|||
print('1.000022 + 0.000050i)^31415 =', |
|||
cpow(complex(1.000022, 0.000050), chns[31415])) |
|||
sq05 = mpf(sqrt(0.5)) |
|||
mat = array([[sq05, 0, sq05, 0, 0, 0], [0, sq05, 0, sq05, 0, 0], [0, sq05, 0, -sq05, 0, 0], |
|||
[-sq05, 0, sq05, 0, 0, 0], [0, 0, 0, 0, 0, 1], [0, 0, 0, 0, 1, 0]]) |
|||
print('matrix A ^ 27182 =') |
|||
print(cpow(mat, chns[27182])) |
|||
print('matrix A ^ 31415 =') |
|||
print(cpow(mat, chns[31415])) |
|||
print('(matrix A ** 27182) ** 31415 =') |
|||
print(cpow(cpow(mat, chns[27182]), chns[31415])) |
|||
</syntaxhighlight>{{out}} |
|||
<pre> |
|||
First one hundred addition chain lengths: |
|||
0 1 2 2 3 3 4 3 4 4 |
|||
5 4 5 5 5 4 5 5 6 5 |
|||
6 6 6 5 6 6 6 6 7 6 |
|||
7 5 6 6 7 6 7 7 7 6 |
|||
7 7 7 7 7 7 8 6 7 7 |
|||
7 7 8 7 8 7 8 8 8 7 |
|||
8 8 8 6 7 7 8 7 8 8 |
|||
9 7 8 8 8 8 8 8 9 7 |
|||
8 8 8 8 8 8 9 8 9 8 |
|||
9 8 9 9 9 7 8 8 8 8 |
|||
Knuth chains for addition chains of 31415 and 27182: |
|||
Exponent: 31415 |
|||
Addition Chain: [1, 2, 3, 5, 7, 14, 28, 56, 61, 122, 244, 488, 976, 1952, 3904, 7808, 15616, 15677, 31293] |
|||
Exponent: 27182 |
|||
Addition Chain: [1, 2, 3, 5, 7, 14, 21, 35, 70, 140, 143, 283, 566, 849, 1698, 3396, 6792, 6799, 13591] |
|||
1.00002206445416^31415 = 1.9999999998924638 |
|||
1.00002550055251^27182 = 1.9999999999792688 |
|||
1.000025 + 0.000058i)^27182 = (-0.01128636963542673+1.9730308496660347j) |
|||
1.000022 + 0.000050i)^31415 = (0.00016144681325535107+1.9960329014194498j) |
|||
matrix A ^ 27182 = |
|||
[[mpf('5.0272320351359433e-4092') 0 mpf('5.0272320351359433e-4092') 0 0 0] |
|||
[0 mpf('5.0272320351359433e-4092') 0 mpf('5.0272320351359433e-4092') 0 0] |
|||
[0 mpf('5.0272320351359433e-4092') 0 mpf('5.0272320351359433e-4092') 0 0] |
|||
[mpf('5.0272320351359433e-4092') 0 mpf('5.0272320351359433e-4092') 0 0 0] |
|||
[0 0 0 0 0 1] |
|||
[0 0 0 0 1 0]] |
|||
matrix A ^ 31415 = |
|||
[[mpf('3.7268602513250562e-4729') 0 mpf('3.7268602513250562e-4729') 0 0 0] |
|||
[0 mpf('3.7268602513250562e-4729') 0 mpf('3.7268602513250562e-4729') 0 0] |
|||
[0 mpf('3.7268602513250562e-4729') 0 mpf('-3.7268602513250562e-4729') 0 |
|||
0] |
|||
[mpf('-3.7268602513250562e-4729') 0 mpf('3.7268602513250562e-4729') 0 0 |
|||
0] |
|||
[0 0 0 0 0 1] |
|||
[0 0 0 0 1 0]] |
|||
(matrix A ** 27182) ** 31415 = |
|||
[[mpf('1.77158544475749e-128528148') 0 mpf('1.77158544475749e-128528148') |
|||
0 0 0] |
|||
[0 mpf('1.77158544475749e-128528148') 0 |
|||
mpf('1.77158544475749e-128528148') 0 0] |
|||
[0 mpf('1.77158544475749e-128528148') 0 |
|||
mpf('1.77158544475749e-128528148') 0 0] |
|||
[mpf('1.77158544475749e-128528148') 0 mpf('1.77158544475749e-128528148') |
|||
0 0 0] |
|||
[0 0 0 0 0 1] |
|||
[0 0 0 0 1 0]] |
|||
</pre> |
|||
=={{header|Racket}}== |
=={{header|Racket}}== |