Started work on builtin LCC

src/CBD/lib/std.py

@@ -1071,6 +1071,8 @@ class Clock(CBD):
 		"""
 		Gets the current time of the clock.
 		"""
+		# FIXME: this will produce wrong signals for all clock inputs and
+		#        all blocks in the clock itself
 		sig = self.getBlockByName("TSum").getSignal("OUT1")
 		if curIt == 0 or len(sig) == 0:
 			return self.__start_time

src/CBD/simulator.py

@@ -499,6 +499,9 @@ class Simulator:
 		# TODO: LCC
 		self.signal("poststep")
 
+	def _lcc_compute(self):
+		self.__computeBlocks(self.__sim_data[1], self.__sim_data[0], self.__sim_data[2])
+
 	def _rewind(self):
 		self.__sim_data[2] -= 1
 		self.model._rewind()

src/CBD/state_events/__init__.py

@@ -0,0 +1,21 @@
+"""
+This package contains the State Event Location logic and
+helper classes.
+"""
+
+from enum import Enum
+
+class Direction(Enum):
+	"""
+	Specifies the direction of the level crossing to check.
+	"""
+
+	ANY = 0
+	"""Any crossing through the level requires an event."""
+
+	FROM_BELOW = 1
+	"""Only a crossing from below through the level causes an event."""
+
+	FROM_ABOVE = 2
+	"""Only a crossing from above through the level causes an event."""
+

src/CBD/state_events/locators.py

@@ -0,0 +1,237 @@
+"""
+This module contains the standard State Event locators.
+"""
+
+from CBD.state_events import Direction
+
+class StateEventLocator:
+	"""
+	Computes the exact level crossing time and locates when a state event must be scheduled.
+
+	Attributes:
+		sim (CBD.simulator.Simulator):  The simulator to which the locator belongs.
+		t_lower (float):                The lower range of the level crossing. It is certain
+										that the crossing happens at a time later than (or
+										equal to) this time.
+	"""
+	def __init__(self):
+		self.sim = None
+		self.t_lower = 0.0
+
+	def setSimulator(self, sim):
+		"""
+		Sets the simulator to the event locator.
+
+		Args:
+			sim (CBD.simulator.Simulator):  The current simulator.
+		"""
+		self.sim = sim
+
+	def detect(self, prev, curr, direction=Direction.ANY):
+		"""
+		Detects that a crossing happened between prev and curr.
+
+		Args:
+			prev (numeric):         The previous value.
+			curr (numeric):         The current value.
+			direction (Direction):  The direction of the crossing to detect.
+									Defaults to :attr:`Direction.ANY`.
+
+		Returns:
+			:code:`True` when the crossing happened, otherwise :code:`False`.
+		"""
+		if direction == Direction.FROM_BELOW:
+			return prev <= 0 <= curr
+		if direction == Direction.FROM_ABOVE:
+			return prev >= 0 >= curr
+		if direction == Direction.ANY:
+			return (prev <= 0 <= curr) or (prev >= 0 >= curr)
+
+		return False
+
+	def detect_signal(self, output_name, level=0.0, direction=Direction.ANY):
+		"""
+		Detects that an output port has a crossing through a specific level.
+
+		Args:
+			output_name (str):      The name of the output port.
+			level (numeric):        The level through which the value must go.
+									Defaults to 0.
+			direction (Direction):  The direction of the crossing to detect.
+									Defaults to :attr:`Direction.ANY`.
+
+		Returns:
+			:code:`True` when the crossing happened, otherwise :code:`False`.
+		"""
+		sig = self.sim.model.getSignal(output_name)
+		if len(sig) < 2:
+			# No crossing possible (yet)
+			return False
+		prev = sig[-2].value - level
+		curr = sig[-1].value - level
+
+		return self.detect(prev, curr, direction)
+
+	def _function(self, output_name, time, level=0.0):
+		"""
+		The internal function. Whenever an algorithm requires the computation of the
+		CBD model at another time, this function can be executed.
+
+		Note:
+			The CBD will remain at the computed time afterwards. Use
+			:meth:`CBD.simulator._rewind` to undo the actions of this
+			function.
+
+		Args:
+			output_name (str):  The output port name for which the crossing point must
+								be computed.
+			time (float):       The time at which the CBD must be computed. Must be
+								larger than the lower bound time.
+			level (float):      The level through which the crossing must be identified.
+								This mainly shifts the signal towards 0, as most algorithms
+								are basically root finders. If the algorithm incorporates
+								the level itself, keep this value at 0 for correct behaviour.
+								Defaults to 0.
+
+		Returns:
+			The signal value of the output at the given time, shifted towards 0.
+		"""
+		assert time >= self.t_lower
+		self.sim._rewind()
+		# TODO: actually update h
+		h = time - self.t_lower
+		self.sim._lcc_compute()
+		return self.sim.model.getSignal(output_name)[-1].value - level
+
+	def run(self, output_name, level=0.0, direction=Direction.ANY):
+		"""
+		Executes the locator for an output.
+
+		Args:
+			output_name (str):      The output port name for which the crossing
+									point must be computed.
+			level (float):          The level through which the crossing must be
+									identified. Defaults to 0.
+			direction (Direction):  The direction of the crossing to detect.
+									Defaults to :attr:`Direction.ANY`.
+
+		Returns:
+			The detected time at which the crossing is suspected to occur.
+		"""
+		sig = self.sim.model.getSignal(output_name)
+		p1 = sig[-2].time, sig[-2].value
+		p2 = sig[-1].time, sig[-1].value
+		self.t_lower = p1[0]
+		t_crossing = self.algorithm(p1, p2, output_name, level, direction)
+		# TODO: reset delta to old value
+		return t_crossing
+
+	# TODO: is the direction even required? Isn't it automatically maintained?
+	def algorithm(self, p1, p2, output_name, level=0.0, direction=Direction.ANY):
+		"""
+		The algorithm that identifies the locator functionality. Must be implemented
+		in sub-classes. Should only ever be called if a crossing exists.
+
+		Args:
+			p1 (tuple):         The (time, value) coordinate before the crossing.
+			p2 (tuple):         The (time, value) coordinate after the crossing.
+			output_name (str):  The output port name for which the crossing point
+								must be computed.
+			level (float):      The level through which the crossing must be
+								identified. Defaults to 0.
+			direction (Direction):  The direction of the crossing to detect.
+								Defaults to :attr:`Direction.ANY`.
+
+		Returns:
+			A suspected time of the crossing.
+		"""
+		raise NotImplementedError()
+
+
+class PreCrossingStateEventLocator(StateEventLocator):
+	"""
+	Assumes that the crossing happens at the start of the interval.
+	Can be used if a precise detection is not a requirement.
+
+	This implementation computes a rough under-estimate.
+	"""
+	def algorithm(self, p1, p2, output_name, level=0.0, direction=Direction.ANY):
+		return p1[0]
+
+
+class PostCrossingStateEventLocator(StateEventLocator):
+	"""
+	Assumes that the crossing happens at the end of the interval.
+	Can be used if a precise detection is not a requirement.
+
+	This implementation computes a rough over-estimate.
+
+	Corresponds to the :code:`if` statement in `Modelica <https://modelica.org/>`_,
+	whereas the other locators can be seen as the :code:`when` statement.
+	"""
+	def algorithm(self, p1, p2, output_name, level=0.0, direction=Direction.ANY):
+		return p2[0]
+
+
+class LinearStateEventLocator(StateEventLocator):
+	"""
+	Uses simple linear interpolation to compute the time of the crossing.
+	This is usually a rough, yet centered estimate.
+
+	This locator should only be used if it is known that the signal is
+	(mostly) linear between the lower and upper bounds.
+	"""
+	def algorithm(self, p1, p2, output_name, level=0.0, direction=Direction.ANY):
+		t1, y1 = p1
+		t2, y2 = p2
+		if y1 == y2:
+			return t1
+		# Use the equation of a line between two points
+		# Formula is easier if x and y axes are swapped.
+		return (t2 - t1) / (y2 - y1) * (level - y1) + t1
+
+
+class RegulaFalsiStateEventLocator(StateEventLocator):
+	"""
+	Implements the Illinois algorithm for finding the root for a crossing problem.
+
+	Args:
+		eps (float):    Half of the upper bound for the relative error.
+						Defaults to 1e-5.
+		n (int):        The maximal amount of iterations to compute. Defaults to
+						5 million iterations.
+
+	See Also:
+		https://en.wikipedia.org/wiki/Regula_falsi
+	"""
+	def __init__(self, eps=1e-5, n=5_000_000):
+		super(RegulaFalsiStateEventLocator, self).__init__()
+
+		self.eps = eps
+		self.n = n
+
+	def algorithm(self, p1, p2, output_name, level=0.0, direction=Direction.ANY):
+		# direction unused, because the algorithm will automatically maintain
+		#   the crossing direction
+		t1, y1 = p1
+		t2, y2 = p2
+		tn, yn = t1, y1
+
+		side = 0
+		for i in range(self.n):
+			if abs(t1 - t2) < self.eps * abs(t1 + t2): break
+			tn = (y1 * t2 - y2 * t1) / (y1 - y2)
+			yn = self._function(output_name, tn, level)
+			if yn * y2 > 0:
+				t2, y2 = tn, yn
+				if side == -1:
+					y1 /= 2
+				side = -1
+			elif yn * y1 > 0:
+				t1, y1 = tn, yn
+				if side == 1:
+					y2 /= 2
+				side = 1
+			else:
+				break
+		return tn