Dual-numbers, auto-differentiation and gradient descent

Gem City Tech AI/ML

Start at the end - Gradient descent

Algorithm to find the minimum of a function

\[ f(x) = x^2 \]

Can be written as python code:

import numpy as np

def f(x):
    return x**2

And plotted:

from matplotlib import pylab as plt

fig,ax = plt.subplots()
x = np.linspace(-10,10,50)
ax.plot(x,f(x))
ax.grid()

Figure 1: \(f(x) = x^2\)

Standard solution using calculus

Minimum \(x_{min}\) occurs where derivative is zero

\[ \left. \frac{d f(x)}{d x}\right\rvert_{x_{min}} = 0 \]

In the example

\[ \frac{d f(x)}{d x} = 2x \]

and

\[ 2 x = 0 \]

when

\[ x=0 \]

What if minimum is hard to find?

• Maybe so complicated that we don’t know how to compute it

Solution:

• create an algorithm to find it - gradient descent

\[ f(x_{k+1}) = f(x_k) - \eta \left. \frac{d f(x)}{d x} \right\rvert_{x_k} \]

Gradient descent has this parameter \(\eta\) which is called the “learning rate” in the machine learning community.

Simple gradient descent algorithm in python

import numpy as np

def gradient_descent(
        x0,              # starting point
        dfdx,            # gradient of f(x)
        eta,             # learning rate
        tolerance = 0.01, # termination condition
        max_steps = 100, # maximum steps
        ):
    
    x = x0 # initialize
    xs = [x0] # list of steps

    for _ in range(max_steps):
        delta = eta * dfdx(x)
        if np.abs(delta) < tolerance:
            break
        x = x - delta
        xs.append(x)

    return x, xs

Gradient descent simple example

Example function from [1]

\[ \begin{align} f(x) &= x^2 - 4 x + 1 \\ \frac{df(x)}{dx} &= 2x - 4 \end{align} \]

Python function definitions for \(f(x)\) and \(\frac{df(x)}{dx}\)

def func(x):
    return x**2 - 4*x + 1

def dfunc(x):
    return 2*x - 4

%%capture
from IPython.display import HTML
import matplotlib.animation as animation

xres,xhist = gradient_descent(8,dfunc,0.1,)

fig, ax = plt.subplots()
# plot underlying function
xvals = np.linspace(-6,10,100)
ax.plot(xvals,func(xvals),'b-')
ax.grid()

def animate(ii,):
    xh = np.array(xhist[0:ii])
    ax.plot(xh,func(xh),'ro-')
    for ii,_x in enumerate(xh):
        ax.annotate(f'{ii}', (_x+0.4, func(_x)-1.5,),color='r')

ani = animation.FuncAnimation(fig,animate,repeat=False,frames=len(xhist)-1,interval=1000)

HTML(ani.to_jshtml())

Once Loop Reflect

Figure 2: Gradient descent for \(\eta=0.1\)

%%capture
from IPython.display import HTML
import matplotlib.animation as animation

xres,xhist = gradient_descent(8,dfunc,0.8,)

fig, ax = plt.subplots()
# plot underlying function
xvals = np.linspace(-6,10,100)
ax.plot(xvals,func(xvals),'b-')
ax.grid()

def animate(ii,):
    xh = np.array(xhist[0:ii])
    ax.plot(xh,func(xh),'ro-')
    for ii,_x in enumerate(xh):
        ax.annotate(f'{ii}', (_x+0.4, func(_x)-1.5,),color='r')

ani = animation.FuncAnimation(fig,animate,repeat=False,frames=len(xhist)-1,interval=1000)

HTML(ani.to_jshtml())

Once Loop Reflect

Figure 3: Gradient descent for \(\eta=0.8\)

Issue

This gradient descent algorithm needs to know \(d f(x) \over d x\). What if \(f(x)\) is so complicated we don’t know how to compute \(d f(x) \over d x\)?

Solution is auto-differentiation.

Auto-differentiation

• How do we automatically compute the derivative?
• With dual-numbers

Dual numbers

Imaginary numbers

You might be familiar with imaginary and complex numbers where \(i\) or \(j\) is the imaginary number. I’m a physicist, so I use \(i\)

\[ i = \sqrt{-1} \]

and the complex number \(z\)

\[ z = x + i y \]

Dual numbers

Dual numbers are similar to imaginary numbers and use the special notation \(\epsilon\) with the properties

\[ \epsilon \ne 0 \]

and

\[ \epsilon^2 = 0 \]

Similar to complex numbers we can construct

\[ z = x + \epsilon y \]

and perform basic math with them

\[ \begin{align} z_1 &= x_1 + \epsilon y_1 \\ z_2 &= x_2 + \epsilon y_2 \end{align} \]

Addition: \[ z_1 + z_2 = (x_1 + x_2) + \epsilon (y_1 + y_2) \]

Multiplication: \[ \begin{align} z_1 z_2 &= x_1 x_2 + \epsilon (x_1 y_2 + x_2 y_1) + \epsilon^2(x_1^2 + y_1^2)\\ &= x_1 x_2 + \epsilon (x_1 y_2 + x_2 y_1) + \cancelto{0}{\epsilon^2}(x_1^2 + y_1^2)\\ &= x_1 x_2 + \epsilon (x_1 y_2 + x_2 y_1) \end{align} \]

Equation of a line

Take a simple function - equation of a line

\[ f(x) = m x + b \]

We know that the derivative of \(f(x)\) is the slope of the line \(m\)

\[ \frac{d f(x)}{dx} = m \]

Suppose that we substitute a dual number for \(x\), let \(x\rightarrow x+\epsilon\). Then

\[ \begin{align} f(x+\epsilon) &= m(x+\epsilon) + b \\ &= mx + b + m\epsilon \\ &= f(x) + \frac{df(x)}{dx} \epsilon \end{align} \]

and we have

• real part is \(mx +b\) the value of the line
• dual part is \(m\) which is the derivative/slope of the line.

This provides a hint of how dual numbers can be used to find both the value and derivative of a function [2] [3].

Fancier functions

\[ f(a+b\epsilon) \]

can be expanded using a Taylor series \[ f(a+b \epsilon) = f(a) + f'(a) b\epsilon + \frac{1}{2} f''(a)(b\epsilon)^2 + \cdots \]

notice that all the terms after \(f'\) are multiplied by \(\epsilon^2\) and are therefore zero, so the series terminates \[ f(a+b\epsilon) = f(a) + f'(a) b \epsilon \]

This is a way to obtain both the value of the function and it’s derivative at the point \(a\).

We can also use this to express other functions as dual numbers [4]

\[ \sin(x+\epsilon) = \sin(x) + \cos(x)\epsilon \]

Chain rule

The chain rule

\[ \frac{d}{dx} g(f(x)) = \frac{d g}{d f} \frac{df}{dx} \]

also applies

\[ \begin{align} g(f(x + \epsilon)) &= g\left(f(x) + \frac{df(x)}{dx}\epsilon\right) \\ &= g\left(f(x)\right) + \frac{dg}{df} \frac{df}{dx} \epsilon \end{align} \]

so we don’t need to know how to symbolically compute the form of any derivative, it’s automatically computed for us [5].

What we need is software that handles the dual-number arithmetic for us.

Revisit gradient descent with auto-diff

Uses the autograd python package [6].

1import autograd.numpy as np
2from autograd import grad

def gradient_descent_autodiff(
        x0,              # starting point
        func,            # f(x)
        eta,             # learning rate
        tolerance = 0.01, # termination condition
        max_steps = 100, # maximum steps
        ):
    
    x = x0 # initialize
    xs = [x0] # list of steps

    dfdx = grad(func)
    
    for _ in range(max_steps):        
        delta = eta * dfdx(x)
        if np.abs(delta) < tolerance:
            break
        x = x - delta
        xs.append(x)

    return x, xs

1

Applies a thin wrapper

2

Creates gradient function

Define function to find the minimum of (same as before).

def func(x):
    return x**2 - 4*x + 1

%%capture
from IPython.display import HTML
import matplotlib.animation as animation

xres,xhist = gradient_descent_autodiff(8.0,func,0.1,)

fig, ax = plt.subplots()
# plot underlying function
xvals = np.linspace(-6,10,100)
ax.plot(xvals,func(xvals),'b-')
ax.grid()

def animate(ii,):
    xh = np.array(xhist[0:ii])
    ax.plot(xh,func(xh),'ro-')
    for ii,_x in enumerate(xh):
        ax.annotate(f'{ii}', (_x+0.4, func(_x)-1.5,),color='r')

ani = animation.FuncAnimation(fig,animate,repeat=False,frames=len(xhist)-1,interval=1000)

HTML(ani.to_jshtml())

Once Loop Reflect

Figure 4: Gradient descent for \(\eta=0.1\) using auto-differentiation

Not surprisingly, the results as the same as before.

Auto-diff software packages - python focused

Notice that many of these software packages have similar names.

• autograd python package [6].
  • \(\rightarrow\) the example already shown
• auto-diff python package [7].
  • \(\rightarrow\) problems with numpy>2 (deprecated features)
• autodiff cython package [8].
  • problems building
• \(xa\partial\) c++ with python bindings [9].

More complicated example using gradient descent to solve

Imagine that you have a large vector in 2-D with some vector on the end that rotates around.

fig,ax = plt.subplots()

bigvec = [50_000,0]

ax.quiver(0,0,*bigvec,angles='xy',scale_units='xy',scale=1,color='r')

angle = np.pi/4

vec = np.array([5000,0])
dcm = np.array( [ [ np.cos(angle), -np.sin(angle)  ],
                  [  np.sin(angle),  np.cos(angle) ] ] )
vecr = dcm @ vec
ax.quiver(bigvec[0],bigvec[1],vecr[0],vecr[1],angles='xy',scale_units='xy',scale=1,color='b')

ax.quiver(0,0,bigvec[0]+vecr[0], bigvec[1]+vecr[1],angles='xy',scale_units='xy',scale=1,color='g')

ax.set_xlim(0,60000)
ax.set_ylim(-6000,6000)

Figure 5: Vector sum (green) of large vector (red) and small vector (blue).

fig, ax = plt.subplots()

bigvec = [50_000*np.cos(np.pi/4), 50_000*np.sin(np.pi/4)]
angles = np.linspace(0,2*np.pi,100)
vec = np.array([5000,0])



def animate(ii):
    angle=angles[ii]

    dcm = np.array( [ [ np.cos(angle), -np.sin(angle)  ],
                      [  np.sin(angle),  np.cos(angle) ] ] )
    vecr = dcm @ vec
    ax.clear()
    ax.set_xlim(0,60000)
    ax.set_ylim(0,60000)
    ax.grid()
    ax.quiver(0,0,*bigvec,angles='xy',scale_units='xy',scale=1,color='r')
    ax.quiver(bigvec[0],bigvec[1],vecr[0],vecr[1],angles='xy',scale_units='xy',scale=1,color='b')
    ax.quiver(0,0,bigvec[0]+vecr[0], bigvec[1]+vecr[1],angles='xy',scale_units='xy',scale=1,color='g')
    
ani = animation.FuncAnimation(fig,animate,repeat=True,frames=len(angles)-1,interval=100)

HTML(ani.to_jshtml())

Once Loop Reflect

Figure 6: Small rotating vector added to large vector

If we don’t know how large the small vector is, can we find it if we can only measure the length of the green vector? (assuming we know the red vector).

Function that computes magnitude of sum

Let’s define a function that computes the green vector magnitude magsum

def magsum(x1,x2):
    s = x1+x2
    return np.linalg.norm(s,axis=-1)

2D Direction cosine matrix for rotations

def dcm(angle):
    """
    Parameters
    ==========
    angle : N,
       array of angles of dimension (N,) for N angles

    Returns
    =======
    dcm : 2,2,N or 2,2
       array of N direction cosing matrices in 2D of shape (2,2,N) if N>1
       array of (2,2) if N=1
    """
    return np.array( [ [ np.cos(angle), -np.sin(angle) ],
                       [ np.sin(angle),  np.cos(angle) ] ] )

bigvec = np.array([[50_000,0]]) # [1x2]
smallvec = np.array([[5000,0]]) # [1x2]
num_angles = 100
angles = np.linspace(0,2*np.pi,num_angles) 

vecr = np.einsum('ijk,mj->ki',dcm(angles),smallvec)
ms = magsum(bigvec,vecr)

ms is the magnitude that we measure.

fig,ax = plt.subplots()
ax.plot(angles,ms)
ax.set_xlabel('Angle [radians]')
ax.set_ylabel('Magnitude')
ax.grid()

Figure 7: Magnitude of vector sum as function of angle

Can we find the components of the small vector from just this rotation data?

Find the minimum of

\[ \min \left\{ |x_\text{tot}| - | x_1 + x_2 | \right\} \]

Define the python function to minimize

from functools import partial

def funcfull(
        magtot,
        x1,
        x2,
        angle,
        ):

    return magtot - magsum(x1,np.einsum('ijk,j->ki',dcm(angles),x2))

func = partial(funcfull,magtot=ms, x1=bigvec, angle=angles)

Test out func

print(f'func(x2=np.array([5000,0])) = {func(x2=np.array([5000,0]))}')

func(x2=np.array([5000,0])) = [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.
 0. 0. 0. 0.]

All zeros as expected for the correct input.

This is not quite the function we want to minimize. We need a function that returns a single value. So, let’s try:

def funcmin(x):
    vals = func(x2=x)
    
    return np.sqrt(np.einsum('i,i->',vals,vals))

funcmin returns a scalar value:

print(f'{funcmin(np.array([6000,0]))}')

7106.247064525054

print(f'{funcmin(np.array([5000,0]))}')

0.0

And now we have a function to minimize.

Run gradient descent on our function

Update the gradient descent function slightly to account for the fact that we are not minimizing with vector instead of scalar input.

import autograd.numpy as np 
from autograd import grad   

def gradient_descent_autodiff(
        x0,              # starting point
        func,            # f(x)
        eta,             # learning rate
        tolerance = 0.01, # termination condition
        max_steps = 100, # maximum steps
        ):
    
    x = x0 # initialize
    xs = [x0] # list of steps

    dfdx = grad(func)
    
    for _ in range(max_steps):        
        delta = eta * dfdx(x)
1        deltamag = np.sqrt(np.einsum('i,i->',delta,delta))
        if np.abs(deltamag) < tolerance:
            break
        x = x - delta
        xs.append(x)

    return x, xs

1

Compute magnitude of vector change to compare with tolerance.

x0 = np.array([6000.0,6000.0])
xmin,xs = gradient_descent_autodiff(x0,funcmin,1.0,max_steps=1000)
print(xmin)
print(len(xs))
print(xs)

[5.00406769e+03 1.10312226e-01]
1001
[array([6000., 6000.]), array([5998.81489656, 5993.0699272 ]), array([5997.62983138, 5986.13982864]), array([5996.44480451, 5979.20970433]), array([5995.25981597, 5972.27955429]), array([5994.07486582, 5965.34937854]), array([5992.88995408, 5958.41917711]), array([5991.7050808 , 5951.48895002]), array([5990.52024601, 5944.55869728]), array([5989.33544975, 5937.62841891]), array([5988.15069206, 5930.69811495]), array([5986.96597297, 5923.7677854 ]), array([5985.78129253, 5916.83743029]), array([5984.59665077, 5909.90704964]), array([5983.41204772, 5902.97664347]), array([5982.22748344, 5896.0462118 ]), array([5981.04295796, 5889.11575465]), array([5979.85847131, 5882.18527204]), array([5978.67402353, 5875.25476399]), array([5977.48961466, 5868.32423053]), array([5976.30524475, 5861.39367167]), array([5975.12091383, 5854.46308743]), array([5973.93662193, 5847.53247784]), array([5972.7523691 , 5840.60184291]), array([5971.56815538, 5833.67118266]), array([5970.3839808 , 5826.74049712]), array([5969.19984541, 5819.80978631]), array([5968.01574924, 5812.87905024]), array([5966.83169233, 5805.94828894]), array([5965.64767473, 5799.01750243]), array([5964.46369647, 5792.08669072]), array([5963.27975758, 5785.15585384]), array([5962.09585812, 5778.22499181]), array([5960.91199812, 5771.29410464]), array([5959.72817762, 5764.36319237]), array([5958.54439665, 5757.432255  ]), array([5957.36065527, 5750.50129256]), array([5956.1769535 , 5743.57030507]), array([5954.9932914 , 5736.63929255]), array([5953.80966899, 5729.70825502]), array([5952.62608633, 5722.77719251]), array([5951.44254344, 5715.84610502]), array([5950.25904037, 5708.91499258]), array([5949.07557716, 5701.98385521]), array([5947.89215385, 5695.05269293]), array([5946.70877049, 5688.12150576]), array([5945.5254271 , 5681.19029373]), array([5944.34212374, 5674.25905684]), array([5943.15886044, 5667.32779513]), array([5941.97563725, 5660.3965086 ]), array([5940.7924542 , 5653.46519729]), array([5939.60931134, 5646.53386121]), array([5938.4262087 , 5639.60250038]), array([5937.24314633, 5632.67111482]), array([5936.06012427, 5625.73970455]), array([5934.87714257, 5618.8082696 ]), array([5933.69420125, 5611.87680997]), array([5932.51130037, 5604.9453257 ]), array([5931.32843997, 5598.01381679]), array([5930.14562008, 5591.08228328]), array([5928.96284075, 5584.15072518]), array([5927.78010202, 5577.21914251]), array([5926.59740393, 5570.28753528]), array([5925.41474653, 5563.35590353]), array([5924.23212986, 5556.42424727]), array([5923.04955395, 5549.49256652]), array([5921.86701886, 5542.5608613 ]), array([5920.68452462, 5535.62913163]), array([5919.50207127, 5528.69737753]), array([5918.31965887, 5521.76559902]), array([5917.13728744, 5514.83379611]), array([5915.95495704, 5507.90196884]), array([5914.77266771, 5500.97011721]), array([5913.59041948, 5494.03824125]), array([5912.40821241, 5487.10634097]), array([5911.22604653, 5480.17441641]), array([5910.04392189, 5473.24246757]), array([5908.86183854, 5466.31049448]), array([5907.67979651, 5459.37849715]), array([5906.49779584, 5452.44647561]), array([5905.31583659, 5445.51442987]), array([5904.1339188 , 5438.58235996]), array([5902.95204251, 5431.65026589]), array([5901.77020775, 5424.71814769]), array([5900.58841459, 5417.78600537]), array([5899.40666306, 5410.85383895]), array([5898.2249532 , 5403.92164846]), array([5897.04328506, 5396.9894339 ]), array([5895.86165868, 5390.05719531]), array([5894.68007412, 5383.12493269]), array([5893.4985314 , 5376.19264608]), array([5892.31703058, 5369.26033548]), array([5891.13557171, 5362.32800093]), array([5889.95415481, 5355.39564243]), array([5888.77277995, 5348.46326   ]), array([5887.59144717, 5341.53085368]), array([5886.41015651, 5334.59842347]), array([5885.22890801, 5327.66596939]), array([5884.04770172, 5320.73349147]), array([5882.86653769, 5313.80098972]), array([5881.68541596, 5306.86846416]), array([5880.50433658, 5299.93591482]), array([5879.32329959, 5293.0033417 ]), array([5878.14230505, 5286.07074484]), array([5876.96135298, 5279.13812424]), array([5875.78044344, 5272.20547993]), array([5874.59957649, 5265.27281193]), array([5873.41875215, 5258.34012026]), array([5872.23797048, 5251.40740493]), array([5871.05723153, 5244.47466596]), array([5869.87653534, 5237.54190338]), array([5868.69588195, 5230.6091172 ]), array([5867.51527142, 5223.67630744]), array([5866.33470379, 5216.74347412]), array([5865.15417911, 5209.81061727]), array([5863.97369742, 5202.87773689]), array([5862.79325878, 5195.944833  ]), array([5861.61286322, 5189.01190563]), array([5860.4325108, 5182.0789548]), array([5859.25220156, 5175.14598052]), array([5858.07193555, 5168.21298281]), array([5856.89171283, 5161.27996169]), array([5855.71153342, 5154.34691719]), array([5854.53139739, 5147.41384931]), array([5853.35130479, 5140.48075808]), array([5852.17125565, 5133.54764351]), array([5850.99125003, 5126.61450563]), array([5849.81128797, 5119.68134446]), array([5848.63136953, 5112.74816   ]), array([5847.45149475, 5105.81495229]), array([5846.27166368, 5098.88172133]), array([5845.09187637, 5091.94846716]), array([5843.91213287, 5085.01518978]), array([5842.73243323, 5078.08188921]), array([5841.55277749, 5071.14856548]), array([5840.37316571, 5064.2152186 ]), array([5839.19359793, 5057.28184859]), array([5838.01407421, 5050.34845547]), array([5836.83459459, 5043.41503926]), array([5835.65515912, 5036.48159998]), array([5834.47576785, 5029.54813763]), array([5833.29642084, 5022.61465226]), array([5832.11711813, 5015.68114386]), array([5830.93785976, 5008.74761246]), array([5829.7586458 , 5001.81405808]), array([5828.5794763 , 4994.88048073]), array([5827.40035129, 4987.94688044]), array([5826.22127084, 4981.01325723]), array([5825.04223499, 4974.0796111 ]), array([5823.8632438 , 4967.14594208]), array([5822.68429731, 4960.21225019]), array([5821.50539558, 4953.27853544]), array([5820.32653865, 4946.34479786]), array([5819.14772658, 4939.41103746]), array([5817.96895942, 4932.47725426]), array([5816.79023722, 4925.54344828]), array([5815.61156004, 4918.60961953]), array([5814.43292792, 4911.67576804]), array([5813.25434091, 4904.74189382]), array([5812.07579907, 4897.80799689]), array([5810.89730246, 4890.87407727]), array([5809.71885111, 4883.94013497]), array([5808.54044509, 4877.00617002]), array([5807.36208444, 4870.07218243]), array([5806.18376923, 4863.13817222]), array([5805.00549949, 4856.2041394 ]), array([5803.82727529, 4849.270084  ]), array([5802.64909668, 4842.33600604]), array([5801.47096371, 4835.40190552]), array([5800.29287643, 4828.46778248]), array([5799.11483491, 4821.53363692]), array([5797.93683918, 4814.59946886]), array([5796.7588893 , 4807.66527833]), array([5795.58098534, 4800.73106533]), array([5794.40312733, 4793.79682989]), array([5793.22531534, 4786.86257203]), array([5792.04754942, 4779.92829176]), array([5790.86982963, 4772.9939891 ]), array([5789.69215601, 4766.05966407]), array([5788.51452863, 4759.12531668]), array([5787.33694753, 4752.19094695]), array([5786.15941278, 4745.25655491]), array([5784.98192443, 4738.32214056]), array([5783.80448252, 4731.38770393]), array([5782.62708713, 4724.45324503]), array([5781.4497383 , 4717.51876387]), array([5780.27243609, 4710.58426049]), array([5779.09518055, 4703.64973489]), array([5777.91797174, 4696.71518709]), array([5776.74080972, 4689.78061711]), array([5775.56369454, 4682.84602497]), array([5774.38662625, 4675.91141068]), array([5773.20960493, 4668.97677426]), array([5772.03263061, 4662.04211572]), array([5770.85570336, 4655.1074351 ]), array([5769.67882324, 4648.17273239]), array([5768.50199029, 4641.23800763]), array([5767.32520459, 4634.30326082]), array([5766.14846618, 4627.36849198]), array([5764.97177512, 4620.43370113]), array([5763.79513148, 4613.49888829]), array([5762.6185353 , 4606.56405348]), array([5761.44198665, 4599.6291967 ]), array([5760.26548558, 4592.69431799]), array([5759.08903216, 4585.75941735]), array([5757.91262644, 4578.8244948 ]), array([5756.73626847, 4571.88955036]), array([5755.55995832, 4564.95458404]), array([5754.38369605, 4558.01959587]), array([5753.20748172, 4551.08458585]), array([5752.03131537, 4544.14955402]), array([5750.85519708, 4537.21450037]), array([5749.67912691, 4530.27942493]), array([5748.5031049 , 4523.34432772]), array([5747.32713113, 4516.40920875]), array([5746.15120565, 4509.47406804]), array([5744.97532852, 4502.53890561]), array([5743.79949981, 4495.60372146]), array([5742.62371956, 4488.66851563]), array([5741.44798785, 4481.73328812]), array([5740.27230473, 4474.79803895]), array([5739.09667027, 4467.86276814]), array([5737.92108452, 4460.92747571]), array([5736.74554755, 4453.99216166]), array([5735.57005941, 4447.05682602]), array([5734.39462018, 4440.12146881]), array([5733.2192299 , 4433.18609003]), array([5732.04388865, 4426.25068971]), array([5730.86859648, 4419.31526787]), array([5729.69335346, 4412.37982451]), array([5728.51815965, 4405.44435965]), array([5727.34301511, 4398.50887331]), array([5726.1679199 , 4391.57336552]), array([5724.99287409, 4384.63783627]), array([5723.81787774, 4377.70228559]), array([5722.64293091, 4370.7667135 ]), array([5721.46803368, 4363.83112001]), array([5720.29318609, 4356.89550513]), array([5719.11838821, 4349.95986889]), array([5717.94364011, 4343.02421129]), array([5716.76894186, 4336.08853236]), array([5715.59429351, 4329.15283211]), array([5714.41969513, 4322.21711055]), array([5713.24514678, 4315.28136771]), array([5712.07064854, 4308.34560359]), array([5710.89620046, 4301.40981822]), array([5709.72180261, 4294.4740116 ]), array([5708.54745506, 4287.53818376]), array([5707.37315787, 4280.60233471]), array([5706.1989111 , 4273.66646446]), array([5705.02471483, 4266.73057304]), array([5703.85056912, 4259.79466045]), array([5702.67647403, 4252.85872671]), array([5701.50242963, 4245.92277184]), array([5700.32843599, 4238.98679585]), array([5699.15449318, 4232.05079877]), array([5697.98060126, 4225.11478059]), array([5696.8067603 , 4218.17874134]), array([5695.63297037, 4211.24268104]), array([5694.45923153, 4204.3065997 ]), array([5693.28554385, 4197.37049733]), array([5692.11190741, 4190.43437395]), array([5690.93832226, 4183.49822957]), array([5689.76478848, 4176.56206422]), array([5688.59130614, 4169.6258779 ]), array([5687.4178753 , 4162.68967062]), array([5686.24449604, 4155.75344242]), array([5685.07116842, 4148.81719329]), array([5683.89789251, 4141.88092326]), array([5682.72466839, 4134.94463234]), array([5681.55149611, 4128.00832054]), array([5680.37837576, 4121.07198788]), array([5679.20530741, 4114.13563438]), array([5678.03229111, 4107.19926004]), array([5676.85932695, 4100.26286489]), array([5675.686415  , 4093.32644894]), array([5674.51355532, 4086.3900122 ]), array([5673.34074799, 4079.45355469]), array([5672.16799308, 4072.51707641]), array([5670.99529066, 4065.5805774 ]), array([5669.8226408 , 4058.64405766]), array([5668.65004358, 4051.7075172 ]), array([5667.47749907, 4044.77095604]), array([5666.30500733, 4037.83437419]), array([5665.13256845, 4030.89777168]), array([5663.9601825 , 4023.96114851]), array([5662.78784955, 4017.02450469]), array([5661.61556967, 4010.08784025]), array([5660.44334294, 4003.15115519]), array([5659.27116943, 3996.21444953]), array([5658.09904921, 3989.27772328]), array([5656.92698237, 3982.34097647]), array([5655.75496897, 3975.40420909]), array([5654.58300909, 3968.46742117]), array([5653.41110281, 3961.53061272]), array([5652.23925019, 3954.59378376]), array([5651.06745133, 3947.65693429]), array([5649.89570628, 3940.72006433]), array([5648.72401514, 3933.7831739 ]), array([5647.55237797, 3926.84626301]), array([5646.38079485, 3919.90933167]), array([5645.20926586, 3912.97237989]), array([5644.03779108, 3906.0354077 ]), array([5642.86637058, 3899.0984151 ]), array([5641.69500444, 3892.1614021 ]), array([5640.52369274, 3885.22436873]), array([5639.35243555, 3878.28731499]), array([5638.18123296, 3871.3502409 ]), array([5637.01008505, 3864.41314647]), array([5635.83899189, 3857.47603172]), array([5634.66795356, 3850.53889665]), array([5633.49697014, 3843.60174128]), array([5632.32604171, 3836.66456563]), array([5631.15516836, 3829.7273697 ]), array([5629.98435015, 3822.79015351]), array([5628.81358718, 3815.85291708]), array([5627.64287952, 3808.91566041]), array([5626.47222725, 3801.97838352]), array([5625.30163046, 3795.04108643]), array([5624.13108923, 3788.10376914]), array([5622.96060363, 3781.16643166]), array([5621.79017375, 3774.22907402]), array([5620.61979968, 3767.29169622]), array([5619.44948149, 3760.35429828]), array([5618.27921927, 3753.41688021]), array([5617.1090131 , 3746.47944202]), array([5615.93886307, 3739.54198372]), array([5614.76876925, 3732.60450533]), array([5613.59873174, 3725.66700685]), array([5612.42875061, 3718.72948831]), array([5611.25882596, 3711.79194971]), array([5610.08895786, 3704.85439107]), array([5608.9191464, 3697.9168124]), array([5607.74939167, 3690.97921371]), array([5606.57969376, 3684.04159501]), array([5605.41005274, 3677.10395631]), array([5604.24046871, 3670.16629764]), array([5603.07094175, 3663.22861899]), array([5601.90147194, 3656.29092038]), array([5600.73205939, 3649.35320183]), array([5599.56270417, 3642.41546334]), array([5598.39340637, 3635.47770493]), array([5597.22416608, 3628.5399266 ]), array([5596.05498338, 3621.60212838]), array([5594.88585838, 3614.66431027]), array([5593.71679115, 3607.72647228]), array([5592.54778178, 3600.78861443]), array([5591.37883037, 3593.85073672]), array([5590.20993701, 3586.91283917]), array([5589.04110178, 3579.97492179]), array([5587.87232478, 3573.0369846 ]), array([5586.7036061 , 3566.09902759]), array([5585.53494582, 3559.16105079]), array([5584.36634404, 3552.2230542 ]), array([5583.19780086, 3545.28503784]), array([5582.02931636, 3538.34700172]), array([5580.86089064, 3531.40894584]), array([5579.69252379, 3524.47087023]), array([5578.5242159 , 3517.53277488]), array([5577.35596707, 3510.59465982]), array([5576.18777739, 3503.65652504]), array([5575.01964695, 3496.71837057]), array([5573.85157586, 3489.78019641]), array([5572.6835642 , 3482.84200258]), array([5571.51561207, 3475.90378908]), array([5570.34771957, 3468.96555593]), array([5569.17988679, 3462.02730313]), array([5568.01211382, 3455.08903069]), array([5566.84440078, 3448.15073864]), array([5565.67674774, 3441.21242697]), array([5564.50915482, 3434.27409569]), array([5563.3416221 , 3427.33574483]), array([5562.17414969, 3420.39737438]), array([5561.00673769, 3413.45898436]), array([5559.83938619, 3406.52057477]), array([5558.67209529, 3399.58214563]), array([5557.50486509, 3392.64369695]), array([5556.3376957 , 3385.70522874]), array([5555.17058721, 3378.766741  ]), array([5554.00353973, 3371.82823375]), array([5552.83655335, 3364.889707  ]), array([5551.66962817, 3357.95116075]), array([5550.50276431, 3351.01259501]), array([5549.33596186, 3344.0740098 ]), array([5548.16922092, 3337.13540513]), array([5547.0025416, 3330.196781 ]), array([5545.835924  , 3323.25813742]), array([5544.66936823, 3316.3194744 ]), array([5543.50287438, 3309.38079195]), array([5542.33644257, 3302.44209008]), array([5541.1700729, 3295.5033688]), array([5540.00376548, 3288.56462812]), array([5538.8375204 , 3281.62586805]), array([5537.67133778, 3274.68708859]), array([5536.50521773, 3267.74828975]), array([5535.33916035, 3260.80947155]), array([5534.17316574, 3253.87063399]), array([5533.00723403, 3246.93177708]), array([5531.84136531, 3239.99290082]), array([5530.6755597 , 3233.05400524]), array([5529.5098173 , 3226.11509032]), array([5528.34413822, 3219.17615609]), array([5527.17852258, 3212.23720256]), array([5526.01297049, 3205.29822972]), array([5524.84748205, 3198.35923759]), array([5523.68205738, 3191.42022617]), array([5522.51669658, 3184.48119548]), array([5521.35139979, 3177.54214552]), array([5520.18616709, 3170.6030763 ]), array([5519.02099862, 3163.66398782]), array([5517.85589448, 3156.7248801 ]), array([5516.69085479, 3149.78575313]), array([5515.52587966, 3142.84660694]), array([5514.36096921, 3135.90744153]), array([5513.19612355, 3128.96825689]), array([5512.0313428 , 3122.02905305]), array([5510.86662708, 3115.08983   ]), array([5509.7019765 , 3108.15058776]), array([5508.53739119, 3101.21132634]), array([5507.37287125, 3094.27204572]), array([5506.20841682, 3087.33274594]), array([5505.044028  , 3080.39342699]), array([5503.87970491, 3073.45408887]), array([5502.71544769, 3066.5147316 ]), array([5501.55125645, 3059.57535518]), array([5500.3871313 , 3052.63595962]), array([5499.22307238, 3045.69654493]), array([5498.0590798, 3038.7571111]), array([5496.89515369, 3031.81765815]), array([5495.73129417, 3024.87818609]), array([5494.56750136, 3017.93869491]), array([5493.40377539, 3010.99918462]), array([5492.24011638, 3004.05965524]), array([5491.07652447, 2997.12010676]), array([5489.91299976, 2990.1805392 ]), array([5488.7495424 , 2983.24095255]), array([5487.58615251, 2976.30134682]), array([5486.42283022, 2969.36172202]), array([5485.25957565, 2962.42207815]), array([5484.09638893, 2955.48241523]), array([5482.9332702 , 2948.54273324]), array([5481.77021958, 2941.6030322 ]), array([5480.6072372 , 2934.66331212]), array([5479.4443232 , 2927.72357299]), array([5478.2814777 , 2920.78381483]), array([5477.11870085, 2913.84403763]), array([5475.95599276, 2906.90424141]), array([5474.79335358, 2899.96442616]), array([5473.63078344, 2893.02459189]), array([5472.46828247, 2886.0847386 ]), array([5471.30585081, 2879.1448663 ]), array([5470.1434886, 2872.204975 ]), array([5468.98119597, 2865.26506469]), array([5467.81897305, 2858.32513539]), array([5466.65682   , 2851.38518708]), array([5465.49473694, 2844.44521979]), array([5464.33272401, 2837.5052335 ]), array([5463.17078135, 2830.56522823]), array([5462.00890911, 2823.62520398]), array([5460.84710743, 2816.68516075]), array([5459.68537644, 2809.74509854]), array([5458.52371629, 2802.80501736]), array([5457.36212713, 2795.86491721]), array([5456.20060908, 2788.92479809]), array([5455.03916231, 2781.98466   ]), array([5453.87778696, 2775.04450296]), array([5452.71648316, 2768.10432695]), array([5451.55525107, 2761.16413199]), array([5450.39409084, 2754.22391807]), array([5449.23300261, 2747.2836852 ]), array([5448.07198653, 2740.34343338]), array([5446.91104274, 2733.40316261]), array([5445.75017141, 2726.46287289]), array([5444.58937268, 2719.52256422]), array([5443.4286467 , 2712.58223661]), array([5442.26799362, 2705.64189006]), array([5441.1074136 , 2698.70152457]), array([5439.94690679, 2691.76114013]), array([5438.78647335, 2684.82073676]), array([5437.62611342, 2677.88031444]), array([5436.46582717, 2670.93987319]), array([5435.30561475, 2663.99941301]), array([5434.14547632, 2657.05893388]), array([5432.98541204, 2650.11843582]), array([5431.82542207, 2643.17791882]), array([5430.66550656, 2636.23738289]), array([5429.50566568, 2629.29682802]), array([5428.34589959, 2622.35625422]), array([5427.18620845, 2615.41566148]), array([5426.02659242, 2608.47504981]), array([5424.86705167, 2601.53441919]), array([5423.70758636, 2594.59376965]), array([5422.54819666, 2587.65310116]), array([5421.38888273, 2580.71241373]), array([5420.22964475, 2573.77170737]), array([5419.07048287, 2566.83098206]), array([5417.91139727, 2559.89023781]), array([5416.75238812, 2552.94947462]), array([5415.59345559, 2546.00869248]), array([5414.43459985, 2539.06789139]), array([5413.27582107, 2532.12707136]), array([5412.11711943, 2525.18623237]), array([5410.9584951 , 2518.24537443]), array([5409.79994825, 2511.30449754]), array([5408.64147906, 2504.36360168]), array([5407.48308771, 2497.42268687]), array([5406.32477438, 2490.4817531 ]), array([5405.16653925, 2483.54080036]), array([5404.00838249, 2476.59982865]), array([5402.85030428, 2469.65883796]), array([5401.69230482, 2462.7178283 ]), array([5400.53438427, 2455.77679967]), array([5399.37654283, 2448.83575205]), array([5398.21878068, 2441.89468544]), array([5397.061098  , 2434.95359984]), array([5395.90349499, 2428.01249525]), array([5394.74597182, 2421.07137166]), array([5393.5885287 , 2414.13022907]), array([5392.4311658 , 2407.18906747]), array([5391.27388332, 2400.24788685]), array([5390.11668145, 2393.30668722]), array([5388.95956039, 2386.36546856]), array([5387.80252032, 2379.42423088]), array([5386.64556145, 2372.48297416]), array([5385.48868397, 2365.5416984 ]), array([5384.33188807, 2358.6004036 ]), array([5383.17517397, 2351.65908975]), array([5382.01854184, 2344.71775684]), array([5380.86199191, 2337.77640486]), array([5379.70552436, 2330.83503382]), array([5378.5491394, 2323.8936437]), array([5377.39283724, 2316.95223449]), array([5376.23661808, 2310.0108062 ]), array([5375.08048213, 2303.0693588 ]), array([5373.92442959, 2296.1278923 ]), array([5372.76846068, 2289.18640669]), array([5371.6125756 , 2282.24490195]), array([5370.45677457, 2275.30337809]), array([5369.30105779, 2268.36183508]), array([5368.14542549, 2261.42027293]), array([5366.98987788, 2254.47869163]), array([5365.83441517, 2247.53709116]), array([5364.67903758, 2240.59547152]), array([5363.52374534, 2233.6538327 ]), array([5362.36853865, 2226.71217469]), array([5361.21341775, 2219.77049748]), array([5360.05838286, 2212.82880105]), array([5358.90343419, 2205.8870854 ]), array([5357.74857199, 2198.94535053]), array([5356.59379646, 2192.00359641]), array([5355.43910785, 2185.06182303]), array([5354.28450638, 2178.1200304 ]), array([5353.12999228, 2171.17821849]), array([5351.97556579, 2164.23638729]), array([5350.82122713, 2157.2945368 ]), array([5349.66697656, 2150.35266699]), array([5348.51281429, 2143.41077786]), array([5347.35874058, 2136.4688694 ]), array([5346.20475566, 2129.52694159]), array([5345.05085977, 2122.58499443]), array([5343.89705316, 2115.64302789]), array([5342.74333606, 2108.70104196]), array([5341.58970874, 2101.75903663]), array([5340.43617143, 2094.81701189]), array([5339.28272438, 2087.87496772]), array([5338.12936785, 2080.93290411]), array([5336.97610208, 2073.99082105]), array([5335.82292734, 2067.04871851]), array([5334.66984387, 2060.10659649]), array([5333.51685194, 2053.16445496]), array([5332.3639518 , 2046.22229392]), array([5331.21114371, 2039.28011334]), array([5330.05842795, 2032.33791322]), array([5328.90580476, 2025.39569352]), array([5327.75327443, 2018.45345425]), array([5326.6008372 , 2011.51119538]), array([5325.44849337, 2004.56891689]), array([5324.29624319, 1997.62661876]), array([5323.14408694, 1990.68430098]), array([5321.9920249 , 1983.74196354]), array([5320.84005734, 1976.7996064 ]), array([5319.68818453, 1969.85722956]), array([5318.53640678, 1962.91483299]), array([5317.38472434, 1955.97241667]), array([5316.23313752, 1949.02998059]), array([5315.08164659, 1942.08752473]), array([5313.93025184, 1935.14504905]), array([5312.77895357, 1928.20255356]), array([5311.62775207, 1921.26003821]), array([5310.47664763, 1914.317503  ]), array([5309.32564055, 1907.37494789]), array([5308.17473112, 1900.43237288]), array([5307.02391966, 1893.48977793]), array([5305.87320646, 1886.54716302]), array([5304.72259182, 1879.60452814]), array([5303.57207606, 1872.66187325]), array([5302.42165948, 1865.71919833]), array([5301.2713424 , 1858.77650337]), array([5300.12112513, 1851.83378833]), array([5298.97100799, 1844.89105319]), array([5297.82099129, 1837.94829793]), array([5296.67107536, 1831.00552252]), array([5295.52126052, 1824.06272693]), array([5294.3715471 , 1817.11991115]), array([5293.22193542, 1810.17707513]), array([5292.07242581, 1803.23421886]), array([5290.92301862, 1796.29134231]), array([5289.77371417, 1789.34844544]), array([5288.6245128 , 1782.40552824]), array([5287.47541485, 1775.46259068]), array([5286.32642067, 1768.51963271]), array([5285.1775306 , 1761.57665433]), array([5284.02874499, 1754.63365549]), array([5282.8800642 , 1747.69063616]), array([5281.73148857, 1740.74759632]), array([5280.58301847, 1733.80453594]), array([5279.43465425, 1726.86145498]), array([5278.28639627, 1719.9183534 ]), array([5277.13824491, 1712.97523119]), array([5275.99020052, 1706.0320883 ]), array([5274.84226349, 1699.08892471]), array([5273.69443417, 1692.14574037]), array([5272.54671296, 1685.20253526]), array([5271.39910023, 1678.25930934]), array([5270.25159637, 1671.31606257]), array([5269.10420176, 1664.37279492]), array([5267.95691678, 1657.42950636]), array([5266.80974184, 1650.48619684]), array([5265.66267733, 1643.54286633]), array([5264.51572364, 1636.59951479]), array([5263.36888118, 1629.65614219]), array([5262.22215036, 1622.71274848]), array([5261.07553159, 1615.76933363]), array([5259.92902527, 1608.82589759]), array([5258.78263182, 1601.88244033]), array([5257.63635167, 1594.9389618 ]), array([5256.49018523, 1587.99546197]), array([5255.34413293, 1581.05194079]), array([5254.1981952 , 1574.10839822]), array([5253.05237249, 1567.16483421]), array([5251.90666521, 1560.22124873]), array([5250.76107382, 1553.27764172]), array([5249.61559877, 1546.33401315]), array([5248.47024049, 1539.39036296]), array([5247.32499945, 1532.44669112]), array([5246.1798761 , 1525.50299757]), array([5245.0348709 , 1518.55928227]), array([5243.88998432, 1511.61554517]), array([5242.74521683, 1504.67178621]), array([5241.6005689 , 1497.72800536]), array([5240.45604101, 1490.78420256]), array([5239.31163364, 1483.84037776]), array([5238.16734729, 1476.89653091]), array([5237.02318243, 1469.95266196]), array([5235.87913958, 1463.00877085]), array([5234.73521922, 1456.06485753]), array([5233.59142186, 1449.12092195]), array([5232.44774802, 1442.17696405]), array([5231.30419821, 1435.23298378]), array([5230.16077294, 1428.28898107]), array([5229.01747275, 1421.34495588]), array([5227.87429816, 1414.40090813]), array([5226.73124971, 1407.45683779]), array([5225.58832793, 1400.51274478]), array([5224.44553338, 1393.56862904]), array([5223.30286661, 1386.62449052]), array([5222.16032816, 1379.68032915]), array([5221.0179186 , 1372.73614486]), array([5219.8756385, 1365.7919376]), array([5218.73348843, 1358.8477073 ]), array([5217.59146897, 1351.9034539 ]), array([5216.44958071, 1344.95917732]), array([5215.30782422, 1338.0148775 ]), array([5214.16620012, 1331.07055437]), array([5213.024709  , 1324.12620786]), array([5211.88335147, 1317.18183791]), array([5210.74212815, 1310.23744443]), array([5209.60103965, 1303.29302737]), array([5208.4600866 , 1296.34858663]), array([5207.31926965, 1289.40412216]), array([5206.17858942, 1282.45963387]), array([5205.03804657, 1275.51512169]), array([5203.89764175, 1268.57058554]), array([5202.75737563, 1261.62602534]), array([5201.61724887, 1254.68144102]), array([5200.47726214, 1247.73683248]), array([5199.33741614, 1240.79219965]), array([5198.19771156, 1233.84754245]), array([5197.05814908, 1226.90286079]), array([5195.91872942, 1219.95815459]), array([5194.7794533 , 1213.01342376]), array([5193.64032143, 1206.0686682 ]), array([5192.50133455, 1199.12388784]), array([5191.36249339, 1192.17908258]), array([5190.2237987 , 1185.23425233]), array([5189.08525124, 1178.28939699]), array([5187.94685178, 1171.34451647]), array([5186.80860108, 1164.39961068]), array([5185.67049993, 1157.45467952]), array([5184.53254912, 1150.50972288]), array([5183.39474945, 1143.56474067]), array([5182.25710174, 1136.61973279]), array([5181.1196068 , 1129.67469913]), array([5179.98226547, 1122.72963958]), array([5178.84507858, 1115.78455405]), array([5177.70804699, 1108.83944242]), array([5176.57117156, 1101.89430459]), array([5175.43445316, 1094.94914043]), array([5174.29789267, 1088.00394985]), array([5173.16149099, 1081.05873272]), array([5172.02524902, 1074.11348892]), array([5170.88916768, 1067.16821835]), array([5169.75324791, 1060.22292087]), array([5168.61749063, 1053.27759637]), array([5167.4818968 , 1046.33224472]), array([5166.34646739, 1039.3868658 ]), array([5165.21120338, 1032.44145948]), array([5164.07610576, 1025.49602563]), array([5162.94117553, 1018.55056412]), array([5161.8064137 , 1011.60507481]), array([5160.67182132, 1004.65955756]), array([5159.53739943,  997.71401225]), array([5158.40314908,  990.76843872]), array([5157.26907135,  983.82283684]), array([5156.13516734,  976.87720646]), array([5155.00143814,  969.93154744]), array([5153.86788487,  962.98585961]), array([5152.73450868,  956.04014284]), array([5151.60131072,  949.09439696]), array([5150.46829215,  942.14862183]), array([5149.33545416,  935.20281727]), array([5148.20279795,  928.25698312]), array([5147.07032476,  921.31111923]), array([5145.9380358 ,  914.36522542]), array([5144.80593235,  907.41930151]), array([5143.67401569,  900.47334735]), array([5142.5422871 ,  893.52736274]), array([5141.4107479 ,  886.58134752]), array([5140.27939943,  879.63530149]), array([5139.14824306,  872.68922447]), array([5138.01728015,  865.74311626]), array([5136.88651211,  858.79697669]), array([5135.75594036,  851.85080554]), array([5134.62556635,  844.90460262]), array([5133.49539154,  837.95836773]), array([5132.36541744,  831.01210066]), array([5131.23564556,  824.06580119]), array([5130.10607744,  817.11946911]), array([5128.97671465,  810.1731042 ]), array([5127.84755879,  803.22670625]), array([5126.71861148,  796.28027501]), array([5125.58987437,  789.33381026]), array([5124.46134914,  782.38731177]), array([5123.33303751,  775.44077929]), array([5122.20494121,  768.49421257]), array([5121.07706201,  761.54761137]), array([5119.94940172,  754.60097543]), array([5118.82196216,  747.65430449]), array([5117.69474522,  740.70759828]), array([5116.56775279,  733.76085654]), array([5115.44098681,  726.81407899]), array([5114.31444926,  719.86726535]), array([5113.18814215,  712.92041532]), array([5112.06206753,  705.97352861]), array([5110.9362275 ,  699.02660494]), array([5109.81062419,  692.07964398]), array([5108.68525978,  685.13264542]), array([5107.56013648,  678.18560895]), array([5106.43525656,  671.23853425]), array([5105.31062233,  664.29142097]), array([5104.18623615,  657.34426878]), array([5103.06210043,  650.39707732]), array([5101.93821762,  643.44984625]), array([5100.81459023,  636.50257519]), array([5099.69122084,  629.55526378]), array([5098.56811205,  622.60791164]), array([5097.44526656,  615.66051836]), array([5096.3226871 ,  608.71308357]), array([5095.20037647,  601.76560684]), array([5094.07833753,  594.81808776]), array([5092.95657322,  587.87052589]), array([5091.83508654,  580.92292081]), array([5090.71388055,  573.97527206]), array([5089.5929584 ,  567.02757918]), array([5088.4723233,  560.0798417]), array([5087.35197857,  553.13205913]), array([5086.23192757,  546.18423097]), array([5085.11217378,  539.23635671]), array([5083.99272075,  532.28843584]), array([5082.87357213,  525.34046781]), array([5081.75473165,  518.39245206]), array([5080.63620316,  511.44438804]), array([5079.51799061,  504.49627515]), array([5078.40009803,  497.5481128 ]), array([5077.28252961,  490.59990038]), array([5076.16528962,  483.65163723]), array([5075.04838245,  476.70332272]), array([5073.93181265,  469.75495617]), array([5072.81558487,  462.80653689]), array([5071.69970391,  455.85806415]), array([5070.58417472,  448.90953723]), array([5069.46900239,  441.96095536]), array([5068.35419218,  435.01231777]), array([5067.23974951,  428.06362364]), array([5066.12567998,  421.11487214]), array([5065.01198937,  414.1660624 ]), array([5063.89868364,  407.21719353]), array([5062.78576897,  400.26826461]), array([5061.67325174,  393.31927468]), array([5060.56113857,  386.37022276]), array([5059.44943629,  379.42110781]), array([5058.33815199,  372.47192878]), array([5057.22729304,  365.52268456]), array([5056.11686705,  358.573374  ]), array([5055.00688197,  351.62399593]), array([5053.897346 ,  344.6745491]), array([5052.78826773,  337.72503222]), array([5051.67965606,  330.77544396]), array([5050.57152027,  323.82578291]), array([5049.46387005,  316.87604764]), array([5048.35671548,  309.92623661]), array([5047.25006711,  302.97634825]), array([5046.14393597,  296.02638089]), array([5045.0383336 ,  289.07633279]), array([5043.93327209,  282.12620214]), array([5042.8287641 ,  275.17598702]), array([5041.72482296,  268.22568544]), array([5040.62146265,  261.27529529]), array([5039.51869792,  254.32481434]), array([5038.41654428,  247.37424026]), array([5037.31501814,  240.42357058]), array([5036.21413684,  233.47280269]), array([5035.11391872,  226.52193382]), array([5034.01438326,  219.57096104]), array([5032.91555115,  212.61988124]), array([5031.81744444,  205.6686911 ]), array([5030.72008662,  198.71738707]), array([5029.62350285,  191.76596539]), array([5028.52772006,  184.81442199]), array([5027.43276722,  177.86275251]), array([5026.33867552,  170.91095228]), array([5025.24547867,  163.95901621]), array([5024.15321322,  157.00693883]), array([5023.06191893,  150.05471416]), array([5021.97163919,  143.1023357 ]), array([5020.88242158,  136.1497963 ]), array([5019.79431851,  129.19708811]), array([5018.70738796,  122.24420245]), array([5017.62169444,  115.29112966]), array([5016.53731021,  108.33785895]), array([5015.45431667,  101.38437815]), array([5014.37280632,   94.43067347]), array([5013.29288516,   87.4767291 ]), array([5012.21467587,   80.52252679]), array([5011.13832213,   73.56804517]), array([5010.06399458,   66.61325894]), array([5008.99189914,   59.65813758]), array([5007.9222893 ,   52.70264361]), array([5006.8554849 ,   45.74672987]), array([5005.79190214,   38.79033511]), array([5004.73210521,   31.83337658]), array([5003.67690267,   24.87573612]), array([5002.62755138,   17.91723089]), array([5001.58627591,   10.95753922]), array([5.0005581e+03, 3.9959411e+00]), array([ 4.99956543e+03, -2.97070309e+00]), array([5.00060415e+03, 3.98935470e+00]), array([ 4.99952964e+03, -2.96537797e+00]), array([5.00065382e+03, 3.98165222e+00]), array([ 4.99949105e+03, -2.95915293e+00]), array([5.00070734e+03, 3.97265258e+00]), array([ 4.99944950e+03, -2.95188296e+00]), array([5.00076495e+03, 3.96214827e+00]), array([ 4.99940479e+03, -2.94340230e+00]), array([5.00082690e+03, 3.94990261e+00]), array([ 4.99935676e+03, -2.93352236e+00]), array([5.00089340e+03, 3.93564706e+00]), array([ 4.99930524e+03, -2.92202990e+00]), array([5.00096468e+03, 3.91907891e+00]), array([ 4.99925007e+03, -2.90868548e+00]), array([5.00104094e+03, 3.89985951e+00]), array([ 4.99919112e+03, -2.89322244e+00]), array([5.00112234e+03, 3.87761347e+00]), array([ 4.99912827e+03, -2.87534683e+00]), array([5.00120901e+03, 3.85192912e+00]), array([ 4.99906145e+03, -2.85473843e+00]), array([5.00130101e+03, 3.82236095e+00]), array([ 4.99899064e+03, -2.83105349e+00]), array([5.00139837e+03, 3.78843446e+00]), array([ 4.99891587e+03, -2.80392949e+00]), array([5.00150099e+03, 3.74965403e+00]), array([ 4.99883722e+03, -2.77299246e+00]), array([5.00160868e+03, 3.70551465e+00]), array([ 4.99875489e+03, -2.73786723e+00]), array([5.00172117e+03, 3.65551760e+00]), array([ 4.99866914e+03, -2.69819098e+00]), array([5.00183802e+03, 3.59919067e+00]), array([ 4.99858035e+03, -2.65362988e+00]), array([5.00195867e+03, 3.53611252e+00]), array([ 4.99848898e+03, -2.60389877e+00]), array([5.00208241e+03, 3.46594043e+00]), array([ 4.99839563e+03, -2.54878280e+00]), array([5.00220840e+03, 3.38844015e+00]), array([ 4.99830098e+03, -2.48815977e+00]), array([5.00233568e+03, 3.30351543e+00]), array([ 4.99820578e+03, -2.42202125e+00]), array([5.00246317e+03, 3.21123458e+00]), array([ 4.99811088e+03, -2.35049013e+00]), array([5.00258972e+03, 3.11185070e+00]), array([ 4.99801714e+03, -2.27383230e+00]), array([5.00271418e+03, 3.00581258e+00]), array([ 4.99792543e+03, -2.19246006e+00]), array([5.00283539e+03, 2.89376352e+00]), array([ 4.99783659e+03, -2.10692614e+00]), array([5.00295226e+03, 2.77652692e+00]), array([ 4.99775139e+03, -2.01790756e+00]), array([5.00306383e+03, 2.65507884e+00]), array([ 4.99767049e+03, -1.92618052e+00]), array([5.00316927e+03, 2.53050943e+00]), array([ 4.99759443e+03, -1.83258828e+00]), array([5.00326796e+03, 2.40397680e+00]), array([ 4.99752362e+03, -1.73800506e+00]), array([5.00335944e+03, 2.27665773e+00]), array([ 4.99745829e+03, -1.64329956e+00]), array([5.00344348e+03, 2.14969989e+00]), array([ 4.99739857e+03, -1.54930118e+00]), array([5.00352002e+03, 2.02417953e+00]), array([ 4.99734441e+03, -1.45677178e+00]), array([5.00358918e+03, 1.90106786e+00]), array([ 4.99729568e+03, -1.36638457e+00]), array([5.00365119e+03, 1.78120740e+00]), array([ 4.99725214e+03, -1.27871075e+00]), array([5.00370642e+03, 1.66529893e+00]), array([ 4.99721350e+03, -1.19421366e+00]), array([5.00375531e+03, 1.55389772e+00]), array([ 4.9971794e+03, -1.1132494e+00]), array([5.00379835e+03, 1.44741782e+00]), array([ 4.99714946e+03, -1.03607246e+00]), array([5.00383604e+03, 1.34614208e+00]), array([ 4.99712330e+03, -9.62844992e-01]), array([5.00386892e+03, 1.25023611e+00]), array([ 4.99710054e+03, -8.93648197e-01]), array([5.00389748e+03, 1.15976432e+00]), array([ 4.99708081e+03, -8.28494573e-01]), array([5.00392220e+03, 1.07470651e+00]), array([ 4.99706375e+03, -7.67340201e-01]), array([5.00394354e+03, 9.94974010e-01]), array([ 4.99704905e+03, -7.10096310e-01]), array([5.00396191e+03, 9.20424610e-01]), array([ 4.99703641e+03, -6.56639766e-01]), array([5.00397769e+03, 8.50875825e-01]), array([ 4.99702556e+03, -6.06822231e-01]), array([5.00399122e+03, 7.86116373e-01]), array([ 4.99701627e+03, -5.60477969e-01]), array([5.00400281e+03, 7.25915816e-01]), array([ 4.99700832e+03, -5.17430307e-01]), array([5.00401271e+03, 6.70032459e-01]), array([ 4.99700153e+03, -4.77496873e-01]), array([5.00402116e+03, 6.18219688e-01]), array([ 4.99699574e+03, -4.40493744e-01]), array([5.00402836e+03, 5.70230933e-01]), array([ 4.99699080e+03, -4.06238648e-01]), array([5.00403450e+03, 5.25823457e-01]), array([ 4.99698660e+03, -3.74553379e-01]), array([5.00403973e+03, 4.84761172e-01]), array([ 4.99698302e+03, -3.45265554e-01]), array([5.00404418e+03, 4.46816650e-01]), array([ 4.99697998e+03, -3.18209838e-01]), array([5.00404796e+03, 4.11772496e-01]), array([ 4.99697739e+03, -2.93228748e-01]), array([5.00405117e+03, 3.79422206e-01]), array([ 4.99697519e+03, -2.70173125e-01]), array([5.00405391e+03, 3.49570625e-01]), array([ 4.99697333e+03, -2.48902346e-01]), array([5.00405622e+03, 3.22034109e-01]), array([ 4.99697174e+03, -2.29284342e-01]), array([5.00405819e+03, 2.96640451e-01]), array([ 4.99697040e+03, -2.11195471e-01]), array([5.00405987e+03, 2.73228638e-01]), array([ 4.99696925e+03, -1.94520291e-01]), array([5.00406128e+03, 2.51648494e-01]), array([ 4.99696828e+03, -1.79151249e-01]), array([5.00406249e+03, 2.31760231e-01]), array([ 4.99696746e+03, -1.64988329e-01]), array([5.00406351e+03, 2.13433957e-01]), array([ 4.99696676e+03, -1.51938671e-01]), array([5.00406438e+03, 1.96549147e-01]), array([ 4.99696617e+03, -1.39916171e-01]), array([5.00406511e+03, 1.80994105e-01]), array([ 4.99696567e+03, -1.28841075e-01]), array([5.00406573e+03, 1.66665417e-01]), array([ 4.99696524e+03, -1.18639581e-01]), array([5.00406626e+03, 1.53467426e-01]), array([ 4.99696488e+03, -1.09243448e-01]), array([5.00406671e+03, 1.41311708e-01]), array([ 4.99696458e+03, -1.00589619e-01]), array([5.00406709e+03, 1.30116580e-01]), array([ 4.99696432e+03, -9.26198601e-02]), array([5.00406741e+03, 1.19806622e-01]), array([ 4.99696410e+03, -8.52804143e-02]), array([5.00406769e+03, 1.10312226e-01])]

Plot the gradient descent process

%%capture

fig,ax = plt.subplots()
ax.plot(5000,0,'bx')
ax.annotate('target',(5000,0))
ax.grid()
xsa = np.array(xs[::10])
def animate(ii,):
    ax.plot(xsa[ii,0],xsa[ii,1],'ro-')

ani = animation.FuncAnimation(fig,animate,repeat=False,frames=len(xsa)-1,interval=100)

HTML(ani.to_jshtml())

Once Loop Reflect

Figure 8: Gradient descent for vector problem

Auto-diff Forward vs Reverse mode

Everything we talked about is called “Forward mode”. Computing derivatives using the chain rule [9].

(Reverse mode is a.k.a. Adjoint mode)

References

[1]

R. Kwiatkowski, “Gradient descent algorithm — a deep dive.” Online, accessed 14 Nov 2024., May 2021. Available: https://towardsdatascience.com/gradient-descent-algorithm-a-deep-dive-cf04e8115f21

[2]

“Dual numbers & automatic differentiation.” Online. Accessed 2024-11-17, Dec. 2014. Available: https://blog.demofox.org/2014/12/30/dual-numbers-automatic-differentiation/

[3]

S. Ritchie, “Dual numbers and automatic differentiation.” Online. Accessed 2024-11-18., Jan. 2021. Available: https://samritchie.io/dual-numbers-and-automatic-differentiation/

[4]

“Dual numbers for first order sensitivity analysis.” Online. Accessed 2024-11-17, 2023. Available: https://ceid.utsa.edu/HYPAD/wp-content/uploads/sites/50/2023/04/3DualNumbers-12.pdf

[5]

F. Peñuñuri, O. Carvente, R. Peón, and C. A. Cruz-Villar, “Dual numbers for algorithmic differentiation.” Online. Accessed 2024-11-17, 2019. Available: https://www.redalyc.org/journal/467/46761359006/html/

[6]

“Autograd.” Online. Accessed 2024-11-17., 2024. Available: https://github.com/HIPS/autograd

[7]

P. Nobel, “Auto_diff: An automatic differentiation package for python.” in Proceedings of the 2020 spring simulation conference (SpringSim ’20). Society for computer simulation international, san diego, CA, USA., Society for Computer Simulation International, 2020, pp. 1–10. Available: https://github.com/PTNobel/autodiff

[8]

“Autodiff.” Online. Accessed 2024-11-18, 2024. Available: https://github.com/ryanirl/autodiff

[9]

“\(xa\partial\) the evolution of automatic differentiation.” Online. Accessed 2024-11-18, 2024. Available: https://auto-differentiation.github.io/