Commit 4bbac678 authored by Phil Jones's avatar Phil Jones

Merge branch 'master' of git.ligo.org:finesse/finesse3 into apidoc-fixes

parents 29bf2d8b 1a654c8f
......@@ -141,19 +141,10 @@ docs/html:
pages:
stage: deploy
needs:
# Everything else must have passed. Grab artifacts only from the doc build job.
- job: test/3.8/functional
artifacts: false
- job: test/3.9/functional
artifacts: false
- job: test/3.8/validation
artifacts: false
- job: test/3.8/validation
artifacts: false
- job: test/coverage
artifacts: false
- job: docs/html
artifacts: true
only:
refs:
- master
script:
- mv docs/build/html public
artifacts:
......
......@@ -318,3 +318,97 @@ class FreeMass(Connector):
ws.owner_id, ws.signal.connections.F_to_PITCH_idx, 0, 0,
) as mat:
mat[:] = -1 / (ws.values.I_pitch * (2*PI*f)**2)
class PendulumMassWorkspace(ConnectorWorkspace):
pass
@model_parameter("mass", "Mass", 1, 0, "kg")
@model_parameter("Qz", "Qz", 1, 0, "")
@model_parameter("fz", "fz", 1, 0, "Hz")
@model_parameter("I_pitch", "Moment of inertia (pitch)", 1, 0, "kg·m^2")
@model_parameter("Qyaw", "Qyaw", 1, 0, "")
@model_parameter("fyaw", "fyaw", 1, 0, "Hz")
@model_parameter("I_yaw", "Moment of inertia (yaw)", 1, 0, "kg·m^2")
@model_parameter("Qpitch", "Qpitch", 1, 0, "")
@model_parameter("fpitch", "fpitch", 1, 0, "Hz")
class Pendulum(Connector):
"""Simple pendulum suspension of an object.
The object being suspended must have a mechanical port with
nodes z, pitch, and yaw and forces F_z, F_pitch, and F_yaw.
"""
def __init__(self, name, mech_port, mass=np.inf, Qz=1000, fz=1, I_yaw=np.inf, Qyaw=1000, fyaw=1, I_pitch=np.inf, Qpitch=1000, fpitch=1):
super().__init__(name)
self.mass = mass
self.Qz = Qz
self.fz = fz
self.Qyaw = Qyaw
self.fyaw = fyaw
self.Qpitch = Qpitch
self.fpitch = fpitch
self.I_yaw = I_yaw
self.I_pitch = I_pitch
# Add motion and force nodes to mech port.
# Here we duplicate the already created mechanical
# nodes in some other connector element
self._add_port("mech", NodeType.MECHANICAL)
self.mech._add_node("z", None, mech_port.z)
self.mech._add_node("yaw", None, mech_port.yaw)
self.mech._add_node("pitch", None, mech_port.pitch)
self.mech._add_node("F_z", None, mech_port.F_z)
self.mech._add_node("F_yaw", None, mech_port.F_yaw)
self.mech._add_node("F_pitch", None, mech_port.F_pitch)
# We just have direct coupling between DOF, no cross-couplings
self._register_node_coupling(
"F_to_Z", self.mech.F_z, self.mech.z,
enabled_check=lambda: self.mass < np.inf and not self.mass.is_changing
)
self._register_node_coupling(
"F_to_YAW", self.mech.F_yaw, self.mech.yaw,
enabled_check=lambda: self.I_yaw < np.inf and not self.I_yaw.is_changing
)
self._register_node_coupling(
"F_to_PITCH", self.mech.F_pitch, self.mech.pitch,
enabled_check=lambda: self.I_pitch < np.inf and not self.I_pitch.is_changing
)
def _get_workspace(self, sim):
if sim.signal:
refill = sim.model.fsig.f.is_changing or any(p.is_changing for p in self.parameters)
ws = PendulumMassWorkspace(self, sim, refill, refill)
ws.signal.set_fill_function(self.fill)
return ws
else:
return None
def fill(self, ws):
cdef:
complex_t s = 2j* PI * ws.sim.model_data.fsig
complex_t num = 1
complex_t den = 1
double omega0
if ws.signal.connections.F_to_Z_idx >= 0:
with ws.sim.signal.component_edge_fill3(
ws.owner_id, ws.signal.connections.F_to_Z_idx, 0, 0,
) as mat:
omega0 = 2 * PI * ws.values.fz
mat[:] = 1 / ws.values.mass * 1/(s**2 + s * omega0/ws.values.Qz + omega0**2)
if ws.signal.connections.F_to_YAW_idx >= 0:
with ws.sim.signal.component_edge_fill3(
ws.owner_id, ws.signal.connections.F_to_YAW_idx, 0, 0,
) as mat:
omega0 = 2 * PI * ws.values.fyaw
mat[:] = 1 / ws.values.I_yaw * 1/(s**2 + s * omega0/ws.values.Qyaw + omega0**2)
if ws.signal.connections.F_to_PITCH_idx >= 0:
with ws.sim.signal.component_edge_fill3(
ws.owner_id, ws.signal.connections.F_to_PITCH_idx, 0, 0,
) as mat:
omega0 = 2 * PI * ws.values.fpitch
mat[:] = 1 / ws.values.I_pitch * 1/(s**2 + s * omega0/ws.values.Qpitch + omega0**2)
\ No newline at end of file
......@@ -442,15 +442,15 @@ cdef void single_yaw_mechanical_frequency_signal_calc (
wx = beam_size(q_P1o.qx, ws.nr1, ws.sim.model_data.lambda0)
# fill_prop_za as off-diagonal -1 is already included in the carrier connection
(<SubCCSView>conn.yaw_P1o[yaw_freq_idx][freq.index]).fill_prop_za_zm (
# factor of 2 because misalignment is 2 * x/ybeta
(<SubCCSView>car_conn.P1i_P1o[carrier_index]), 0, 2*ws.nr1 * wx * a_2_o_factor, ws.K_yaw_sig.data_view, False
# factor of 2 because misalignment is 2 * x/ybeta, but 0.5 factor from upper/lower SB gain
(<SubCCSView>car_conn.P1i_P1o[carrier_index]), 0, 2/2*ws.nr1 * wx * a_2_o_factor, ws.K_yaw_sig.data_view, False
)
if conn.yaw_P2o[yaw_freq_idx][freq.index]:
wx = beam_size(q_P2o.qx, ws.nr2, ws.sim.model_data.lambda0)
# fill_prop_za as off-diagonal -1 is already included in the carrier connection
(<SubCCSView>conn.yaw_P2o[yaw_freq_idx][freq.index]).fill_prop_za_zm (
# factor of 2 because misalignment is 2 * x/ybeta
(<SubCCSView>car_conn.P2i_P2o[carrier_index]), 0, 2*ws.nr2 * wx * a_2_o_factor, ws.K_yaw_sig.data_view, False
# factor of 2 because misalignment is 2 * x/ybeta, but 0.5 factor from upper/lower SB gain
(<SubCCSView>car_conn.P2i_P2o[carrier_index]), 0, 2/2*ws.nr2 * wx * a_2_o_factor, ws.K_yaw_sig.data_view, False
)
# TODO ddb - need to add in transmission signal effects
......@@ -524,8 +524,8 @@ cdef void single_pitch_mechanical_frequency_signal_calc (
wy = beam_size(q_P1o.qy, ws.nr1, ws.sim.model_data.lambda0)
# fill_prop_za as off-diagonal -1 is already included in the carrier connection
(<SubCCSView>conn.pitch_P1o[pitch_freq_idx][freq.index]).fill_prop_za_zm (
# factor of 2 because misalignment is 2 * x/ybeta
(<SubCCSView>car_conn.P1i_P1o[carrier_index]), 0, 2 * ws.nr1 * wy * a_2_o_factor, ws.K_pitch_sig.data_view, False
# factor of 2 because misalignment is 2 * x/ybeta, but 0.5 factor from upper/lower SB gain
(<SubCCSView>car_conn.P1i_P1o[carrier_index]), 0, -2/2 * ws.nr1 * wy * a_2_o_factor, ws.K_pitch_sig.data_view, False
)
if conn.pitch_P2o[pitch_freq_idx][freq.index]:
......@@ -533,9 +533,9 @@ cdef void single_pitch_mechanical_frequency_signal_calc (
# fill_prop_za as off-diagonal -1 is already included in the carrier connection
# Extra minus sign here because of coordinate system from back side of of tilted mirror
(<SubCCSView>conn.pitch_P2o[pitch_freq_idx][freq.index]).fill_prop_za_zm (
# factor of 2 because misalignment is 2 * x/ybeta
# factor of 2 because misalignment is 2 * x/ybeta, but 0.5 factor from upper/lower SB gain
# minus because side 2 pitch sends beam upwards
(<SubCCSView>car_conn.P2i_P2o[carrier_index]), 0, -2 * ws.nr2 * wy * a_2_o_factor, ws.K_pitch_sig.data_view, False
(<SubCCSView>car_conn.P2i_P2o[carrier_index]), 0, 2/2 * ws.nr2 * wy * a_2_o_factor, ws.K_pitch_sig.data_view, False
)
# -------------------------------------------------
# Optical to mechanical connections
......
......@@ -77,7 +77,7 @@ class ModelElement:
return "<'{}' @ {} ({})>".format(
self.name, hex(id(self)), self.__class__.__name__
)
def __deepcopy__(self, memo):
new = object.__new__(type(self))
memo[id(self)] = new
......@@ -450,8 +450,8 @@ cdef class ElementWorkspace:
swaps = []
for i, p in enumerate(symbolic_params):
for dependent_p in p.value.parameters():
if dependent_p in symbolic_params:
j = symbolic_params.index(dependent_p)
if dependent_p.parameter in symbolic_params:
j = symbolic_params.index(dependent_p.parameter)
if j > i:
swaps.append((i, j))
for i, j in swaps:
......
......@@ -1272,6 +1272,7 @@ class KatSpec(metaclass=Singleton):
("nothing", components.Nothing),
# Mechanics.
("free_mass", mechanical.FreeMass),
("pendulum", mechanical.Pendulum),
("ligo_triple", ligo.LIGOTripleSuspension),
# Lock.
("lock", locks.Lock),
......
......@@ -123,10 +123,18 @@ cdef class CarrierSignalMatrixSimulation:
public list variable_workspaces
readonly list gouy_phase_workspaces
# Tracing stuff
### Tracing stuff ###
public dict cavity_workspaces
NodeBeamParam* trace
# The TraceForest of geometrically changing branches. This is an
# empty forest for any simulation in which geometric parameters
# are not changing.
readonly TraceForest trace_forest
# A dict of {<tuple of newly unstable cavities> : <contingent TraceForest>}
# required for when a scan results in a geometrically changing cavity becoming
# unstable -> invalidating self.trace_forest temporarily for that data point
dict contingent_trace_forests
bint needs_reflag_changing_q # Used when exiting from unstable cavity regions
bint retrace
# List of workspaces for components which scatter modes
......@@ -151,10 +159,14 @@ cdef class CarrierSignalMatrixSimulation:
cpdef run_signal(self)
# Methods to construct the changing TraceForest for the simulation
cdef void _determine_changing_beam_params(self)
cdef void _setup_trace_forest(self)
cdef void _setup_single_trace_tree(self, TraceTree tree)
cdef void _determine_changing_beam_params(self, TraceForest forest=?, bint set_tree_node_ids=?)
cdef void _setup_trace_forest(self, TraceForest forest=?, bint set_tree_node_ids=?)
cdef void _setup_single_trace_tree(self, TraceTree tree, bint set_tree_node_ids=?)
# Find the newly unstable cavity instances from the changing forest
cdef tuple _find_new_unstable_cavities(self)
cdef TraceForest _initialise_contingent_forest(self, tuple unstable_cavities)
# Perform the beam trace on the changing TraceForest
cdef void _propagate_trace(self, TraceTree tree)
cdef void _propagate_trace(self, TraceTree tree, bint symmetric)
cdef bint trace_beam(self)
......@@ -285,12 +285,12 @@ cdef class MatrixSystemSolver:
nio = []
for _ in owner._registered_connections[name]:
nio.append(owner.nodes[_])
enabled_check = owner._enabled_checks.get(name, None)
if enabled_check:
enabled_check = enabled_check()
else:
enabled_check = True
enabled_check = True
nio = tuple(nio)
......@@ -865,6 +865,8 @@ cdef class CarrierSignalMatrixSimulation:
self.model = model
self.name = name
self.trace = NULL
self.contingent_trace_forests = {}
self.needs_reflag_changing_q = False
self.do_matrix_solving = needs_matrix
def build(self):
......@@ -1006,14 +1008,16 @@ cdef class CarrierSignalMatrixSimulation:
# and TEM Gouy phases at lasers
return self.set_gouy_phases()
cdef void _determine_changing_beam_params(self):
cdef void _determine_changing_beam_params(
self, TraceForest forest=None, bint set_tree_node_ids=True,
):
cdef:
Py_ssize_t i
Py_ssize_t num_nodes = self.carrier.num_nodes
ConnectorWorkspace ws
KnmConnectorWorkspace kws
# re-set all beam parameter changing flags to false initially
# Re-set all beam parameter changing flags to false initially
for i in range(num_nodes):
self.trace[i].is_changing = False
......@@ -1021,33 +1025,42 @@ cdef class CarrierSignalMatrixSimulation:
LOGGER.info("Flagging changing beam parameters.")
# Prepare the forest for simulation by setting all the node_id attributes
# and flag the corresponding self.trace entries as changing
self._setup_trace_forest()
self._setup_trace_forest(forest, set_tree_node_ids)
# now tell each knm workspace whether it is changing or not
# Now tell each knm workspace whether it is changing or not
# so that only changing scattering matrices get recomputed
# from here on
for kws in self.to_scatter_matrix_compute:
kws.flag_changing_knm_workspaces()
cdef void _setup_trace_forest(self):
cdef void _setup_trace_forest(self, TraceForest forest=None, bint set_tree_node_ids=True):
cdef:
Py_ssize_t tree_idx
TraceTree tree
for tree_idx in range(self.trace_forest.N_trees):
tree = self.trace_forest.forest[tree_idx]
self._setup_single_trace_tree(tree)
if forest is None:
forest = self.trace_forest
for tree_idx in range(forest.N_trees):
tree = forest.forest[tree_idx]
self._setup_single_trace_tree(tree, set_tree_node_ids)
if self.trace_forest.N_trees:
LOGGER.info("Determined changing trace trees:%s", self.trace_forest.draw())
LOGGER.info(
"Determined changing trace trees:%s", self.trace_forest.draw_by_order()
)
cdef void _setup_single_trace_tree(self, TraceTree tree):
cdef void _setup_single_trace_tree(self, TraceTree tree, bint set_tree_node_ids=True):
cdef:
TraceTree ltree = tree.left
TraceTree rtree = tree.right
tree.node_id = self.carrier.node_id(tree.node)
tree.opp_node_id = self.carrier.node_id(tree.node.opposite)
# Only ever need to do this once, so avoid repeating when reflagging
# changing beam params after exiting unstable cavity regions
if set_tree_node_ids:
tree.node_id = self.carrier.node_id(tree.node)
tree.opp_node_id = self.carrier.node_id(tree.node.opposite)
self.trace[tree.node_id].is_changing = True
self.trace[tree.opp_node_id].is_changing = True
......@@ -1056,8 +1069,75 @@ cdef class CarrierSignalMatrixSimulation:
if rtree is not None:
self._setup_single_trace_tree(rtree)
cdef tuple _find_new_unstable_cavities(self):
cdef:
Py_ssize_t tree_idx
TraceTree tree
CavityWorkspace cav_ws
bint source_is_cav
double gx, gy
list ch_unstable_cavities = []
for tree_idx in range(self.trace_forest.N_trees):
tree = self.trace_forest.forest[tree_idx]
if tree.is_source:
cav_ws = self.cavity_workspaces.get(tree.dependency)
source_is_cav = cav_ws is not None
if source_is_cav: # Tree is an internal cavity tree
# The geometrically changing cavity has become unstable
# so inform that this is the case
if not cav_ws.is_stable:
ch_unstable_cavities.append(cav_ws.owner)
gx = cav_ws.gx
gy = cav_ws.gy
if float_eq(gx, gy):
LOGGER.warning(
"Cavity %s is unstable with g = %s",
tree.dependency.name,
gx,
)
else:
LOGGER.warning(
"Cavity %s is unstable with gx = %s, gy = %s",
tree.dependency.name,
gx,
gy,
)
return tuple(ch_unstable_cavities)
cdef TraceForest _initialise_contingent_forest(self, tuple unstable_cavities):
cdef TraceForest contingent_forest = TraceForest(self.model)
# FIXME (sjr) Currently not taking the model.sim_trace_config into account,
# need to replace order accordingly if there are entries in this
# config dict
cdef list order = self.model.trace_order
for uc in unstable_cavities:
order.remove(uc)
# If there are no dependencies left after disabling the
# unstable cavities then a beam trace cannot be performed
# at this data point so no forest can be planted
if not order:
LOGGER.warning(
"Cannot build a contingent trace forest as the simulation "
"is in a regime with no stable cavities nor Gauss objects."
)
return None
contingent_forest.symmetric = self.trace_forest.symmetric
contingent_forest.plant(order)
self._determine_changing_beam_params(contingent_forest)
return contingent_forest
@cython.initializedcheck(False)
cdef void _propagate_trace(self, TraceTree tree):
cdef void _propagate_trace(self, TraceTree tree, bint symmetric):
cdef:
TraceTree ltree = tree.left
TraceTree rtree = tree.right
......@@ -1073,7 +1153,7 @@ cdef class CarrierSignalMatrixSimulation:
# other dependency trees. Note these are only performed
# on the left tree; see TraceForest._add_backwards_nonsymm_trees
# for details.
if self.trace_forest.symmetric or (not tree.do_nonsymm_reverse and not tree.do_inv_transform):
if symmetric or (not tree.do_nonsymm_reverse and not tree.do_inv_transform):
qx2 = transform_q(tree.left_abcd_x, qx1, tree.nr, ltree.nr)
qy2 = transform_q(tree.left_abcd_y, qy1, tree.nr, ltree.nr)
elif tree.do_inv_transform:
......@@ -1089,11 +1169,11 @@ cdef class CarrierSignalMatrixSimulation:
self.trace[ltree.node_id].qx = qx2
self.trace[ltree.node_id].qy = qy2
if self.trace_forest.symmetric:
if symmetric:
self.trace[ltree.opp_node_id].qx = -conj(qx2)
self.trace[ltree.opp_node_id].qy = -conj(qy2)
self._propagate_trace(ltree)
self._propagate_trace(ltree, symmetric)
if rtree is not None:
qx2 = transform_q(tree.right_abcd_x, qx1, tree.nr, rtree.nr)
......@@ -1101,11 +1181,11 @@ cdef class CarrierSignalMatrixSimulation:
self.trace[rtree.node_id].qx = qx2
self.trace[rtree.node_id].qy = qy2
if self.trace_forest.symmetric:
if symmetric:
self.trace[rtree.opp_node_id].qx = -conj(qx2)
self.trace[rtree.opp_node_id].qy = -conj(qy2)
self._propagate_trace(rtree)
self._propagate_trace(rtree, symmetric)
cdef bint trace_beam(self):
cdef:
......@@ -1114,72 +1194,116 @@ cdef class CarrierSignalMatrixSimulation:
CavityWorkspace cav_ws
bint source_is_cav
double gx, gy
complex_t qx_src, qy_src
Py_ssize_t source_node_id, source_opp_node_id
bint valid = not self.trace_forest.num_source_trees
# The actual trace forest which gets traced. In most circumstances
# this will be self.trace_forest but if changing cavities enter an
# unstable regime then this will be temporarily swapped out for the
# contingent forest for this data point (see below).
TraceForest trace_forest
# Objects necessary for dealing with newly unstable cavities
tuple ch_unstable_cavities
TraceForest contingent_forest
# No changing beam parameters, do nothing
if not self.trace_forest.N_trees:
return True
# First we loop over the source trees and find any changing
# cavities which have become unstable
ch_unstable_cavities = self._find_new_unstable_cavities()
# If we did find any newly unstable cavities then the current trace_forest
# is invalidated at the current data point so we must build a new forest with
# the unstable cavities disabled
# NOTE (sjr) Don't worry about increased Python interaction in
# this block as this will rarely be executed anyway
if ch_unstable_cavities:
LOGGER.info(
"Attempting to use a contingent trace forest due "
"to the presence of unstable cavities"
)
# Look-up the combination of unstable cavities to see if a
# contingent forest was already built from this
contingent_forest = self.contingent_trace_forests.get(ch_unstable_cavities)
# No previous forest built from the given combination of disabled
# unstable cavities so need to build one here
if contingent_forest is None:
LOGGER.debug(
"For unstable cavity combination %s no cached contingent "
"trace forest found, now attempting to build a new one...",
[uc.name for uc in ch_unstable_cavities],
)
contingent_forest = self._initialise_contingent_forest(ch_unstable_cavities)
# If there are no dependencies left after disabling the
# unstable cavities then a beam trace cannot be performed
# at this data point so inform of this on return
if contingent_forest is None:
return False
# Cache the contingent forest for this combination of unstable
# cavities as these typically occur in blocks (or across strides)
# of data points so we don't want to keep rebuilding the same
# contingency forests for identical unstable cavity combos
self.contingent_trace_forests[ch_unstable_cavities] = contingent_forest
else:
LOGGER.debug(
"For unstable cavity combination %s found and using "
"cached contingent trace forest:%s",
[uc.name for uc in ch_unstable_cavities],
contingent_forest.draw_by_order()
)
for tree_idx in range(self.trace_forest.N_trees):
tree = self.trace_forest.forest[tree_idx]
# Make sure only the correctly changing beam parameters, according
# to self.trace_forest, get reflagged when exiting from the unstable
# region again
self.needs_reflag_changing_q = True
if tree.is_source:
cav_ws = self.cavity_workspaces.get(tree.dependency)
source_is_cav = cav_ws is not None
if source_is_cav: # cavity source
if not cav_ws.is_stable:
gx = cav_ws.gx
gy = cav_ws.gy
if float_eq(gx, gy):
LOGGER.warning(
"Skipping beam traces from cavity %s as it is unstable "
"with g = %s",
tree.dependency.name,
gx,
)
else:
LOGGER.warning(
"Skipping beam traces from cavity %s as it is unstable "
"with gx = %s, gy = %s",
tree.dependency.name,
gx,
gy,
)
# Use the contingent forest for this data point
trace_forest = contingent_forest
continue
# Otherwise we just use the standard changing trace forest of the simulation
else:
trace_forest = self.trace_forest
# TODO (sjr) If the above statement is triggered then I need to
# notify the dependent sub-trees that their dependency
# now has an invalid beam parameter
# -> not too sure what will happen currently but I suspect
# that the previous data point trace data will just be
# used for these external cavity tree traces; at worst
# though we may have UB occurring
# If we've just exited an unstable cavity region where a contingent trace
# forest was being used, then we need to reflag the beam parameters which
# are changing
if self.needs_reflag_changing_q:
self._determine_changing_beam_params(forest=None, set_tree_node_ids=False)
self.needs_reflag_changing_q = False
# If tree is a source and beam parameter is not None then
# we have at least one valid dependency so the beam trace
# is ok
valid = True
# Now do the actual beam tracing by simply traversing the forest
# and propagating the beam through each tree
for tree_idx in range(trace_forest.N_trees):
tree = trace_forest.forest[tree_idx]
if tree.is_source:
cav_ws = self.cavity_workspaces.get(tree.dependency)
source_is_cav = cav_ws is not None
source_node_id = tree.node_id
source_opp_node_id = tree.opp_node_id
if not source_is_cav: # gauss command
if not source_is_cav: # Source tree is from a Gauss
qx_src = tree.dependency.qx
qy_src = tree.dependency.qy
else:
else: # Source tree is from a Cavity
qx_src = cav_ws.qx
qy_src = cav_ws.qy
self.trace[source_node_id].qx = qx_src
self.trace[source_node_id].qy = qy_src
if self.trace_forest.symmetric:
if trace_forest.symmetric:
self.trace[source_opp_node_id].qx = -conj(qx_src)