triflow.plugins package¶
Submodules¶
triflow.plugins.displays module¶
This module regroups different displays: function and coroutine written in order to give extra informationto the user during the simulation (plot, post-processing...)
-
class
triflow.plugins.displays.
bokeh_fields_update
(simul, keys=None, line_kwargs={}, fig_kwargs={}, notebook=True)[source]¶ Bases:
object
Display fields data in a interactive Bokeh plot displayed in a jupyter notebook.
Parameters: - keys (None, optional) – list of the dependant variables to be displayed
- line_kwargs (dict of dict) – dictionnary with vars as key and a dictionnary of keywords arguments passed to the lines plots
- fig_kwargs (dict of dict) – dictionnary with vars as key and a dictionnary of keywords arguments passed to the figs plots
- init_notebook (True, optional) – if True, initialize the javascript component needed for bokeh.
-
class
triflow.plugins.displays.
bokeh_probes_update
(simul, probes, line_kwargs={}, fig_kwargs={}, notebook=True)[source]¶ Bases:
object
Display custom probes in a interactive Bokeh plot displayed in a jupyter notebook.
Parameters: - probes (dictionnary of callable) – Dictionnary with {name: callable} used to plot the probes. The signature is the same as in the hooks and return the value we want to plot.
- line_kwargs (dict of dict) – dictionnary with vars as key and a dictionnary of keywords arguments passed to the lines plots
- fig_kwargs (dict of dict) – dictionnary with vars as key and a dictionnary of keywords arguments passed to the figs plots
- init_notebook (True, optional) – if True, initialize the javascript component needed for bokeh.
triflow.plugins.schemes module¶
This module regroups all the implemented temporal schemes. They are written as callable class which take the model and some control arguments at the init, and perform a computation step every time they are called.
- The following solvers are implemented:
- Backward and Forward Euler, Crank-Nicolson method (with the Theta class)
- Some Rosenbrock Wanner schemes (up to the 6th order) with time controler
- All the scipy.integrate.ode integrators with the scipy_ode class.
-
class
triflow.plugins.schemes.
RODASPR
(model, tol=0.01, time_stepping=True, max_iter=None, dt_min=None)[source]¶ Bases:
triflow.plugins.schemes.ROW_general
6th order Rosenbrock scheme, with time stepping
Parameters: - model (triflow.Model) – triflow Model
- tol (float, optional, default 1E-2) – tolerance factor for the time stepping. The time step will adapt to ensure that the maximum relative error on all fields stay under that value.
- time_stepping (bool, optional, default True) – allow a variable internal time-step to ensure good agreement between computing performance and accuracy.
- max_iter (float or None, optional, default None) – maximum internal iteration allowed
- dt_min (float or None, optional, default None) – minimum internal time step allowed
-
class
triflow.plugins.schemes.
ROS2
(model)[source]¶ Bases:
triflow.plugins.schemes.ROW_general
Second order Rosenbrock scheme, without time stepping
Parameters: model (triflow.Model) – triflow Model
-
class
triflow.plugins.schemes.
ROS3PRL
(model, tol=0.01, time_stepping=True, max_iter=None, dt_min=None)[source]¶ Bases:
triflow.plugins.schemes.ROW_general
4th order Rosenbrock scheme, with time stepping
Parameters: - model (triflow.Model) – triflow Model
- tol (float, optional, default 1E-2) – tolerance factor for the time stepping. The time step will adapt to ensure that the maximum relative error on all fields stay under that value.
- time_stepping (bool, optional, default True) – allow a variable internal time-step to ensure good agreement between computing performance and accuracy.
- max_iter (float or None, optional, default None) – maximum internal iteration allowed
- dt_min (float or None, optional, default None) – minimum internal time step allowed
-
class
triflow.plugins.schemes.
ROS3PRw
(model, tol=0.01, time_stepping=True, max_iter=None, dt_min=None)[source]¶ Bases:
triflow.plugins.schemes.ROW_general
Third order Rosenbrock scheme, with time stepping
Parameters: - model (triflow.Model) – triflow Model
- tol (float, optional, default 1E-2) – tolerance factor for the time stepping. The time step will adapt to ensure that the maximum relative error on all fields stay under that value.
- time_stepping (bool, optional, default True) – allow a variable internal time-step to ensure good agreement between computing performance and accuracy.
- max_iter (float or None, optional, default None) – maximum internal iteration allowed
- dt_min (float or None, optional, default None) – minimum internal time step allowed
-
class
triflow.plugins.schemes.
ROW_general
(model, alpha, gamma, b, b_pred=None, time_stepping=False, tol=None, max_iter=None, dt_min=None)[source]¶ Bases:
object
Rosenbrock Wanner class of temporal solvers
The implementation and the different parameters can be found in http://www.digibib.tu-bs.de/?docid=00055262
-
__call__
(t, fields, dt, pars, hook=<function ROW_general.<lambda>>)[source]¶ Perform a step of the solver: took a time and a system state as a triflow Fields container and return the next time step with updated container.
Parameters: - t (float) – actual time step
- fields (triflow.Fields) – actual system state in a triflow Fields
- dt (float) – temporal step-size
- pars (dict) – physical parameters of the model
- hook (callable, optional) – any callable taking the actual time, fields and parameters and return modified fields and parameters. Will be called every internal time step and can be used to include time dependent or conditionnal parameters, boundary conditions...
- container –
Returns: tuple – updated time and fields container
Return type: Raises: NotImplementedError
– raised if a time stepping is requested but the scheme do not provide the b predictor coefficients.ValueError
– raised if time_stepping is True and tol is not provided.
-
-
class
triflow.plugins.schemes.
Theta
(model, theta=1, solver=<function spsolve>)[source]¶ Bases:
object
- Simple theta-based scheme where theta is a weight
- if theta = 0, the scheme is a forward-euler scheme if theta = 1, the scheme is a backward-euler scheme if theta = 0.5, the scheme is called a Crank-Nicolson scheme
Parameters: - model (triflow.Model) – triflow Model
- theta (int, optional, default 1) – weight of the theta-scheme
- solver (callable, optional, default scipy.sparse.linalg.spsolve) – method able to solve a Ax = b linear equation with A a sparse matrix. Take A and b as argument and return x.
-
__call__
(t, fields, dt, pars, hook=<function Theta.<lambda>>)[source]¶ Perform a step of the solver: took a time and a system state as a triflow Fields container and return the next time step with updated container.
Parameters: - t (float) – actual time step
- fields (triflow.Fields) – actual system state in a triflow Fields container
- dt (float) – temporal step-size
- pars (dict) – physical parameters of the model
- hook (callable, optional) – any callable taking the actual time, fields and parameters and return modified fields and parameters. Will be called every internal time step and can be used to include time dependent or conditionnal parameters, boundary conditions...
Returns: tuple – updated time and fields container
Return type:
-
class
triflow.plugins.schemes.
scipy_ode
(model, integrator='vode', **integrator_kwargs)[source]¶ Bases:
object
Proxy written around the scipy.integrate.ode class. Give access to all the scpy integrators.
Parameters: - model (triflow.Model) – triflow Model
- integrator (str, optional, default 'vode') – name of the chosen scipy integration scheme.
- **integrator_kwargs – extra arguments provided to the scipy integration scheme.
-
__call__
(t, fields, dt, pars, hook=<function scipy_ode.<lambda>>)[source]¶ Perform a step of the solver: took a time and a system state as a triflow Fields container and return the next time step with updated container.
Parameters: - t (float) – actual time step
- fields (triflow.Fields) – actual system state in a triflow Fields
- dt (float) – temporal step-size
- pars (dict) – physical parameters of the model
- hook (callable, optional) – any callable taking the actual time, fields and parameters and return modified fields and parameters. Will be called every internal time step and can be used to include time dependent or conditionnal parameters, boundary conditions...
- container –
Returns: tuple – updated time and fields container
Return type: Raises: RuntimeError
– Description
triflow.plugins.signals module¶
This module provides a Signal class usefull for time variable boundary conditions.
- Available signals:
- ConstantSignal: just an offset signal
- ForcedSignal: a sinusoidal wave
- WhiteNoise: a noisy signal
- BrownNoise: a noisy signal with Fourrier modes set to 0 for a fraction of the available modes.
-
class
triflow.plugins.signals.
AdditiveSignal
(signal_a, signal_b, op)[source]¶ Bases:
triflow.plugins.signals.Signal
Additive signal. Proxy the different signals and sum their interpolation functions.
-
class
triflow.plugins.signals.
ConstantSignal
(offset: float)[source]¶ Bases:
triflow.plugins.signals.Signal
Offset signal
Parameters: offset (float) – Value of the offset Examples
Set an offset to a white noise:
>>> import numpy as np >>> from triflow.plugins import signals >>> t = np.linspace(0, 1, 1000) >>> noisy = signals.GaussianWhiteNoise(frequency_sampling=200) >>> offset = signals.ConstantSignal(offset=1) >>> noise_with_offset = noisy + offset >>> noise_without_offset = noise_with_offset - offset >>> np.isclose(noisy(t), noise_without_offset(t)).all() True
-
class
triflow.plugins.signals.
GaussianBrownNoise
(frequency_sampling: float, frequency_cut: float, mean=0, std=0.2, n=1000, filter_std=25, filter_window_size=500, seed=None)[source]¶ Bases:
triflow.plugins.signals.GaussianWhiteNoise
Gaussian Brown noise signal, seeded with numpy.random.randn and with frequency cut at specific value.
Parameters: - frequency_sampling (float) – Frequency sampling (the highest frequency of the white signal spectrum).
- frequency_cut (float) – Filter frequency cut.
- mean (int, optional, default 1000) – sampling number of the signal.
- n (int, optional, default 1000) – sampling number of the signal.
- seed (int or None, optional) – pseudo-random number generator seed. Signals with same signal_period, sampling number and seed will be similar.
Examples
White signal with frequencies cut after 50 Hz:
>>> from triflow.plugins import signals >>> brown_signal = signals.GaussianBrownNoise(frequency_sampling=200, ... frequency_cut=50) >>> spectrum_frequencies, spectrum_density = brown_signal.fourrier_spectrum >>> np.isclose(spectrum_density[spectrum_frequencies > 50], 0, atol=1E-3).all() True
-
class
triflow.plugins.signals.
GaussianWhiteNoise
(frequency_sampling, mean=0, std=0.2, n=1000, seed=None)[source]¶ Bases:
triflow.plugins.signals.Signal
Gaussian White noise signal, seeded with numpy.random.randn
Parameters: - frequency_sampling (float) – Frequency sampling (the highest frequency of the white signal spectrum).
- mean (int, optional, default 1000) – sampling number of the signal.
- n (int, optional, default 1000) – sampling number of the signal.
- seed (int or None, optional) – pseudo-random number generator seed. Signals with same signal_period, sampling number and seed will be similar.
Examples
Simple white signal with a 200 Hz frequency sampling:
>>> import numpy as np >>> from triflow.plugins import signals >>> signal = signals.GaussianWhiteNoise(frequency_sampling=200)
Forcing the pseudo-random-number generator seed for reproductible signal:
>>> import numpy as np >>> from triflow.plugins import signals >>> t = np.linspace(0, 10, 1000) >>> signal1 = signals.GaussianWhiteNoise(frequency_sampling=200, seed=50) >>> signal2 = signals.GaussianWhiteNoise(frequency_sampling=200, seed=50) >>> np.isclose(signal1(t) - signal2(t), 0).all() True >>> signal3 = signals.GaussianWhiteNoise(frequency_sampling=200) >>> np.isclose(signal1(t) - signal3(t), 0).all() False
-
class
triflow.plugins.signals.
Signal
(signal_period: float, n: int = 1000, **kwargs)[source]¶ Bases:
object
Base class for signal object.
Parameters: - signal_period (float) – period of the signal.
- n (int, optional) – sampling number of the signal.
- **kwargs – extra arguments provided to the custom signal template.
Raises: NotImplementedError
- this class is not supposed to be used by the user. If initialized directly, it will raise a NotImplementedError.
-
__add__
(other_signal)[source]¶ Parameters: other_signal (triflow.plugins.signal.Signal) – second signal to which this one will be added. Returns: added signal. Return type: triflow.plugins.signal.AdditiveSignal
-
__call__
(t: float)[source]¶ return the signal value for the time t.
Parameters: t (float) – time where the signal is evaluated. Returns: amplitude of the signal a the time t. Return type: float
-
class
triflow.plugins.signals.
SinusoidalSignal
(frequency: float, amplitude: float, phase: float = 0, n: int = 1000)[source]¶ Bases:
triflow.plugins.signals.Signal
Simple sinusoidal signal
Parameters: - frequency (float) – frequency of the signal
- amplitude (float) – amplitude of the signal
- phase (float, optional) – phase of the signam
- n (int, optional) – number of sample for 2 period of the signal
Examples
We generate a 5 Hz signal with an amplitude of 0.5, and we check the signal.
>>> from triflow.plugins import signals >>> t = np.linspace(0, 100, 10000) >>> signal = signals.SinusoidalSignal(frequency=5, amplitude=.5) >>> spectrum_frequencies, spectrum_density = signal.fourrier_spectrum >>> print(f"Amplitude: {signal(t).max():g}") Amplitude: 0.499999 >>> print(f"main mode: {spectrum_frequencies[spectrum_density.argmax()]:g} Hz") main mode: 4.995 Hz
The class give a warning if the number of sample is too low and lead to aliasing. You can try it with
>>> import logging >>> logger = logging.getLogger('triflow.plugins.signals') >>> logger.handlers = [] >>> logger.addHandler(logging.StreamHandler()) >>> logger.setLevel('INFO') >>> aliased_signal = signals.SinusoidalSignal(5, 1, n=30)