Scikit Finite Elements#
- class pybamm.ScikitFiniteElement(options=None)[source]#
A class which implements the steps specific to the finite element method during discretisation. The class uses scikit-fem to discretise the problem to obtain the mass and stiffness matrices. At present, this class is only used for solving the Poisson problem -grad^2 u = f in the y-z plane (i.e. not the through-cell direction).
For broadcast, we follow the default behaviour from SpatialMethod.
Extends:
pybamm.spatial_methods.spatial_method.SpatialMethod
- assemble_mass_form(symbol, boundary_conditions, region='interior')[source]#
Assembles the form of the finite element mass matrix over the domain interior or boundary.
- Parameters:
symbol (
pybamm.Variable
) – The variable corresponding to the equation for which we are calculating the mass matrix.boundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
region (str, optional) – The domain over which to assemble the mass matrix form. Can be “interior” (default) or “boundary”.
- Returns:
The (sparse) mass matrix for the spatial method.
- Return type:
- bc_apply(M, boundary, zero=False)[source]#
Adjusts the assembled finite element matrices to account for boundary conditions.
- Parameters:
M (
scipy.sparse.coo_matrix
) – The assembled finite element matrix to adjust.boundary (
numpy.array
) – Array of the indices which correspond to the boundary.zero (bool, optional) – If True, the rows of M given by the indicies in boundary are set to zero. If False, the diagonal element is set to one. default is False.
- boundary_integral(child, discretised_child, region)[source]#
Implementation of the boundary integral operator. See
pybamm.SpatialMethod.boundary_integral()
- boundary_integral_vector(domain, region)[source]#
A node in the expression tree representing an integral operator over the boundary of a domain
\[I = \int_{\partial a}\!f(u)\,du,\]where \(\partial a\) is the boundary of the domain, and \(u\in\text{domain boundary}\).
- Parameters:
domain (list) – The domain(s) of the variable in the integrand
region (str) – The region of the boundary over which to integrate. If region is entire the integration is carried out over the entire boundary. If region is negative tab or positive tab then the integration is only carried out over the appropriate part of the boundary corresponding to the tab.
- Returns:
The finite element integral vector for the domain
- Return type:
- boundary_mass_matrix(symbol, boundary_conditions)[source]#
Calculates the mass matrix for the finite element method assembled over the boundary.
- Parameters:
symbol (
pybamm.Variable
) – The variable corresponding to the equation for which we are calculating the mass matrix.boundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
- Returns:
The (sparse) mass matrix for the spatial method.
- Return type:
- boundary_value_or_flux(symbol, discretised_child, bcs=None)[source]#
Returns the average value of the symbol over the negative tab (“negative tab”) or the positive tab (“positive tab”) in the Finite Element Method.
Overwrites the default
pybamm.SpatialMethod.boundary_value()
- definite_integral_matrix(child, vector_type='row')[source]#
Matrix for finite-element implementation of the definite integral over the entire domain
\[I = \int_{\Omega}\!f(s)\,dx\]for where \(\Omega\) is the domain.
- Parameters:
child (
pybamm.Symbol
) – The symbol being integratedvector_type (str, optional) – Whether to return a row or column vector (default is row)
- Returns:
The finite element integral vector for the domain
- Return type:
- divergence(symbol, discretised_symbol, boundary_conditions)[source]#
Matrix-vector multiplication to implement the divergence operator. See
pybamm.SpatialMethod.divergence()
- gradient(symbol, discretised_symbol, boundary_conditions)[source]#
Matrix-vector multiplication to implement the gradient operator. The gradient w of the function u is approximated by the finite element method using the same function space as u, i.e. we solve w = grad(u), which corresponds to the weak form w*v*dx = grad(u)*v*dx, where v is a suitable test function.
- Parameters:
symbol (
pybamm.Symbol
) – The symbol that we will take the Laplacian of.discretised_symbol (
pybamm.Symbol
) – The discretised symbol of the correct sizeboundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
- Returns:
A concatenation that contains the result of acting the discretised gradient on the child discretised_symbol. The first column corresponds to the y-component of the gradient and the second column corresponds to the z component of the gradient.
- Return type:
class: pybamm.Concatenation
- gradient_matrix(symbol, boundary_conditions)[source]#
Gradient matrix for finite elements in the appropriate domain.
- Parameters:
symbol (
pybamm.Symbol
) – The symbol for which we want to calculate the gradient matrixboundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
- Returns:
The (sparse) finite element gradient matrix for the domain
- Return type:
- gradient_squared(symbol, discretised_symbol, boundary_conditions)[source]#
Multiplication to implement the inner product of the gradient operator with itself. See
pybamm.SpatialMethod.gradient_squared()
- indefinite_integral(child, discretised_child, direction)[source]#
Implementation of the indefinite integral operator. The input discretised child must be defined on the internal mesh edges. See
pybamm.SpatialMethod.indefinite_integral()
- integral(child, discretised_child, integration_dimension)[source]#
Vector-vector dot product to implement the integral operator. See
pybamm.SpatialMethod.integral()
- laplacian(symbol, discretised_symbol, boundary_conditions)[source]#
Matrix-vector multiplication to implement the Laplacian operator.
- Parameters:
symbol (
pybamm.Symbol
) – The symbol that we will take the Laplacian of.discretised_symbol (
pybamm.Symbol
) – The discretised symbol of the correct sizeboundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
- Returns:
Contains the result of acting the discretised gradient on the child discretised_symbol
- Return type:
class: pybamm.Array
- mass_matrix(symbol, boundary_conditions)[source]#
Calculates the mass matrix for the finite element method.
- Parameters:
symbol (
pybamm.Variable
) – The variable corresponding to the equation for which we are calculating the mass matrix.boundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
- Returns:
The (sparse) mass matrix for the spatial method.
- Return type:
- spatial_variable(symbol)[source]#
Creates a discretised spatial variable compatible with the FiniteElement method.
- Parameters:
symbol (
pybamm.SpatialVariable
) – The spatial variable to be discretised.- Returns:
Contains the discretised spatial variable
- Return type:
- stiffness_matrix(symbol, boundary_conditions)[source]#
Laplacian (stiffness) matrix for finite elements in the appropriate domain.
- Parameters:
symbol (
pybamm.Symbol
) – The symbol for which we want to calculate the Laplacian matrixboundary_conditions (dict) – The boundary conditions of the model ({symbol: {“negative tab”: neg. tab bc, “positive tab”: pos. tab bc}})
- Returns:
The (sparse) finite element stiffness matrix for the domain
- Return type: