1*49fe348cSAndroid Build Coastguard Worker#!/usr/bin/env python3 2*49fe348cSAndroid Build Coastguard Worker# 3*49fe348cSAndroid Build Coastguard Worker# Copyright 2022 Google LLC 4*49fe348cSAndroid Build Coastguard Worker# 5*49fe348cSAndroid Build Coastguard Worker# Licensed under the Apache License, Version 2.0 (the "License"); 6*49fe348cSAndroid Build Coastguard Worker# you may not use this file except in compliance with the License. 7*49fe348cSAndroid Build Coastguard Worker# You may obtain a copy of the License at 8*49fe348cSAndroid Build Coastguard Worker# 9*49fe348cSAndroid Build Coastguard Worker# http://www.apache.org/licenses/LICENSE-2.0 10*49fe348cSAndroid Build Coastguard Worker# 11*49fe348cSAndroid Build Coastguard Worker# Unless required by applicable law or agreed to in writing, software 12*49fe348cSAndroid Build Coastguard Worker# distributed under the License is distributed on an "AS IS" BASIS, 13*49fe348cSAndroid Build Coastguard Worker# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14*49fe348cSAndroid Build Coastguard Worker# See the License for the specific language governing permissions and 15*49fe348cSAndroid Build Coastguard Worker# limitations under the License. 16*49fe348cSAndroid Build Coastguard Worker# 17*49fe348cSAndroid Build Coastguard Worker 18*49fe348cSAndroid Build Coastguard Workerimport numpy as np 19*49fe348cSAndroid Build Coastguard Workerimport matplotlib.pyplot as plt 20*49fe348cSAndroid Build Coastguard Worker 21*49fe348cSAndroid Build Coastguard Worker 22*49fe348cSAndroid Build Coastguard Workerdef fast_exp2(x, t, p): 23*49fe348cSAndroid Build Coastguard Worker 24*49fe348cSAndroid Build Coastguard Worker p = p.astype(np.float32) 25*49fe348cSAndroid Build Coastguard Worker x = x.astype(np.float32) 26*49fe348cSAndroid Build Coastguard Worker 27*49fe348cSAndroid Build Coastguard Worker m = ((x + 0.5/8) % (1/8)) - (0.5/8) 28*49fe348cSAndroid Build Coastguard Worker e = int((x - m) * 8) 29*49fe348cSAndroid Build Coastguard Worker 30*49fe348cSAndroid Build Coastguard Worker y = ((((p[0]*m) + p[1])*m + p[2])*m + p[3])*m + p[4] 31*49fe348cSAndroid Build Coastguard Worker y = y * 2**(e // 8) * t[e % 8] 32*49fe348cSAndroid Build Coastguard Worker 33*49fe348cSAndroid Build Coastguard Worker return y.astype(np.float32) 34*49fe348cSAndroid Build Coastguard Worker 35*49fe348cSAndroid Build Coastguard Workerdef approx_exp2(): 36*49fe348cSAndroid Build Coastguard Worker 37*49fe348cSAndroid Build Coastguard Worker x = np.arange(0, 1/8, step=1e-6) 38*49fe348cSAndroid Build Coastguard Worker p = np.polyfit(x, 2 ** x, 4) 39*49fe348cSAndroid Build Coastguard Worker t = [ 2**(i/8) for i in range(8) ] 40*49fe348cSAndroid Build Coastguard Worker 41*49fe348cSAndroid Build Coastguard Worker x = np.arange(-10, 10, step=1e-3) 42*49fe348cSAndroid Build Coastguard Worker y = [ fast_exp2(x[i], t, p) for i in range(len(x)) ] 43*49fe348cSAndroid Build Coastguard Worker 44*49fe348cSAndroid Build Coastguard Worker e = np.abs(y - 2**x) / (2 ** x) 45*49fe348cSAndroid Build Coastguard Worker 46*49fe348cSAndroid Build Coastguard Worker print('{{ {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e}, \n' 47*49fe348cSAndroid Build Coastguard Worker ' {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e}, '.format(*t)) 48*49fe348cSAndroid Build Coastguard Worker print('{{ {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e} }}'.format(*p)) 49*49fe348cSAndroid Build Coastguard Worker print('Max relative error: ', np.max(e)) 50*49fe348cSAndroid Build Coastguard Worker print('Max RMS error: ', np.sqrt(np.mean(e ** 2))) 51*49fe348cSAndroid Build Coastguard Worker 52*49fe348cSAndroid Build Coastguard Worker if False: 53*49fe348cSAndroid Build Coastguard Worker fig, (ax1, ax2) = plt.subplots(2) 54*49fe348cSAndroid Build Coastguard Worker 55*49fe348cSAndroid Build Coastguard Worker ax1.plot(x, 2**x, label='Reference') 56*49fe348cSAndroid Build Coastguard Worker ax1.plot(x, y, label='Approximation') 57*49fe348cSAndroid Build Coastguard Worker ax1.legend() 58*49fe348cSAndroid Build Coastguard Worker 59*49fe348cSAndroid Build Coastguard Worker ax2.plot(x, e, label='Relative Error') 60*49fe348cSAndroid Build Coastguard Worker ax2.legend() 61*49fe348cSAndroid Build Coastguard Worker 62*49fe348cSAndroid Build Coastguard Worker plt.show() 63*49fe348cSAndroid Build Coastguard Worker 64*49fe348cSAndroid Build Coastguard Worker 65*49fe348cSAndroid Build Coastguard Workerdef fast_log2(x, p): 66*49fe348cSAndroid Build Coastguard Worker 67*49fe348cSAndroid Build Coastguard Worker p = p.astype(np.float32) 68*49fe348cSAndroid Build Coastguard Worker x = x.astype(np.float32) 69*49fe348cSAndroid Build Coastguard Worker 70*49fe348cSAndroid Build Coastguard Worker (x, e) = np.frexp(x) 71*49fe348cSAndroid Build Coastguard Worker 72*49fe348cSAndroid Build Coastguard Worker y = ((((p[0]*x) + p[1])*x + p[2])*x + p[3])*x + p[4] 73*49fe348cSAndroid Build Coastguard Worker 74*49fe348cSAndroid Build Coastguard Worker return (e ) + y.astype(np.float32) 75*49fe348cSAndroid Build Coastguard Worker 76*49fe348cSAndroid Build Coastguard Workerdef approx_log2(): 77*49fe348cSAndroid Build Coastguard Worker 78*49fe348cSAndroid Build Coastguard Worker x = np.logspace(-1, 0, base=2, num=100) 79*49fe348cSAndroid Build Coastguard Worker p = np.polyfit(x, np.log2(x), 4) 80*49fe348cSAndroid Build Coastguard Worker 81*49fe348cSAndroid Build Coastguard Worker x = np.logspace(-2, 5, num=10000) 82*49fe348cSAndroid Build Coastguard Worker y = [ fast_log2(x[i], p) for i in range(len(x)) ] 83*49fe348cSAndroid Build Coastguard Worker e = np.abs(y - np.log2(x)) 84*49fe348cSAndroid Build Coastguard Worker 85*49fe348cSAndroid Build Coastguard Worker print('{{ {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e}, {:14.8e} }}' 86*49fe348cSAndroid Build Coastguard Worker .format(p[0], p[1], p[2], p[3], p[4])) 87*49fe348cSAndroid Build Coastguard Worker print('Max absolute error: ', np.max(e)) 88*49fe348cSAndroid Build Coastguard Worker print('Max RMS error: ', np.sqrt(np.mean(e ** 2))) 89*49fe348cSAndroid Build Coastguard Worker 90*49fe348cSAndroid Build Coastguard Worker if False: 91*49fe348cSAndroid Build Coastguard Worker fig, (ax1, ax2) = plt.subplots(2) 92*49fe348cSAndroid Build Coastguard Worker 93*49fe348cSAndroid Build Coastguard Worker ax1.plot(x, np.log2(x), label='Reference') 94*49fe348cSAndroid Build Coastguard Worker ax1.plot(x, y, label='Approximation') 95*49fe348cSAndroid Build Coastguard Worker ax1.legend() 96*49fe348cSAndroid Build Coastguard Worker 97*49fe348cSAndroid Build Coastguard Worker ax2.plot(x, e, label = 'Absolute error') 98*49fe348cSAndroid Build Coastguard Worker ax2.legend() 99*49fe348cSAndroid Build Coastguard Worker 100*49fe348cSAndroid Build Coastguard Worker plt.show() 101*49fe348cSAndroid Build Coastguard Worker 102*49fe348cSAndroid Build Coastguard Worker 103*49fe348cSAndroid Build Coastguard Workerdef table_db_q16(): 104*49fe348cSAndroid Build Coastguard Worker 105*49fe348cSAndroid Build Coastguard Worker k = 10 * np.log10(2); 106*49fe348cSAndroid Build Coastguard Worker 107*49fe348cSAndroid Build Coastguard Worker for i in range(32): 108*49fe348cSAndroid Build Coastguard Worker a = k * np.log2(np.ldexp(32 + i , -5)) - (i // 16) * (k/2); 109*49fe348cSAndroid Build Coastguard Worker b = k * np.log2(np.ldexp(32 + i+1, -5)) - (i // 16) * (k/2); 110*49fe348cSAndroid Build Coastguard Worker 111*49fe348cSAndroid Build Coastguard Worker an = np.ldexp(a, 15) + 0.5 112*49fe348cSAndroid Build Coastguard Worker bn = np.ldexp(b - a, 15) + 0.5 113*49fe348cSAndroid Build Coastguard Worker print('{{ {:5d}, {:4d} }},' 114*49fe348cSAndroid Build Coastguard Worker .format(int(np.ldexp(a, 15) + 0.5), 115*49fe348cSAndroid Build Coastguard Worker int(np.ldexp(b - a, 15) + 0.5)), 116*49fe348cSAndroid Build Coastguard Worker end = ' ' if i % 4 < 3 else '\n') 117*49fe348cSAndroid Build Coastguard Worker 118*49fe348cSAndroid Build Coastguard Worker 119*49fe348cSAndroid Build Coastguard Workerif __name__ == '__main__': 120*49fe348cSAndroid Build Coastguard Worker 121*49fe348cSAndroid Build Coastguard Worker print('\n--- Approximation of 2^n ---') 122*49fe348cSAndroid Build Coastguard Worker approx_exp2() 123*49fe348cSAndroid Build Coastguard Worker 124*49fe348cSAndroid Build Coastguard Worker print('\n--- Approximation of log2(n) ---') 125*49fe348cSAndroid Build Coastguard Worker approx_log2() 126*49fe348cSAndroid Build Coastguard Worker 127*49fe348cSAndroid Build Coastguard Worker print('\n--- Table of fixed Q16 dB ---') 128*49fe348cSAndroid Build Coastguard Worker table_db_q16() 129*49fe348cSAndroid Build Coastguard Worker 130*49fe348cSAndroid Build Coastguard Worker print('') 131