Time Domain Simulations

Author: Kenny Moc (kmoc@ualberta.ca)

Overview

This notebook demonstrates how to perform a time-domain simulation of the 40m PRMI. The PRM and BS optics are subject to noise forces measured at the 40m. The ITMs are assumed to be stationary.

Key Points to Note

  • Currently, the pendulums start at rest. We plan to add randomized velocity. (WIP)

  • Noise forces repeat every 16 seconds due to nperseg being set to 16 seconds in the Welch method.

  • It looks like the pendulums are not as velocity damped as they should be. Example on how to change velocity damping is shown below.

Validity of Time Domain models

../_images/time_domain.drawio.png

Finesse is frequency domain model. This time-domain simulation invokes a frequency domain model to compute the optical steady state conditions. We can then compute DC electrical control signals from these optical fields and use these to advance the mechanical model states. The mechanical model can then be used to set the mirror positions for the next steady state optical calculation. See above diagram for more details.

For this to be valid, the lowest frequencies in our optical model must be much faster than the aliasing frequency of our mechanical model. In this model, the update rate is 16 kHz. Therefore, any optical effects with frequencies less than about 100 kHz, such as optomechanical interactions, can’t be trivially modelled. For the 40m, this is true when locking at low power.

Imports

from finesse.detectors import PowerDetector
import finesse_40m
from finesse_40m.time_simulation.prmi import Simulation
from finesse_40m.time_simulation.prmi.factory import make_model

import numpy as np
import matplotlib.pyplot as plt

Vanilla Simulation

Creates the model, randomizes the initial PRCL and MICH values, then creates the simulation

model = make_model()

model.PRCL.DC = np.random.randint(-360, 360)
model.MICH.DC = np.random.randint(-360, 360)

sim = Simulation(model,runtime=1) # dummy simulation
print(f"Sample rate is fixed at {1e-3*sim.fs:.3f} kHz "
      f"and sample time is: {1e-3*sim.dt:.3f} ms")
sim = Simulation(model, runtime=.2*sim.fs) # Overwrite with desired parameters
print(f"Total runtime (s): {sim.duration}")

Writing data to /root/finesse-data/finesse-40m

Writing data to /root/finesse-data/finesse-40m

Sample rate is fixed at 16.384 kHz and sample time is: 0.000 ms
Total runtime (s): 0.2

results = sim.run()

A dedicated list exists with all the PDs specific to this class, or you can access them by interrogating the model.

from finesse_40m.time_simulation.prmi.factory import pd_list, pdh_list
print(pd_list + pdh_list)

# Note one can also do
# Get pds
print(
    [detector for detector in model.detectors if isinstance(detector, PowerDetector)]
)

['POP', 'AS_pd', 'REFL', 'Pinx', 'Piny', 'POP11_I', 'POP11_Q', 'POP22_I', 'POP22_Q', 'POP55_I', 'POP55_Q', 'POP110_I', 'POP110_Q', 'REFL11_I', 'REFL11_Q', 'REFL22_I', 'REFL22_Q', 'REFL55_I', 'REFL55_Q', 'REFL110_I', 'REFL110_Q', 'REFL33_I', 'REFL33_Q', 'REFL165_I', 'REFL165_Q', 'AS11_I', 'AS11_Q', 'AS22_I', 'AS22_Q', 'AS55_I', 'AS55_Q', 'AS110_I', 'AS110_Q']
[<'Px' @ 0x7f58ba8ea390 (PowerDetector)>, <'Pinx' @ 0x7f58bb2317c0 (PowerDetector)>, <'Py' @ 0x7f58ba97fb90 (PowerDetector)>, <'Piny' @ 0x7f58ba989c40 (PowerDetector)>, <'Pprc' @ 0x7f58ba8b6810 (PowerDetector)>, <'Psrc' @ 0x7f58ba8b4710 (PowerDetector)>, <'Pin' @ 0x7f58ba17c740 (PowerDetector)>, <'PreflPRM' @ 0x7f58ba17e870 (PowerDetector)>, <'Prefl' @ 0x7f58ba17eed0 (PowerDetector)>, <'Ppop' @ 0x7f58ba17e300 (PowerDetector)>, <'Pas' @ 0x7f58ba17e480 (PowerDetector)>, <'POP' @ 0x7f58ba0663f0 (PowerDetector)>, <'REFL' @ 0x7f58ba57d310 (PowerDetector)>, <'AS_pd' @ 0x7f58ba6aa480 (PowerDetector)>]

We can now plot the results.

def doplot(sim, results):
    fig, ax = plt.subplots()
    t = sim.t
    ax.plot(t, results['AS_pd'], label='AS')
    ax.plot(t, results['REFL'], label='REFL')
    ax.set_ylabel('Normalized Power')
    ax.legend(title='Port')
    ax.set_xlabel('Time (s)')
    return fig, ax
doplot(sim, results)

(<Figure size 640x480 with 1 Axes>,
 <Axes: xlabel='Time (s)', ylabel='Normalized Power'>)
../_images/time_domain_4_1.png

Adding a controller

To add a controller, we need to overwrite the controller method in the simulation environment.

NOTE: The done flag is not yet implemented. The simulation will run for the full runtime regardless of the done flag. However, it still needs to be returned by the function

# Controllers can read any signal from the pdh_list and pd_list and apply forces to the pendulums accordingly.
def kenny_new_controller(observation_dict):
    """A simple integral controller"""
    done = False
    if observation_dict['REFL'] > 0:
        apply_prm_force = 1
        apply_mich_force = -1
    if observation_dict['AS_pd'] < 0.1:
        done = True
    return apply_prm_force, apply_mich_force, done

sim.env.controller = kenny_new_controller
results = sim.run()
doplot(sim, results)

(<Figure size 640x480 with 1 Axes>,
 <Axes: xlabel='Time (s)', ylabel='Normalized Power'>)
../_images/time_domain_5_2.png

Changing Physical Attributes of Pendulum

# Changes the beamsplitter pendulum damping coefficient
sim.env.bs_pendulum.damping_coefficient = 1.2

# Changes local gravity for prm pendulum
sim.env.prm_pendulum.g = 10

results = sim.run()
doplot(sim, results)

(<Figure size 640x480 with 1 Axes>,
 <Axes: xlabel='Time (s)', ylabel='Normalized Power'>)
../_images/time_domain_6_2.png