General¶
-
class
ConstantFunction
(constant, **kwargs)¶ Bases:
pyinduct.core.Function
A
Function
that returns a constant value.This function can be differentiated without limits.
- Parameters
constant (number) – value to return
- Keyword Arguments
**kwargs – All other kwargs get passed to
Function
.
-
derive
(self, order=1)¶ Spatially derive this
Function
.This is done by neglecting order derivative handles and to select handle as the new evaluation_handle.
- Parameters
order (int) – the amount of derivations to perform
- Raises
TypeError – If order is not of type int.
ValueError – If the requested derivative order is higher than the provided one.
- Returns
Function
the derived function.
-
class
FieldVariable
(function_label, order=0, 0, weight_label=None, location=None, exponent=1, raised_spatially=False)¶ Bases:
pyinduct.placeholder.Placeholder
Class that represents terms of the systems field variable .
- Parameters
function_label (str) – Label of shapefunctions to use for approximation, see
register_base()
for more information about how to register an approximation basis.tuple of int (order) – Tuple of temporal_order and spatial_order derivation order.
weight_label (str) – Label of weights for which coefficients are to be calculated (defaults to function_label).
location – Where the expression is to be evaluated.
exponent – Exponent of the term.
Examples
Assuming some shapefunctions have been registered under the label
"phi"
the following expressions hold:>>> x_dt_dzz = FieldVariable("phi", order=(1, 2))
>>> x_dtt_at_3 = FieldVariable("phi", order=(2, 0), location=3)
-
derive
(self, *, temp_order=0, spat_order=0)¶ Derive the expression to the specified order.
- Parameters
temp_order – Temporal derivative order.
spat_order – Spatial derivative order.
- Returns
The derived expression.
- Return type
Note
This method uses keyword only arguments, which means that a call will fail if the arguments are passed by order.
-
class
Function
(eval_handle, domain=- np.inf, np.inf, nonzero=- np.inf, np.inf, derivative_handles=None)¶ Bases:
pyinduct.core.BaseFraction
Most common instance of a
BaseFraction
. This class handles all tasks concerning derivation and evaluation of functions. It is used broad across the toolbox and therefore incorporates some very specific attributes. For example, to ensure the accurateness of numerical handling functions may only evaluated in areas where they provide nonzero return values. Also their domain has to be taken into account. Therefore the attributes domain and nonzero are provided.To save implementation time, ready to go version like
LagrangeFirstOrder
are provided in thepyinduct.simulation
module.For the implementation of new shape functions subclass this implementation or directly provide a callable eval_handle and callable derivative_handles if spatial derivatives are required for the application.
- Parameters
eval_handle (callable) – Callable object that can be evaluated.
domain ((list of) tuples) – Domain on which the eval_handle is defined.
nonzero (tuple) – Region in which the eval_handle will return
output. Must be a subset of domain (nonzero) –
derivative_handles (list) – List of callable(s) that contain
of eval_handle (derivatives) –
-
add_neutral_element
(self)¶ Return the neutral element of addition for this object.
In other words: self + ret_val == self.
-
derivative_handles
(self)¶
-
derive
(self, order=1)¶ Spatially derive this
Function
.This is done by neglecting order derivative handles and to select handle as the new evaluation_handle.
- Parameters
order (int) – the amount of derivations to perform
- Raises
TypeError – If order is not of type int.
ValueError – If the requested derivative order is higher than the provided one.
- Returns
Function
the derived function.
-
static
from_data
(x, y, **kwargs)¶ Create a
Function
based on discrete data by interpolating.The interpolation is done by using
interp1d
from scipy, the kwargs will be passed.
-
function_handle
(self)¶
-
function_space_hint
(self)¶ Return the hint that this function is an element of the an scalar product space which is uniquely defined by the scalar product
scalar_product_hint()
.Note
If you are working on different function spaces, you have to overwrite this hint in order to provide more properties which characterize your specific function space. For example the domain of the functions.
-
get_member
(self, idx)¶ Implementation of the abstract parent method.
Since the
Function
has only one member (itself) the parameter idx is ignored and self is returned.- Parameters
idx – ignored.
- Returns
self
-
mul_neutral_element
(self)¶ Return the neutral element of multiplication for this object.
In other words: self * ret_val == self.
-
raise_to
(self, power)¶ Raises the function to the given power.
Warning
Derivatives are lost after this action is performed.
- Parameters
power (
numbers.Number
) – power to raise the function to- Returns
raised function
-
scalar_product_hint
(self)¶ Return the hint that the
_dot_product_l2()
has to calculated to gain the scalar product.
-
class
Input
(function_handle, index=0, order=0, exponent=1)¶ Bases:
pyinduct.placeholder.Placeholder
Class that works as a placeholder for an input of the system.
- Parameters
function_handle (callable) – Handle that will be called by the simulation unit.
index (int) – If the system’s input is vectorial, specify the element to be used.
order (int) – temporal derivative order of this term (See
Placeholder
).exponent (numbers.Number) – See
FieldVariable
.
Note
if order is nonzero, the callable is expected to return the temporal derivatives of the input signal by returning an array of
len(order)+1
.
-
class
IntegralTerm
(integrand, limits, scale=1.0)¶ Bases:
pyinduct.placeholder.EquationTerm
Class that represents an integral term in a weak equation.
- Parameters
integrand –
limits (tuple) –
scale –
-
class
Product
(a, b=None)¶ Bases:
object
Represents a product.
- Parameters
a –
b –
-
get_arg_by_class
(self, cls)¶ Extract element from product that is an instance of cls.
- Parameters
cls –
- Returns
- Return type
list
-
class
ScalarFunction
(function_label, order=0, location=None)¶ Bases:
pyinduct.placeholder.SpatialPlaceholder
Class that works as a placeholder for spatial functions in an equation. An example could be spatial dependent coefficients.
- Parameters
function_label (str) – label under which the function is registered
order (int) – spatial derivative order to use
location – location to evaluate at
- Warns
There seems to be a problem when this function is used in combination
with the :py:class:`.Product` class. Make sure to provide this class as
first argument to any product you define.
Todo
see warning.
-
static
from_scalar
(scalar, label, **kwargs)¶ create a
ScalarFunction
from scalar values.- Parameters
scalar (array like) – Input that is used to generate the placeholder. If a number is given, a constant function will be created, if it is callable it will be wrapped in a
Function
and registered.label (string) – Label to register the created base.
**kwargs – All kwargs that are not mentioned below will be passed to
Function
.
- Keyword Arguments
order (int) – See constructor.
location (int) – See constructor.
overwrite (bool) – See
register_base()
- Returns
Placeholder object that can be used in a weak formulation.
- Return type
-
class
ScalarTerm
(argument, scale=1.0)¶ Bases:
pyinduct.placeholder.EquationTerm
Class that represents a scalar term in a weak equation.
- Parameters
argument –
scale –
-
class
SecondOrderOperator
(a2=0, a1=0, a0=0, alpha1=0, alpha0=0, beta1=0, beta0=0, domain=- np.inf, np.inf)¶ Interface class to collect all important parameters that describe a second order ordinary differential equation.
- Parameters
a2 (Number or callable) – coefficient .
a1 (Number or callable) – coefficient .
a0 (Number or callable) – coefficient .
alpha1 (Number) – coefficient .
alpha0 (Number) – coefficient .
beta1 (Number) – coefficient .
beta0 (Number) – coefficient .
-
static
from_dict
(param_dict, domain=None)¶
-
static
from_list
(param_list, domain=None)¶
-
get_adjoint_problem
(self)¶ Return the parameters of the operator describing the the problem
where the are constant and whose boundary conditions are given by
The following mapping is used:
- Returns
Parameter set describing .
- Return type
-
class
TestFunction
(function_label, order=0, location=None, approx_label=None)¶ Bases:
pyinduct.placeholder.SpatialPlaceholder
Class that works as a placeholder for test functions in an equation.
- Parameters
function_label (str) – Label of the function test base.
order (int) – Spatial derivative order.
location (Number) – Point of evaluation / argument of the function.
approx_label (str) – Label of the approximation test base.
-
class
WeakFormulation
(terms, name, dominant_lbl=None)¶ Bases:
object
This class represents the weak formulation of a spatial problem. It can be initialized with several terms (see children of
EquationTerm
). The equation is interpreted as- Parameters
terms (list) – List of object(s) of type EquationTerm.
name (string) – Name of this weak form.
dominant_lbl (string) – Name of the variable that dominates this weak form.
-
compute_rad_robin_eigenfrequencies
(param, l, n_roots=10, show_plot=False)¶ Return the first
n_roots
eigenfrequencies (and eigenvalues )to the eigenvalue problem
- Parameters
param (array_like) –
l (numbers.Number) – Right boundary value of the domain .
n_roots (int) – Amount of eigenfrequencies to be compute.
show_plot (bool) – A plot window of the characteristic equation appears if it is
True
.
- Returns
- Return type
tuple –> two numpy.ndarrays of length
nroots
-
eliminate_advection_term
(param, domain_end)¶ This method performs a transformation
on the system, which eliminates the advection term from a reaction-advection-diffusion equation of the type:
The boundary can be given by robin
dirichlet
or mixed boundary conditions.
- Parameters
param (array_like) –
domain_end (float) – upper bound of the spatial domain
- Raises
TypeError – If is callable but no derivative handle is
defined for it. –
- Returns
Parameters
the transformed system
and the corresponding boundary conditions ( and/or set to None by dirichlet boundary condition).
- Return type
SecondOrderOperator or tuple
-
find_roots
(function, grid, n_roots=None, rtol=1e-05, atol=1e-08, cmplx=False, sort_mode='norm')¶ Searches n_roots roots of the function on the given grid and checks them for uniqueness with aid of rtol.
In Detail
scipy.optimize.root()
is used to find initial candidates for roots of . If a root satisfies the criteria given by atol and rtol it is added. If it is already in the list, a comprehension between the already present entries’ error and the current error is performed. If the newly calculated root comes with a smaller error it supersedes the present entry.- Raises
ValueError – If the demanded amount of roots can’t be found.
- Parameters
function (callable) – Function handle for math:f(boldsymbol{x}) whose roots shall be found.
grid (list) – Grid to use as starting point for root detection. The th element of this list provides sample points for the th parameter of .
n_roots (int) – Number of roots to find. If none is given, return all roots that could be found in the given area.
rtol – Tolerance to be exceeded for the difference of two roots to be unique: .
atol – Absolute tolerance to zero: .
cmplx (bool) – Set to True if the given function is complex valued.
sort_mode (str) – Specify tho order in which the extracted roots shall be sorted. Default “norm” sorts entries by their norm, while “component” will sort them in increasing order by every component.
- Returns
numpy.ndarray of roots; sorted in the order they are returned by .
-
get_in_domain_transformation_matrix
(k1, k2, mode='n_plus_1')¶ Returns the transformation matrix M.
M is one part of a transformation
where x is the field variable of an interior point controlled parabolic system and y is the field variable of an boundary controlled parabolic system. T is a (Fredholm-) integral transformation (which can be approximated with M).
- Parameters
k1 –
k2 –
mode –
Available modes
n_plus_1: M.shape =
2n: M.shape = (2n,2n),
- Returns
Transformation matrix M.
- Return type
numpy.array
-
get_parabolic_dirichlet_weak_form
(init_func_label, test_func_label, input_handle, param, spatial_domain)¶ Return the weak formulation of a parabolic 2nd order system, using an inhomogeneous dirichlet boundary at both sides.
- Parameters
init_func_label (str) – Label of shape base to use.
test_func_label (str) – Label of test base to use.
input_handle (
SimulationInput
) – Input.param (tuple) – Parameters of the spatial operator.
spatial_domain (#) – Spatial domain of the problem.
spatial_domain – Spatial domain of the
problem. (#) –
- Returns
Weak form of the system.
- Return type
-
get_parabolic_robin_weak_form
(shape_base_label, test_base_label, input_handle, param, spatial_domain, actuation_type_point=None)¶ Provide the weak formulation for the diffusion system with advection term, reaction term, robin boundary condition and robin actuation.
- Parameters
shape_base_label (str) – State space base label
test_base_label (str) – Test base label
input_handle (
SimulationInput
) – System inputparam (array-like) –
List of parameters:
(numbers.Number) ~ diffusion coefficient
(callable) ~ advection coefficient
(callable) ~ reaction coefficient
(numbers.Number) ~ constants for robin boundary conditions
spatial_domain (tuple) – Limits of the spatial domain
actuation_type_point (numbers.number) – Here you can shift the point of actuation from to a other point in the spatial domain.
- Returns
strings for the created base lables for the advection and reaction coefficient
- Return type
tuple