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
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
|
#!/usr/bin/env python3
#
# Copyright 2022 Google LLC
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
#
import numpy as np
import matplotlib.pyplot as plt
def fast_exp2(x, p):
p = p.astype(np.float32)
x = x.astype(np.float32)
y = (((((p[0]*x) + p[1])*x + p[2])*x + p[3])*x + p[4])*x + 1
return np.power(y.astype(np.float32), 16)
def approx_exp2():
x = np.arange(-8, 8, step=1e-3)
p = np.polyfit(x, ((2 ** (x/16)) - 1) / x, 4)
y = [ fast_exp2(x[i], p) for i in range(len(x)) ]
e = np.abs(y - 2**x) / (2 ** x)
print('{{ {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e} }}'
.format(p[0], p[1], p[2], p[3], p[4]))
print('Max relative error: ', np.max(e))
print('Max RMS error: ', np.sqrt(np.mean(e ** 2)))
if False:
fig, (ax1, ax2) = plt.subplots(2)
ax1.plot(x, 2**x, label='Reference')
ax1.plot(x, y, label='Approximation')
ax1.legend()
ax2.plot(x, e, label='Relative Error')
ax2.legend()
plt.show()
def fast_log2(x, p):
p = p.astype(np.float32)
x = x.astype(np.float32)
(x, e) = np.frexp(x)
y = ((((p[0]*x) + p[1])*x + p[2])*x + p[3])*x + p[4]
return (e ) + y.astype(np.float32)
def approx_log2():
x = np.logspace(-1, 0, base=2, num=100)
p = np.polyfit(x, np.log2(x), 4)
x = np.logspace(-2, 5, num=10000)
y = [ fast_log2(x[i], p) for i in range(len(x)) ]
e = np.abs(y - np.log2(x))
print('{{ {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e} }}'
.format(p[0], p[1], p[2], p[3], p[4]))
print('Max absolute error: ', np.max(e))
print('Max RMS error: ', np.sqrt(np.mean(e ** 2)))
if False:
fig, (ax1, ax2) = plt.subplots(2)
ax1.plot(x, np.log2(x), label='Reference')
ax1.plot(x, y, label='Approximation')
ax1.legend()
ax2.plot(x, e, label = 'Absolute error')
ax2.legend()
plt.show()
def table_db_q16():
k = 10 * np.log10(2);
for i in range(32):
a = k * np.log2(np.ldexp(32 + i , -5)) - (i // 16) * (k/2);
b = k * np.log2(np.ldexp(32 + i+1, -5)) - (i // 16) * (k/2);
an = np.ldexp(a, 15) + 0.5
bn = np.ldexp(b - a, 15) + 0.5
print('{{ {:5d}, {:4d} }},'
.format(int(np.ldexp(a, 15) + 0.5),
int(np.ldexp(b - a, 15) + 0.5)),
end = ' ' if i % 4 < 3 else '\n')
if __name__ == '__main__':
print('\n--- Approximation of 2^n ---')
approx_exp2()
print('\n--- Approximation of log2(n) ---')
approx_log2()
print('\n--- Table of fixed Q16 dB ---')
table_db_q16()
print('')
|
