# Example notebook to show how adaptive frequency oscillators (AFOs) work

In the first part of the notebook, we show how to implement an AFO in native python. Then we show how to use the c++ library to make faster simulation in the second part of the notebook. An example with an input with a time varying frequency is shown last


```python
# standard imports
%matplotlib widget
import numpy as np
import matplotlib.pylab as plt
```

## Adaptive Frequency Oscillator in Python
We wish to implement an AFO whose equation is
$$\begin{align} \dot{\phi} & = \lambda \omega - K \sin(\phi) F(t) \\ \dot{\omega} &= - K \sin(\phi) F(t) \end{align}$$
where $$F(t) = \sin(\omega_F t)$$


```python
def integrate_afo(phi0, omega0, omegaF, t_init, t_end, K=10000, lamb=1., save_dt=0.001, dt=0.00001):
    """ this function integrates an AFO from t_init to t_end 
    it starts with initial conditions phi0 and omega0
    returns: numpy vectors t, phi and omega
    """
    # hoe many integration steps in the internal loop
    internal_step = int(round(save_dt/dt))

    # how many steps till t_end
    num_steps = int(round((t_end - t_init)/save_dt))

    # we preallocate memory
    t = np.zeros(num_steps+1)
    phi = np.zeros_like(t)
    omega = np.zeros_like(t)
    
    # we set initial conditions
    t[0] = t_init
    omega[0] = omega0
    phi[0] = phi0
    
    # our temp variables
    phi_temp = phi0
    omega_temp = omega0
    t_temp = t_init
    
    # the main integration loop
    for i in range(num_steps):
        #internal integration loop
        for j in range(internal_step):
            pert = -K * np.sin(omegaF * t_temp) * np.sin(phi_temp)
            phi_temp += (lamb * omega_temp + pert)*dt
            omega_temp += pert*dt
            t_temp += dt
        
        # save data
        t[i+1] = t_temp
        phi[i+1] = phi_temp
        omega[i+1] = omega_temp
    
    return t, phi, omega
```


```python
# we define the coupling stength, input frequency, lambda
K =  1000.
omegaF = 30
lamb = 1.

# we define the time resolution (save_dt) and the internal integration step dt
dt = 0.00001
save_dt = 0.001

# duration of the integration
t_end = 20.
t_init = 0.

phi0 = 0
omega0 = 100

# we integrate
t, phi, omega = integrate_afo(phi0, omega0, omegaF, t_init, t_end, K=10000, lamb=1., save_dt=0.001, dt=0.00001)

# now we plot the results
plt.figure()
plt.subplot(2,1,1)
plt.plot(t, phi)
plt.ylabel(r'$\phi$')
plt.subplot(2,1,2)
plt.plot(t, lamb * omega)
plt.plot([t_init, t_end],[omegaF, omegaF], '--k')
plt.ylabel(r'$\omega$')
plt.xlabel('Time [s]')
```


    Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …





    Text(0.5, 0, 'Time [s]')



## Adaptive Frequency Oscillator using C++ library and python bindings
Now we do the same using the c++ library (it needs to be compiled and installed - cf. README)


```python
import pyafos # we import the library
```


```python
# we define the coupling stength, input frequency, lambda
K =  1000.
omegaF = 30
lamb = 1.

# we define the time resolution (save_dt) and the internal integration step dt
dt = 0.00001
save_dt = 0.001

# duration of the integration
t_end = 20.
t_init = 0.

phi0 = 0
omega0 = 100

#run an integration
oscill = pyafos.PhaseAFO()
oscill.initialize(K, lamb)
oscill.input().vec_of_sines(np.array([omegaF]),np.array([1.0]),np.array([0.0]))
pyafos.integrate(oscill, t_init,t_end,np.array([phi0,omega0]),dt,save_dt)

#get the data to be plotted    
t = oscill.t()
phi = oscill.y()[0,:]
omega = oscill.y()[1,:]

# now we plot the results
plt.figure()
plt.subplot(2,1,1)
plt.plot(t, phi)
plt.ylabel(r'$\phi$')
plt.subplot(2,1,2)
plt.plot(t, lamb * omega)
plt.plot([t_init, t_end],[omegaF, omegaF], '--k')
plt.ylabel(r'$\omega$')
plt.xlabel('Time [s]')
```


    Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …





    Text(0.5, 0, 'Time [s]')



## Tracking a time varying frequency

We want to track a frequency that is changing over time $F(t) = \sin\left(\omega_F t + \frac{1}{\omega_C}\sin(\omega_C t)\right)$


```python
# we define the coupling stength, input frequency, lambda
K =  1000.
omegaF = 50
lamb = 1.
omegaC = 0.5

# we define the time resolution (save_dt) and the internal integration step dt
dt = 0.00001
save_dt = 0.001

# duration of the integration
t_end = 50.
t_init = 0.

phi0 = 0
omega0 = 40

#run an integration
oscill = pyafos.PhaseAFO()
oscill.initialize(K, lamb)
oscill.input().frequency_changing_sine(omegaF, omegaC)
pyafos.integrate(oscill, t_init,t_end,np.array([phi0,omega0]),dt,save_dt)

#get the data to be plotted    
t = oscill.t()
phi = oscill.y()[0,:]
omega = oscill.y()[1,:]

# now we plot the results
plt.figure()
plt.subplot(2,1,1)
plt.plot(t, phi)
plt.ylabel(r'$\phi$')
plt.subplot(2,1,2)
plt.plot(t, lamb * omega)
plt.plot(t, omegaF + np.sin(omegaC*t))
plt.ylim([40,60])
plt.ylabel(r'$\omega$')
plt.xlabel('Time [s]')
```


    Canvas(toolbar=Toolbar(toolitems=[('Home', 'Reset original view', 'home', 'home'), ('Back', 'Back to previous …





    Text(0.5, 0, 'Time [s]')




```python

```

987c5dbe84d361d0f44cfa4707d63d49862f2433
ipynb
Jupyter Notebook
python/adaptive frequency oscillators example.ipynb
righetti/AFOs
79e23f6ef278621e66eb775ff2fb5156eaac1261
BSD-3-Clause
2021-03-31T06:41:31.000Z
2021-12-21T05:24:55.000Z
python/adaptive frequency oscillators example.ipynb
righetti/AFOs
79e23f6ef278621e66eb775ff2fb5156eaac1261
BSD-3-Clause
python/adaptive frequency oscillators example.ipynb
righetti/AFOs
79e23f6ef278621e66eb775ff2fb5156eaac1261
BSD-3-Clause
Qwen/Qwen-72B
1. YES
2. YES