Discretisation#

class pybamm.Discretisation(mesh=None, spatial_methods=None, check_model=True, remove_independent_variables_from_rhs=False)[source]#

The discretisation class, with methods to process a model and replace Spatial Operators with Matrices and Variables with StateVectors

Parameters:
  • mesh (pybamm.Mesh) – contains all submeshes to be used on each domain

  • spatial_methods (dict) – a dictionary of the spatial methods to be used on each domain. The keys correspond to the model domains and the values to the spatial method.

  • check_model (bool, optional) – If True, model checks are performed after discretisation. For large systems these checks can be slow, so can be skipped by setting this option to False. When developing, testing or debugging it is recommended to leave this option as True as it may help to identify any errors. Default is True.

  • remove_independent_variables_from_rhs (bool, optional) – If True, model checks to see whether any variables from the RHS are used in any other equation. If a variable meets all of the following criteria (not used anywhere in the model, len(rhs)>1), then the variable is moved to be explicitly integrated when called by the solution object. Default is False.

check_model(model)[source]#

Perform some basic checks to make sure the discretised model makes sense.

check_tab_conditions(symbol, bcs)[source]#

Check any boundary conditions applied on “negative tab”, “positive tab” and “no tab”. For 1D current collector meshes, these conditions are converted into boundary conditions on “left” (tab at z=0) or “right” (tab at z=l_z) depending on the tab location stored in the mesh. For 2D current collector meshes, the boundary conditions can be applied on the tabs directly.

Parameters:
Returns:

The dictionary of boundary conditions, with the keys changed to “left” and “right” where necessary.

Return type:

dict

check_variables(model)[source]#

Check variables in variable list against rhs. Be lenient with size check if the variable in model.variables is broadcasted, or a concatenation (if broadcasted, variable is a multiplication with a vector of ones)

create_mass_matrix(model)[source]#

Creates mass matrix of the discretised model. Note that the model is assumed to be of the form M*y_dot = f(t,y), where M is the (possibly singular) mass matrix.

Parameters:

model (pybamm.BaseModel) – Discretised model. Must have attributes rhs, initial_conditions and boundary_conditions (all dicts of {variable: equation})

Returns:

  • pybamm.Matrix – The mass matrix

  • pybamm.Matrix – The inverse of the ode part of the mass matrix (required by solvers which only accept the ODEs in explicit form)

process_boundary_conditions(model)[source]#

Discretise model boundary_conditions, also converting keys to ids

Parameters:

model (pybamm.BaseModel) – Model to dicretise. Must have attributes rhs, initial_conditions and boundary_conditions (all dicts of {variable: equation})

Returns:

Dictionary of processed boundary conditions

Return type:

dict

process_dict(var_eqn_dict, ics=False)[source]#

Discretise a dictionary of {variable: equation}, broadcasting if necessary (can be model.rhs, model.algebraic, model.initial_conditions or model.variables).

Parameters:
  • var_eqn_dict (dict) – Equations ({variable: equation} dict) to dicretise (can be model.rhs, model.algebraic, model.initial_conditions or model.variables)

  • ics (bool, optional) – Whether the equations are initial conditions. If True, the equations are scaled by the reference value of the variable, if given

Returns:

new_var_eqn_dict – Discretised equations

Return type:

dict

process_initial_conditions(model)[source]#

Discretise model initial_conditions.

Parameters:

model (pybamm.BaseModel) – Model to dicretise. Must have attributes rhs, initial_conditions and boundary_conditions (all dicts of {variable: equation})

Returns:

Tuple of processed_initial_conditions (dict of initial conditions) and concatenated_initial_conditions (numpy array of concatenated initial conditions)

Return type:

tuple

process_model(model, inplace=True)[source]#

Discretise a model. Currently inplace, could be changed to return a new model.

Parameters:
  • model (pybamm.BaseModel) – Model to dicretise. Must have attributes rhs, initial_conditions and boundary_conditions (all dicts of {variable: equation})

  • inplace (bool, optional) – If True, discretise the model in place. Otherwise, return a new discretised model. Default is True.

Returns:

model_disc – The discretised model. Note that if inplace is True, model will have also been discretised in place so model == model_disc. If inplace is False, model != model_disc

Return type:

pybamm.BaseModel

Raises:

pybamm.ModelError – If an empty model is passed (model.rhs = {} and model.algebraic = {} and model.variables = {})

process_rhs_and_algebraic(model)[source]#

Discretise model equations - differential (‘rhs’) and algebraic.

Parameters:

model (pybamm.BaseModel) – Model to dicretise. Must have attributes rhs, initial_conditions and boundary_conditions (all dicts of {variable: equation})

Returns:

Tuple of processed_rhs (dict of processed differential equations), processed_concatenated_rhs, processed_algebraic (dict of processed algebraic equations) and processed_concatenated_algebraic

Return type:

tuple

process_symbol(symbol)[source]#

Discretise operators in model equations. If a symbol has already been discretised, the stored value is returned.

Parameters:

symbol (pybamm.expression_tree.symbol.Symbol) – Symbol to discretise

Returns:

Discretised symbol

Return type:

pybamm.expression_tree.symbol.Symbol

set_internal_boundary_conditions(model)[source]#

A method to set the internal boundary conditions for the submodel. These are required to properly calculate the gradient. Note: this method modifies the state of self.boundary_conditions.

set_variable_slices(variables)[source]#

Sets the slicing for variables.

Parameters:

variables (iterable of pybamm.Variables) – The variables for which to set slices