triflow.core package

Submodules

triflow.core.fields module

class triflow.core.fields.BaseFields(**inputs)[source]

Bases: object

Specialized container which expose the data as a structured numpy array, give access to the dependants variables and the herlpers function as attributes (as a numpy rec array) and is able to give access to a flat view of the dependent variables only (which is needed by the ode solvers for all the linear algebra manipulation).

Parameters:**inputs (numpy.array) – named argument providing x, the dependent variables and the helper functions. All of these are mendatory and a KeyError will be raised if a data is missing.
array

numpy.array – vanilla numpy array containing the data

size

int – Number of discretisation nodes

copy()[source]
static factory(dependent_variables, helper_functions)[source]

Fields factory generating specialized container build around a triflow Model.

Parameters:
  • dependent_variables (iterable of str) – name of the dependent variables
  • helper_functions (iterable of str) – name of the helper functions
Returns:

Specialized container which expose the data as a structured numpy array

Return type:

triflow.BaseFields

fill(flat_array)[source]

take a flat numpy array and update inplace the dependent variables of the container

Parameters:flat_array (numpy.ndarray) – flat array which will be put in the dependant variable flat array.
flat

numpy.ndarray.view – flat view of the main numpy array

structured

numpy.ndarray.view – structured view of the main numpy array

uarray

numpy.ndarray.view – view of the dependent variables of the main numpy array

uflat

return a flatten copy of the main numpy array with only the dependant variables.

Be carefull, modification of these data will not be reflected on the main array!

triflow.core.model module

class triflow.core.model.Model(differential_equations, dependent_variables, parameters=None, help_functions=None)[source]

Bases: object

Contain finite difference approximation and routine of the dynamical system

Take a mathematical form as input, use Sympy to transform it as a symbolic expression, perform the finite difference approximation and expose theano optimized routine for both the right hand side of the dynamical system and Jacobian matrix approximation.

Parameters:
  • differential_equations (iterable of str or str) – the right hand sides of the partial differential equations written as \(\frac{\partial U}{\partial t} = F(U)\), where the spatial derivative can be written as dxxU or dx(U, 2) or with the sympy notation Derivative(U, x, x)
  • dependent_variables (iterable of str or str) – the dependent variables with the same order as the temporal derivative of the previous arg.
  • parameters (iterable of str or str, optional, default None) – list of the parameters. Can be feed with a scalar of an array with the same size
  • help_functions (None, optional) – All fields which have not to be solved with the time derivative but will be derived in space.
F

triflow.F_Routine – Callable used to compute the right hand side of the dynamical system

F_array

numpy.ndarray of sympy.Expr – Symbolic expressions of the right hand side of the dynamical system

J

triflow.J_Routine – Callable used to compute the Jacobian of the dynamical system

J_array

numpy.ndarray of sympy.Expr – Symbolic expressions of the Jacobian side of the dynamical system

Properties
----------
fields_template

Model specific Fields container used to store and access to the model variables in an efficient way.

save: Save a binary of the Model with pre-optimized F and J routines

Examples

A simple diffusion equation:

>>> from triflow import Model
>>> model = Model("k * dxxU", "U", "k")

A coupled system of convection-diffusion equation:

>>> from triflow import Model
>>> model = Model(["k1 * dxxU - c1 * dxV",
...                "k2 * dxxV - c2 * dxU",],
...                ["U", "V"], ["k1", "k2", "c1", "c2"])
fields_template
static load(filename)[source]

load a pre-compiled triflow model. The internal of theano allow a caching of the model. Will be slow if it is the first time the model is loaded on the system.

Parameters:filename (str) – path of the pre-compiled model
Returns:triflow pre-compiled model
Return type:triflow.core.Model
save(filename)[source]

Save the model as a binary pickle file.

Parameters:filename (str) – name of the file where the model is saved.
Returns:
Return type:None

triflow.core.routines module

class triflow.core.routines.F_Routine(matrix, args, pars, ufunc, reduced=False)[source]

Bases: triflow.core.routines.ModelRoutine

Compute the right hand side of the dynamical system \(\frac{\partial U}{\partial t} = F(U)\)

Parameters:
  • fields (triflow.Fields) – triflow fields container generated by a triflow.Model containing the actual state of the dependent variables and helper functions.
  • pars (dict) – dictionnary with the different physical parameters of the model and the ‘periodic’ key.
Returns:

flat array containing the right hand side of the dynamical system.

Return type:

numpy.ndarray

diff_approx(fields, pars, eps=1e-08)[source]
class triflow.core.routines.J_Routine(matrix, args, pars, ufunc, reduced=False)[source]

Bases: triflow.core.routines.ModelRoutine

Compute the right hand side of the dynamical system \(\frac{\partial U}{\partial t} = F(U)\)

Parameters:
  • fields (triflow.Fields) – triflow fields container generated by a triflow.Model containing the actual state of the dependent variables and helper functions.
  • pars (dict) – dictionnary with the different physical parameters of the model and the ‘periodic’ key.
  • sparse (bool, optional, default True) – whether should the matrix returned as dense or sparse form.
Returns:

scipy.sparse.CSC or numpy.ndarray

Return type:

sparse or dense form (depending of the sparse argument) of the Jacobian approximation of the dynamical system right hand side.

class triflow.core.routines.ModelRoutine(matrix, args, pars, ufunc, reduced=False)[source]

Bases: object

triflow.core.simulation module

class triflow.core.simulation.Simulation(model, t, fields, physical_parameters, dt, id=None, hook=<function Simulation.<lambda>>, scheme=<class 'triflow.plugins.schemes.RODASPR'>, tmax=None, **kwargs)[source]

Bases: object

High level container used to run simulation build on triflow Model. This object is an iterable which will yield every time step until the parameters ‘tmax’ is reached if provided. By default, the solver use a 6th order ROW solver, an implicit method with integrated time-stepping.

Parameters:
  • model (triflow.Model) – Contain finite difference approximation and routine of the dynamical system
  • t (float) – initial time
  • fields (triflow.Fields) – triflow container filled with initial conditions
  • physical_parameters (dict) – physical parameters of the simulation
  • id (None, optional) – name of the simulation. A 2 word slug will be generated if not provided.
  • 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...
  • scheme (callable, optional, default triflow.schemes.RODASPR) – an callable object which take the simulation state and return the next step. Its signature is scheme.__call__(fields, t, dt, pars, hook) and it should return the next time and the updated fields. It take the model and extra positional and named arguments.
  • **kwargs (*args,) –

    extra arguments passed to the scheme.

  • **kwargs

    extra arguments passed to the scheme.

dt

float – output time step

fields

triflow.Fields – triflow container filled with actual data

i

int – actual iteration

id

str – name of the simulation

model

triflow.Model – triflow Model used in the simulation

physical_parameters

dict – physical parameters of the simulation

status

str – status of the simulation, one of the following one: (‘created’, ‘running’, ‘finished’, ‘failed’)

t

float – actual time

tmax

float or None, default None – stopping time of the simulation. Not stopping if set to None.

Examples

>>> import numpy as np
>>> import triflow
>>> model = triflow.Model(["k1 * dxxU",
...                        "k2 * dxxV"],
...                       ["U", "V"],
...                       ["k1", "k2"])
>>> x = np.linspace(0, 100, 1000, endpoint=False)
>>> U = np.cos(x * 2 * np.pi / 100)
>>> V = np.sin(x * 2 * np.pi / 100)
>>> fields = model.fields_template(x=x, U=U, V=V)
>>> pars = {'k1': 1, 'k2': 1, 'periodic': True}
>>> simulation = triflow.Simulation(model, 0, fields,
...                                 pars, dt=5, tmax=50)
>>> for t, fields in simulation:
...    pass
>>> print(t)
50
add_display(display, *display_args, **display_kwargs)[source]

add a display for the simulation.

Parameters:
  • display (callable) – a display as the one available in triflow.displays
  • *display_args – positional arguments for the display function (other than the simulation itself)
  • **display_kwargs – named arguments for the display function
compute()[source]

Generator which yield the actual state of the system every dt.

Yields:tuple (t, fields) – Actual time and updated fields container.

Module contents