## "Fossies" - the Fresh Open Source Software Archive

### Member "NZMATH-1.2.0/nzmath/algorithm.py" (19 Nov 2012, 6063 Bytes) of package /linux/misc/old/NZMATH-1.2.0.tar.gz:

As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) Python source code syntax highlighting (style: standard) with prefixed line numbers. Alternatively you can here view or download the uninterpreted source code file. For more information about "algorithm.py" see the Fossies "Dox" file reference documentation.

```    1 """
2 Abstract (higher order) functions
3 """
4
5 import functools
6
7 def digital_method(coefficients, val, add, mul, act, power, zero, one):
8     """
9     Evaluate a univariate polynomial at val.
10
11     The polynomial is given as an iterator 'coefficients' which yields
12     (degree, coefficient) pair in descending order of degree.
13
14     If the polynomial is of R-coefficients, then val should be in an
15     R-algebra D.
16
17     All operations 'add', 'mul', 'act', 'power', 'zero', 'one' should
18     be explicitly given, where:
20       'mul' multiplication (D x D -> D),
21       'act' action of R (R x D -> D),
22       'power' powering (D x Z -> D) and
23       'zero', 'one' constants in D.
24     """
25     result = zero
26     val_pow = {0:one, 1:val}
27     d = None
28     for d_next, c_next in coefficients:
29         if d:
30             diff = d - d_next
31             if diff not in val_pow:
32                 val_pow[diff] = power(val, diff)
33             result = add(mul(result, val_pow[diff]), act(c_next, one))
34         else:
35             result = act(c_next, one)
36         d = d_next
37     if d:
38         return mul(result, val_pow.get(d, power(val, d)))
39     else:
40         return result
41
42 def digital_method_func(add, mul, act, power, zero, one):
43     """
44     Return a function which evaluates polynomial (binary function)
45     constructed from given functions.
46     """
48                                act=act, power=power, zero=zero, one=one)
49     result.___doc__ = """
50     Evaluate a univariate polynomial at val.
51
52     The polynomial is given as an iterator 'coefficients' which yields
53     (degree, coefficient) pair in descending order of degree.
54
55     If the polynomial is of R-coefficients, then val should be in an
56     R-algebra D.
57     """
58     return result
59
60 def rl_binary_powering(element, index, mul, square=None, one=None):
61     """
62      powering by using right-left binary method
63      (assume that index >= 0)
64     """
65     square = _set_square(mul, square)
66     sol = one
67     mul_part = element
68     while True:
69         if index & 1:
70             try:
71                 sol = mul(sol, mul_part)
72             except: # one is None!
73                 sol = mul_part
74         index >>= 1
75         if not index:
76             return sol
77         else:
78             mul_part = square(mul_part)
79
80 def lr_binary_powering(element, index, mul, square=None, one=None):
81     """
82     powering by using left-right binary method.
83     (assume that index >= 0)
84     """
85     import math
86     square = _set_square(mul, square)
87     if not index:
88         return one
89     spot = 1 << long(math.log(index, 2))
90     sol = one
91     while True:
92         if spot & index:
93             try:
94                 sol = mul(sol, element)
95             except: # one is None!
96                 sol = element
97         spot >>= 1
98         if not spot:
99             return sol
100         sol = square(sol)
101
102 def window_powering(element, index, mul, square=None, one=None):
103     """
104     powering by using small-window method
105     window size is selected by average analystic optimization
106     (assume that index >= 0)
107     """
108     import math
109     square = _set_square(mul, square)
110     if not index:
111         return one
112     log_n = long(math.log(index, 2))
113     # Find the proper window size
114     size = 2
115     pow_size = 2
116     while log_n > (size + 1) * (size + 2) * (pow_size - 1):
117         pow_size <<= 1
118         size += 1
119     # Precomputation
120     sqr = square(element)
121     pre_table = [element]
122     pow_size -= 1
123     for i in range(pow_size):
124         pre_table.append(mul(pre_table[-1], sqr))
125
126     near_n = 1 << log_n
127     spot = near_n
128     while spot:
129         if not(spot & index):
130             sol = square(sol)
131             spot >>= 1
132         else:
133             # Find the window
134             f_spot, e_spot = spot, spot >> (size - 1)
135             t_size = size
136             if not(e_spot):
137                 e_spot = 1
138                 t_size = int(math.log(f_spot, 2)) + 1
139             while True:
140                 if e_spot & index:
141                     spot = (e_spot >> 1)
142                     window = 0
143                     sqr = 1
144                     while f_spot != e_spot: # Compute value of window
145                         if index & e_spot:
146                             window += sqr
147                         sqr <<= 1
148                         e_spot <<= 1
149                     window += sqr
150                     break
151                 t_size -= 1
152                 e_spot <<= 1
153             # Compute a part of powering
154             if f_spot == near_n:
155                 sol = pre_table[(window - 1) >> 1]
156             else:
157                 for i in range(t_size):
158                     sol = square(sol)
159                 sol = mul(sol, pre_table[(window - 1) >> 1])
160     return sol
161
162 def _set_square(mul, square):
163     """
164     Set square if square is None
165     (for the initialization of methods for powering)
166     """
167     if square == None:
168         return lambda x : mul(x, x)
169     else:
170         return square
171
172 def powering_func(mul, square=None, one=None, type=0):
173     """
174     Return a function which computes powering (binary function)
175     constructed from given functions.
176     """
177     if type == 0:
178         result = functools.partial(
179             rl_binary_powering, mul=mul, one=one, square=square)
180         result.___doc__ = """
181         compute powering, that is, element ** index.
182
183         This function uses right-left binary method.
184         (assume that index >= 0)
185         """
186     elif type == 1:
187         result = functools.partial(
188             lr_binary_powering, mul=mul, one=one, square=square)
189         result.___doc__ = """
190         compute powering, that is, element ** index.
191
192         This function uses left-right binary method.
193         (assume that index >= 0)
194         """
195     elif type == 2:
196         result = functools.partial(
197             window_powering, mul=mul, one=one, square=square)
198         result.___doc__ = """
199         compute powering, that is, element ** index.
200
201         This function uses window method.
202         (assume that index >= 0)
203         """
204     return result
```