Coverage for physioblocks / computing / quantities.py: 100%
73 statements
« prev ^ index » next coverage.py v7.13.1, created at 2026-01-09 16:40 +0100
1# SPDX-FileCopyrightText: Copyright INRIA
2#
3# SPDX-License-Identifier: LGPL-3.0-only
4#
5# Copyright INRIA
6#
7# This file is part of PhysioBlocks, a library mostly developed by the
8# [Ananke project-team](https://team.inria.fr/ananke) at INRIA.
9#
10# Authors:
11# - Colin Drieu
12# - Dominique Chapelle
13# - François Kimmig
14# - Philippe Moireau
15#
16# PhysioBlocks is free software: you can redistribute it and/or modify it under the
17# terms of the GNU Lesser General Public License as published by the Free Software
18# Foundation, version 3 of the License.
19#
20# PhysioBlocks is distributed in the hope that it will be useful, but WITHOUT ANY
21# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
22# PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.
23#
24# You should have received a copy of the GNU Lesser General Public License along with
25# PhysioBlocks. If not, see <https://www.gnu.org/licenses/>.
27"""
28Define **Quantities** to represent discretized numeric values along with some
29convenience functions to handle them.
30"""
32from collections.abc import Iterable
33from typing import Any, Generic, TypeVar
35import numpy as np
36from numpy.typing import NDArray
38from physioblocks.base.operators import AbstractBaseOperators
40T = TypeVar("T", np.float64, NDArray[np.float64])
43class Quantity(AbstractBaseOperators, Generic[T]):
44 """
45 Represent discretized values by holding their current and new numeric values.
47 Examples
48 ^^^^^^^^
50 >>> scalar = Quantity(0.5)
51 >>> scalar.current == scalar.new # 0.5
53 >>> vector = Quantity([0.1, 0.2, 0.3])
54 >>> vector.current == vector.new # [0.1, 0.2, 0.3]
56 It can also be initialized with numpy arrays
58 >>> Quantity(np.ones(3)) # [1, 1, 1]
60 Matrixes are flattened when initializing a quantity: it can only
61 be a scalar or vector.
63 >>> Quantity(np.ones((2, 2))) # [1, 1, 1, 1]
64 """
66 _current: T
67 """Current value (at step n) of the quantity"""
69 _new: T
70 """New value (at step n + 1) of the quantity"""
72 def __init__(self, value: Any):
73 """
74 Quantity Constructor
76 :param value: Initialization value for the Quantity
77 :type value: Any
78 """
79 if isinstance(value, Quantity):
80 self.initialize(value.current)
81 self.update(value.new)
82 else:
83 self.initialize(value)
85 @property
86 def current(self) -> T:
87 """
88 Get the current value (at step n) of the quantity
90 :return: The current value
91 :rtype: np.float64 | NDArray[np.float64]
92 """
93 return self._current
95 @property
96 def new(self) -> T:
97 """
98 Get the new value (at step n + 1) of the quantity
100 :return: The new value
101 :rtype: np.float64 | NDArray[np.float64]
102 """
103 return self._new
105 @property
106 def size(self) -> int:
107 """
108 Get the vector size of the quantity value (or 1 for scalars)
110 :return: The vector size of the quantity
111 :rtype: int
112 """
113 return self.current.size
115 def update(self, new: Any) -> None:
116 """
117 Set the new value (at step n + 1) of the quantity
119 :param new: The new value
120 :type new: Any
122 :raise ValueError: Exception raised if the new value and the quantity
123 value are not the same size
125 Example
126 ^^^^^^^
128 .. code:: python
130 >>> q = Quantity(0.5)
131 >>> q.current # 0.5
132 >>> q.new # 0.5
134 >>> q.update(0.6)
135 >>> q.current # 0.5
136 >>> q.new # 0.6
137 """
138 new_value = np.array(new) if isinstance(new, Iterable) else np.float64(new)
140 if new_value.size != self.size:
141 raise ValueError("New value should be the same size as the current value")
143 self._new = new_value # type: ignore
145 def initialize(self, value: Any) -> None:
146 """
147 Set the current and new value of the quantity
149 :param value: The value to set
150 :type value: Any
152 Example
153 ^^^^^^^
155 .. code:: python
157 >>> q = Quantity(0.5)
158 >>> q.update(0.6)
159 >>> q.current # 0.5
160 >>> q.new # 0.6
162 >>> q.initialize(0.1)
163 >>> q.current # 0.1
164 >>> q.new # 0.1
165 """
166 if isinstance(value, Iterable):
167 array_value = np.array(value, np.float64).flatten()
168 if array_value.size != 0:
169 self._current = array_value # type: ignore
170 else:
171 raise ValueError(
172 str.format("Quantity size can not be 0. Got: {0}", array_value)
173 )
174 else:
175 self._current = np.float64(value) # type: ignore
177 self._new = self._current
179 # Base operators
181 @property
182 def operation_value(self) -> T:
183 """
184 Value used with base operators.
186 Direct operation are allowed and performed on the ``current`` value of the
187 quantity.
189 :return: the current value
190 :rtype: Any
193 Example
194 ^^^^^^^
196 .. code:: python
198 >>> q1 = Quantity(1.0)
199 >>> q2 = Quantity(2.0)
200 >>> q1.update(1.5)
202 # Quantity sum is performed on current value:
203 >>> q1 + q2 # 3.0
205 >>> q1.initialise(1.5)
206 >>> q1 + q2 # 3.5
208 In place operators call ``initialise`` on the **Quantity**
210 Example
211 ^^^^^^^
213 .. code:: python
215 >>> q1 = Quantity(1.0)
216 >>> q1 += 0.5
217 >>> q1.current # 1.5
218 >>> q1.new # 1.5
220 """
221 return self.current
223 # In place Operations
225 def __iadd__(self, other: Any) -> "Quantity[Any]":
226 self.initialize(self.current + other)
227 return self
229 def __isub__(self, other: Any) -> "Quantity[Any]":
230 self.initialize(self.current - other)
231 return self
233 def __imul__(self, other: Any) -> "Quantity[Any]":
234 self.initialize(self.current * other)
235 return self
237 def __imatmul__(self, other: Any) -> "Quantity[Any]":
238 self.initialize(self.current @ other)
239 return self
241 def __itruediv__(self, other: Any) -> "Quantity[Any]":
242 self.initialize(self.current / other)
243 return self
245 def __ifloordiv__(self, other: Any) -> "Quantity[Any]":
246 self.initialize(self.current // other)
247 return self
249 def __imod__(self, other: Any) -> "Quantity[Any]":
250 self.initialize(self.current % other)
251 return self
253 def __ipow__(self, other: Any) -> "Quantity[Any]":
254 self.initialize(self.current**other)
255 return self
258def diff(q: Quantity[T]) -> T:
259 """
260 Compute the difference between the new and the current value of the given
261 **Quantity**
263 :param q: The quantity
264 :type q: Quantity
266 :return: The difference betwen new and current for the quantity
267 :rtype: np.float64 | NDArray[np.float64]
269 Examples
270 ^^^^^^^^
272 .. code:: python
274 >>> q = Quantity(1.0)
275 >>> diff(q) # 0.0
277 >>> q.update(1.1)
278 >>> diff(q) # 0.1
280 """
281 return q.new - q.current
284def mid_alpha(q: Quantity[Any], scheme_time_shift: Any) -> Any:
285 """
286 Compute the discretized value of the Quantity for a shifted time scheme.
288 .. math::
290 x^{n+\\alpha} = (\\frac{1}{2} - \\alpha) x^{n} + (\\frac{1}{2} +
291 \\alpha) x^{n + 1}
293 :param q: The quantity
294 :type value: Quantity
296 :param scheme_time_shift: the time shift of the scheme.
297 :type scheme_time_shift: float
299 :raise ValueError: If scheme_time_shift is not a scalar,
300 exception raised if scheme_time_shift and q are not the same size.
302 :return: the discretized value of the quantity
303 :rtype: np.float64 | NDArray[np.float64]
305 Example
306 ^^^^^^^
308 .. code:: python
310 >>> q = Quantity(1.0)
311 >>> q.update(2.0)
312 >>> mid_alpha(q, 0.25) # 1.75
313 """
314 return (0.5 - scheme_time_shift) * q.current + (0.5 + scheme_time_shift) * q.new
317def mid_point(q: Quantity[Any]) -> Any:
318 """
319 Compute the discretized value of the Quantity for a mid point time scheme
321 .. math::
323 x^{n + \\frac{1}{2}} = \\frac{1}{2} (x^{n} + x^{n + 1})
325 :param q: Quantity
326 :type q: Quantity
328 :return: Discretized value of the quantity
329 :rtype: np.float64 | NDArray[np.float64]
331 Example
332 ^^^^^^^
334 .. code:: python
336 >>> q = Quantity(1.0)
337 >>> q.update(2.0)
338 >>> mid_point(q) # 1.5
339 """
340 return 0.5 * (q.current + q.new)
343def sign(q: Quantity[np.float64]) -> float:
344 """
345 Get 1 if the new value is superior to current value,
346 -1 othwerwise
348 :param q: the quantity
349 :type q: Quantity
351 :return: 1 or -1
352 :rtype: np.float64
353 """
354 return -1.0 if diff(q) < 0.0 else 1.0