#numerical library
import numpy as np
#plot library
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
from matplotlib.ticker import LinearLocator, FormatStrFormatter
import matplotlib.pyplot as plt
from pprint import pprint

"""
Explanation: Norm approximation on restricted quantized domain
Fast approximation of the norm over the value of a 10bit unsigned int 
using batch gradient descent to minimize the squared error
Approximation formula:
$$\gamma (x+y) \approx \sqrt{x^2+y^2} \ x,y \in [0,1023] \subseteq \mathbb{N}$$
where $$\gamma = 0.7531854654594905$$
import needed modules
End of explanation
"""


x,y = np.indices([1024,1024])

"""
Explanation: inizialization of the variable
End of explanation
"""


%timeit (x**2+y**2)**0.5

%timeit 0.7531854654594905*(x+y)

"""
Explanation: Guessing the speed gain if this work
End of explanation
"""


epoch_of_training = 1000
learning_rate = 1e-8
gamma = 1

"""
Explanation: So if this approximation work we can stimate the norm over x7 time faster

hyper parameters initialization
End of explanation
"""


init_sq_err = np.sum(0.5*((gamma*(x+y) - (x**2+y**2)**0.5)**2))
init_sq_err

"""
Explanation: Initial Error
Error formula:
$$E = \frac{1}{2}\sum_{1}^{n} (\gamma (x+y)-\sqrt{x^2+y^2})^2$$
End of explanation
"""


for i in range(epoch_of_training):
    gamma -= learning_rate * np.mean((gamma*(x+y)-(x**2+y**2)**0.5)*(x+y))

"""
Explanation: the partial derivative of the error
Partial derivative formula:
$$\frac{\partial }{\partial \gamma}E =\sum_{1}^{n} (\gamma (x+y)-\sqrt{x^2+y^2})(x+y)$$
The stocastic grdient descent
$$\gamma^{i} = \gamma^{i-1} - \eta \frac{\partial }{\partial \gamma}E$$
End of explanation
"""


gamma

fin_sq_err = np.sum(0.5*((gamma*(x+y) - (x**2+y**2)**0.5)**2))
fin_sq_err

delta_sq_err = init_sq_err - fin_sq_err
delta_sq_err

Error = abs(gamma*(x+y) - (x**2+y**2)**0.5)
print(np.max(Error))
print(np.mean(Error))
print(np.min(Error))

"""
Explanation: Results
End of explanation
"""


fig = plt.figure(figsize=[12,12])
ax = fig.gca(projection='3d')

X = np.arange(1, 1024, 8)
Y = np.arange(1, 1024, 8)
X, Y = np.meshgrid(X, Y)
Z = np.sqrt(X**2 + Y**2)
F = abs(gamma*(X+Y)-Z)/(Z)

surf = ax.plot_surface(X, Y, F, rstride=2, cstride=2, cmap=cm.jet,linewidth=1)

ax.set_xlabel('X')
ax.set_xlim(-10, 1034)
ax.set_ylabel('Y')
ax.set_ylim(-10, 1034)
ax.set_zlabel('Z')
ax.set_zlim(0, 0.30)

ax.invert_yaxis()

ax.zaxis.set_major_locator(LinearLocator(10))
ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))

fig.colorbar(surf)

plt.show()


"""
Explanation: Percentual Error Plot
End of explanation
"""


