diff --git a/CHANGELOG.md b/CHANGELOG.md index 1b8e1a869b..5204e0bc82 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,8 @@ ## Features +- Added a new unary operator, `EvaluateAt`, that evaluates a spatial variable at a given position ([#3573](https://github.com/pybamm-team/PyBaMM/pull/3573)) +- Added a method, `insert_reference_electrode`, to `pybamm.lithium_ion.BaseModel` that insert a reference electrode to measure the electrolyte potential at a given position in space and adds new variables that mimic a 3E cell setup. ([#3573](https://github.com/pybamm-team/PyBaMM/pull/3573)) - Serialisation added so models can be written to/read from JSON ([#3397](https://github.com/pybamm-team/PyBaMM/pull/3397)) - Mechanical parameters are now a function of stoichiometry and temperature ([#3576](https://github.com/pybamm-team/PyBaMM/pull/3576)) diff --git a/docs/source/api/expression_tree/unary_operator.rst b/docs/source/api/expression_tree/unary_operator.rst index ad5bb0a48f..e6a3cbe554 100644 --- a/docs/source/api/expression_tree/unary_operator.rst +++ b/docs/source/api/expression_tree/unary_operator.rst @@ -34,6 +34,9 @@ Unary Operators .. autoclass:: pybamm.Mass :members: +.. autoclass:: pybamm.BoundaryMass + :members: + .. autoclass:: pybamm.Integral :members: @@ -58,6 +61,9 @@ Unary Operators .. autoclass:: pybamm.BoundaryGradient :members: +.. autoclass:: pybamm.EvaluateAt + :members: + .. autoclass:: pybamm.UpwindDownwind :members: diff --git a/docs/source/examples/index.rst b/docs/source/examples/index.rst index 4afaa6eeeb..e025ea71b4 100644 --- a/docs/source/examples/index.rst +++ b/docs/source/examples/index.rst @@ -65,6 +65,7 @@ The notebooks are organised into subfolders, and can be viewed in the galleries notebooks/models/rate-capability.ipynb notebooks/models/saving_models.ipynb notebooks/models/SEI-on-cracks.ipynb + notebooks/models/simulate-3E-cell.ipynb notebooks/models/simulating-ORegan-2022-parameter-set.ipynb notebooks/models/SPM.ipynb notebooks/models/SPMe.ipynb diff --git a/docs/source/examples/notebooks/models/simulate-3E-cell.ipynb b/docs/source/examples/notebooks/models/simulate-3E-cell.ipynb new file mode 100644 index 0000000000..c92cb53465 --- /dev/null +++ b/docs/source/examples/notebooks/models/simulate-3E-cell.ipynb @@ -0,0 +1,210 @@ +{ + "cells": [ + { + "attachments": {}, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Simulating a 3E cell\n", + "\n", + "In this notebook we show how to insert a reference electrode to mimic a 3E cell." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Note: you may need to restart the kernel to use updated packages.\n" + ] + } + ], + "source": [ + "%pip install \"pybamm[plot,cite]\" -q # install PyBaMM if it is not installed\n", + "import pybamm" + ] + }, + { + "attachments": {}, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We first load a model" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "model = pybamm.lithium_ion.DFN()" + ] + }, + { + "attachments": {}, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Next we use the helper function `insert_reference_electrode` to insert a reference electrode into the model. This function takes the position of the reference electrode as an optional argument. If no position is given, the reference electrode is inserted at the midpoint of the separator. The helper function adds the new variables \"Reference electrode potential [V]\", \"Negative electrode 3E potential [V]\" and \"Positive electrode 3E potential [V]\" to the model.\n", + "\n", + "In this example we will explicitly pass a position to show how it is done" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "L_n = model.param.n.L # Negative electrode thickness [m]\n", + "L_s = model.param.s.L # Separator thickness [m]\n", + "L_ref = L_n + L_s / 2 # Reference electrode position [m]\n", + "\n", + "model.insert_reference_electrode(L_ref)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Next we can set up a simulation and solve the model as usual" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "sim = pybamm.Simulation(model)\n", + "sim.solve([0, 3600])" + ] + }, + { + "attachments": {}, + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's plot a comparison of the 3E potentials and the potential difference between the solid and electrolyte phases at the electrode/separator interfaces" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "application/vnd.jupyter.widget-view+json": { + "model_id": "d1f4e1ed03764660b87ee56b135b24a7", + "version_major": 2, + "version_minor": 0 + }, + "text/plain": [ + "interactive(children=(FloatSlider(value=0.0, description='t', max=1.0, step=0.01), Output()), _dom_classes=('w…" + ] + }, + "metadata": {}, + "output_type": "display_data" + }, + { + "data": { + "text/plain": [ + "" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "sim.plot(\n", + " [\n", + " [\n", + " \"Negative electrode surface potential difference at separator interface [V]\",\n", + " \"Negative electrode 3E potential [V]\",\n", + " ],\n", + " [\n", + " \"Positive electrode surface potential difference at separator interface [V]\",\n", + " \"Positive electrode 3E potential [V]\",\n", + " ],\n", + " \"Voltage [V]\",\n", + " ]\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "[1] Joel A. E. Andersson, Joris Gillis, Greg Horn, James B. Rawlings, and Moritz Diehl. CasADi – A software framework for nonlinear optimization and optimal control. Mathematical Programming Computation, 11(1):1–36, 2019. doi:10.1007/s12532-018-0139-4.\n", + "[2] Marc Doyle, Thomas F. Fuller, and John Newman. Modeling of galvanostatic charge and discharge of the lithium/polymer/insertion cell. Journal of the Electrochemical society, 140(6):1526–1533, 1993. doi:10.1149/1.2221597.\n", + "[3] Charles R. Harris, K. Jarrod Millman, Stéfan J. van der Walt, Ralf Gommers, Pauli Virtanen, David Cournapeau, Eric Wieser, Julian Taylor, Sebastian Berg, Nathaniel J. Smith, and others. Array programming with NumPy. Nature, 585(7825):357–362, 2020. doi:10.1038/s41586-020-2649-2.\n", + "[4] Scott G. Marquis, Valentin Sulzer, Robert Timms, Colin P. Please, and S. Jon Chapman. An asymptotic derivation of a single particle model with electrolyte. Journal of The Electrochemical Society, 166(15):A3693–A3706, 2019. doi:10.1149/2.0341915jes.\n", + "[5] Valentin Sulzer, Scott G. Marquis, Robert Timms, Martin Robinson, and S. Jon Chapman. Python Battery Mathematical Modelling (PyBaMM). Journal of Open Research Software, 9(1):14, 2021. doi:10.5334/jors.309.\n", + "\n" + ] + } + ], + "source": [ + "pybamm.print_citations()" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "venv", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.6" + }, + "orig_nbformat": 4, + "vscode": { + "interpreter": { + "hash": "9ff3d0c7e37de5f5aa47f4f719e4c84fc6cba7b39c571a05173422444e82fa58" + } + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/docs/source/examples/notebooks/models/simulating-ORegan-2022-parameter-set.ipynb b/docs/source/examples/notebooks/models/simulating-ORegan-2022-parameter-set.ipynb index 7eb647fc97..f20f385601 100644 --- a/docs/source/examples/notebooks/models/simulating-ORegan-2022-parameter-set.ipynb +++ b/docs/source/examples/notebooks/models/simulating-ORegan-2022-parameter-set.ipynb @@ -163,7 +163,7 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3.7.4 ('dev': venv)", + "display_name": "dev", "language": "python", "name": "python3" }, @@ -177,7 +177,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.4" + "version": "3.9.16" }, "toc": { "base_numbering": 1, @@ -194,7 +194,7 @@ }, "vscode": { "interpreter": { - "hash": "0f0e5a277ebcf03e91e138edc3d4774b5dee64e7d6640c0d876f99a9f6b2a4dc" + "hash": "bca2b99bfac80e18288b793d52fa0653ab9b5fe5d22e7b211c44eb982a41c00c" } } }, diff --git a/pybamm/__init__.py b/pybamm/__init__.py index cab7914cd9..f9b809d121 100644 --- a/pybamm/__init__.py +++ b/pybamm/__init__.py @@ -54,6 +54,7 @@ from .logger import logger, set_logging_level, get_new_logger from .settings import settings from .citations import Citations, citations, print_citations + # # Classes for the Expression Tree # diff --git a/pybamm/discretisations/discretisation.py b/pybamm/discretisations/discretisation.py index bb6e678f4c..62110b1676 100644 --- a/pybamm/discretisations/discretisation.py +++ b/pybamm/discretisations/discretisation.py @@ -865,6 +865,10 @@ def _process_symbol(self, symbol): return child_spatial_method.boundary_value_or_flux( symbol, disc_child, self.bcs ) + elif isinstance(symbol, pybamm.EvaluateAt): + return child_spatial_method.evaluate_at( + symbol, disc_child, symbol.position + ) elif isinstance(symbol, pybamm.UpwindDownwind): direction = symbol.name # upwind or downwind return spatial_method.upwind_or_downwind( diff --git a/pybamm/expression_tree/unary_operators.py b/pybamm/expression_tree/unary_operators.py index 2705685814..319429183c 100644 --- a/pybamm/expression_tree/unary_operators.py +++ b/pybamm/expression_tree/unary_operators.py @@ -325,7 +325,14 @@ def _unary_jac(self, child_jac): def set_id(self): """See :meth:`pybamm.Symbol.set_id()`""" self._id = hash( - (self.__class__, self.name, self.slice.start, self.slice.stop, self.children[0].id, *tuple(self.domain)) + ( + self.__class__, + self.name, + self.slice.start, + self.slice.stop, + self.children[0].id, + *tuple(self.domain), + ) ) def _unary_evaluate(self, child): @@ -465,7 +472,9 @@ def _unary_new_copy(self, child): def _sympy_operator(self, child): """Override :meth:`pybamm.UnaryOperator._sympy_operator`""" - sympy_Divergence = have_optional_dependency("sympy.vector.operators", "Divergence") + sympy_Divergence = have_optional_dependency( + "sympy.vector.operators", "Divergence" + ) return sympy_Divergence(child) @@ -616,7 +625,18 @@ def integration_variable(self): def set_id(self): """See :meth:`pybamm.Symbol.set_id()`""" self._id = hash( - (self.__class__, self.name, *tuple([integration_variable.id for integration_variable in self.integration_variable]), self.children[0].id, *tuple(self.domain)) + ( + self.__class__, + self.name, + *tuple( + [ + integration_variable.id + for integration_variable in self.integration_variable + ] + ), + self.children[0].id, + *tuple(self.domain), + ) ) def _unary_new_copy(self, child): @@ -757,7 +777,13 @@ def __init__(self, child, vector_type="row"): def set_id(self): """See :meth:`pybamm.Symbol.set_id()`""" self._id = hash( - (self.__class__, self.name, self.vector_type, self.children[0].id, *tuple(self.domain)) + ( + self.__class__, + self.name, + self.vector_type, + self.children[0].id, + *tuple(self.domain), + ) ) def _unary_new_copy(self, child): @@ -851,7 +877,13 @@ def __init__(self, child, side, domain): def set_id(self): """See :meth:`pybamm.Symbol.set_id()`""" self._id = hash( - (self.__class__, self.name, self.side, self.children[0].id, *tuple([(k, tuple(v)) for k, v in self.domains.items()])) + ( + self.__class__, + self.name, + self.side, + self.children[0].id, + *tuple([(k, tuple(v)) for k, v in self.domains.items()]), + ) ) def _evaluates_on_edges(self, dimension): @@ -872,7 +904,8 @@ def evaluate_for_shape(self): class BoundaryOperator(SpatialOperator): """ - A node in the expression tree which gets the boundary value of a variable. + A node in the expression tree which gets the boundary value of a variable on its + primary domain. Parameters ---------- @@ -909,7 +942,13 @@ def __init__(self, name, child, side): def set_id(self): """See :meth:`pybamm.Symbol.set_id()`""" self._id = hash( - (self.__class__, self.name, self.side, self.children[0].id, *tuple([(k, tuple(v)) for k, v in self.domains.items()])) + ( + self.__class__, + self.name, + self.side, + self.children[0].id, + *tuple([(k, tuple(v)) for k, v in self.domains.items()]), + ) ) def _unary_new_copy(self, child): @@ -923,7 +962,8 @@ def _evaluate_for_shape(self): class BoundaryValue(BoundaryOperator): """ - A node in the expression tree which gets the boundary value of a variable. + A node in the expression tree which gets the boundary value of a variable on its + primary domain. Parameters ---------- @@ -998,7 +1038,8 @@ def to_json(self): class BoundaryGradient(BoundaryOperator): """ - A node in the expression tree which gets the boundary flux of a variable. + A node in the expression tree which gets the boundary flux of a variable on its + primary domain. Parameters ---------- @@ -1012,6 +1053,59 @@ def __init__(self, child, side): super().__init__("boundary flux", child, side) +class EvaluateAt(SpatialOperator): + """ + A node in the expression tree which evaluates a symbol at a given position in space + in its primary domain. Currently this is only implemented for 1D primary domains. + + Parameters + ---------- + child : :class:`pybamm.Symbol` + The variable to evaluate + position : :class:`pybamm.Symbol` + The position in space on the symbol's primary domain at which to evaluate + the symbol. + """ + + def __init__(self, child, position): + self.position = position + + # "evaluate at" of a child takes the primary domain from secondary domain + # of the child + # tertiary auxiliary domain shift down to secondary, quarternary to tertiary + domains = { + "primary": child.domains["secondary"], + "secondary": child.domains["tertiary"], + "tertiary": child.domains["quaternary"], + } + + super().__init__("evaluate", child, domains) + + def set_id(self): + """See :meth:`pybamm.Symbol.set_id()`""" + self._id = hash( + ( + self.__class__, + self.name, + self.position, + self.children[0].id, + *tuple([(k, tuple(v)) for k, v in self.domains.items()]), + ) + ) + + def _unary_jac(self, child_jac): + """See :meth:`pybamm.UnaryOperator._unary_jac()`.""" + return pybamm.Scalar(0) + + def _unary_new_copy(self, child): + """See :meth:`UnaryOperator._unary_new_copy()`.""" + return self.__class__(child, self.position) + + def _evaluate_for_shape(self): + """See :meth:`pybamm.Symbol.evaluate_for_shape_using_domain()`""" + return pybamm.evaluate_for_shape_using_domain(self.domains) + + class UpwindDownwind(SpatialOperator): """ A node in the expression tree representing an upwinding or downwinding operator. diff --git a/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py b/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py index cc736d6d04..fbe19b0d42 100644 --- a/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py +++ b/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py @@ -462,3 +462,51 @@ def set_convection_submodel(self): self.submodels[ "through-cell convection" ] = pybamm.convection.through_cell.NoConvection(self.param, self.options) + + def insert_reference_electrode(self, position=None): + """ + Insert a reference electrode to measure the electrolyte potential at a given + position in space. Adds model variables for the electrolyte potential at the + reference electrode and for the potential difference between the electrode + potentials measured at the electrode/current collector interface and the + reference electrode. Only implemented for 1D models (i.e. where the + 'dimensionality' option is 0). + + Parameters + ---------- + position : :class:`pybamm.Symbol`, optional + The position in space at which to measure the electrolyte potential. If + None, defaults to the mid-point of the separator. + """ + if self.options["dimensionality"] != 0: + raise NotImplementedError( + "Reference electrode can only be inserted for models where " + "'dimensionality' is 0. For other models, please add a reference " + "electrode manually." + ) + + param = self.param + if position is None: + position = param.n.L + param.s.L / 2 + + phi_e_ref = pybamm.EvaluateAt( + self.variables["Electrolyte potential [V]"], position + ) + phi_p = pybamm.boundary_value( + self.variables["Positive electrode potential [V]"], "right" + ) + variables = { + "Positive electrode 3E potential [V]": phi_p - phi_e_ref, + "Reference electrode potential [V]": phi_e_ref, + } + + if self.options["working electrode"] == "both": + phi_n = pybamm.boundary_value( + self.variables["Negative electrode potential [V]"], "left" + ) + variables.update( + { + "Negative electrode 3E potential [V]": phi_n - phi_e_ref, + } + ) + self.variables.update(variables) diff --git a/pybamm/parameters/parameter_values.py b/pybamm/parameters/parameter_values.py index d5f12f362f..049910ae9e 100644 --- a/pybamm/parameters/parameter_values.py +++ b/pybamm/parameters/parameter_values.py @@ -721,6 +721,15 @@ def _process_symbol(self, symbol): # f_a_dist in the size average needs to be processed if isinstance(new_symbol, pybamm.SizeAverage): new_symbol.f_a_dist = self.process_symbol(new_symbol.f_a_dist) + # position in evaluate at needs to be processed, and should be a Scalar + if isinstance(new_symbol, pybamm.EvaluateAt): + new_symbol_position = self.process_symbol(new_symbol.position) + if not isinstance(new_symbol_position, pybamm.Scalar): + raise ValueError( + "'position' in 'EvaluateAt' must evaluate to a scalar" + ) + else: + new_symbol.position = new_symbol_position return new_symbol # Functions diff --git a/pybamm/spatial_methods/finite_volume.py b/pybamm/spatial_methods/finite_volume.py index 5c32e5a2c0..636243f829 100644 --- a/pybamm/spatial_methods/finite_volume.py +++ b/pybamm/spatial_methods/finite_volume.py @@ -1023,6 +1023,46 @@ def boundary_value_or_flux(self, symbol, discretised_child, bcs=None): return boundary_value + def evaluate_at(self, symbol, discretised_child, position): + """ + Returns the symbol evaluated at a given position in space. + + Parameters + ---------- + symbol: :class:`pybamm.Symbol` + The boundary value or flux symbol + discretised_child : :class:`pybamm.StateVector` + The discretised variable from which to calculate the boundary value + position : :class:`pybamm.Scalar` + The point in one-dimensional space at which to evaluate the symbol. + + Returns + ------- + :class:`pybamm.MatrixMultiplication` + The variable representing the value at the given point. + """ + # Get mesh nodes + domain = discretised_child.domain + mesh = self.mesh[domain] + nodes = mesh.nodes + repeats = self._get_auxiliary_domain_repeats(discretised_child.domains) + + # Find the index of the node closest to the value + index = np.argmin(np.abs(nodes - position.value)) + + # Create a sparse matrix with a 1 at the index + sub_matrix = csr_matrix(([1], ([0], [index])), shape=(1, mesh.npts)) + # repeat across auxiliary domains + matrix = csr_matrix(kron(eye(repeats), sub_matrix)) + + # Index into the discretised child + out = pybamm.Matrix(matrix) @ discretised_child + + # `EvaluateAt` removes domain + out.clear_domains() + + return out + def process_binary_operators(self, bin_op, left, right, disc_left, disc_right): """Discretise binary operators in model equations. Performs appropriate averaging of diffusivities if one of the children is a gradient operator, so diff --git a/pybamm/spatial_methods/spatial_method.py b/pybamm/spatial_methods/spatial_method.py index acb0227bc2..a461d6c150 100644 --- a/pybamm/spatial_methods/spatial_method.py +++ b/pybamm/spatial_methods/spatial_method.py @@ -331,7 +331,7 @@ def internal_neumann_condition( def boundary_value_or_flux(self, symbol, discretised_child, bcs=None): """ - Returns the boundary value or flux using the approriate expression for the + Returns the boundary value or flux using the appropriate expression for the spatial method. To do this, we create a sparse vector 'bv_vector' that extracts either the first (for side="left") or last (for side="right") point from 'discretised_child'. @@ -377,6 +377,26 @@ def boundary_value_or_flux(self, symbol, discretised_child, bcs=None): out.clear_domains() return out + def evaluate_at(self, symbol, discretised_child, position): + """ + Returns the symbol evaluated at a given position in space. + + Parameters + ---------- + symbol: :class:`pybamm.Symbol` + The boundary value or flux symbol + discretised_child : :class:`pybamm.StateVector` + The discretised variable from which to calculate the boundary value + position : :class:`pybamm.Scalar` + The point in one-dimensional space at which to evaluate the symbol. + + Returns + ------- + :class:`pybamm.MatrixMultiplication` + The variable representing the value at the given point. + """ + raise NotImplementedError + def mass_matrix(self, symbol, boundary_conditions): """ Calculates the mass matrix for a spatial method. diff --git a/tests/unit/test_expression_tree/test_operations/test_copy.py b/tests/unit/test_expression_tree/test_operations/test_copy.py index 0340e56bb1..6800f9092f 100644 --- a/tests/unit/test_expression_tree/test_operations/test_copy.py +++ b/tests/unit/test_expression_tree/test_operations/test_copy.py @@ -60,6 +60,7 @@ def test_symbol_new_copy(self): pybamm.maximum(a, b), pybamm.SparseStack(mat, mat), pybamm.Equality(a, b), + pybamm.EvaluateAt(a, 0), ]: self.assertEqual(symbol, symbol.new_copy()) diff --git a/tests/unit/test_expression_tree/test_operations/test_jac.py b/tests/unit/test_expression_tree/test_operations/test_jac.py index c6e04d331f..503e7321ea 100644 --- a/tests/unit/test_expression_tree/test_operations/test_jac.py +++ b/tests/unit/test_expression_tree/test_operations/test_jac.py @@ -236,6 +236,12 @@ def test_index(self): jac = ind.jac(vec).evaluate(y=np.linspace(0, 2, 5)).toarray() np.testing.assert_array_equal(jac, np.array([[0, 0, 0, 0, 0]])) + def test_evaluate_at(self): + y = pybamm.StateVector(slice(0, 4)) + expr = pybamm.EvaluateAt(y, 2) + jac = expr.jac(y).evaluate(y=np.linspace(0, 2, 4)) + np.testing.assert_array_equal(jac, 0) + def test_jac_of_number(self): """Jacobian of a number should be zero""" a = pybamm.Scalar(1) diff --git a/tests/unit/test_expression_tree/test_unary_operators.py b/tests/unit/test_expression_tree/test_unary_operators.py index 3c34db7dcd..6ae6b62d05 100644 --- a/tests/unit/test_expression_tree/test_unary_operators.py +++ b/tests/unit/test_expression_tree/test_unary_operators.py @@ -454,6 +454,11 @@ def test_index(self): pybamm.Index(vec, 5) pybamm.settings.debug_mode = debug_mode + def test_evaluate_at(self): + a = pybamm.Symbol("a", domain=["negative electrode"]) + f = pybamm.EvaluateAt(a, 1) + self.assertEqual(f.position, 1) + def test_upwind_downwind(self): # upwind of scalar symbol should fail a = pybamm.Symbol("a") @@ -673,9 +678,10 @@ def test_not_constant(self): self.assertFalse((2 * a).is_constant()) def test_to_equation(self): - sympy = have_optional_dependency("sympy") - sympy_Divergence = have_optional_dependency("sympy.vector.operators", "Divergence") + sympy_Divergence = have_optional_dependency( + "sympy.vector.operators", "Divergence" + ) sympy_Gradient = have_optional_dependency("sympy.vector.operators", "Gradient") a = pybamm.Symbol("a", domain="negative particle") diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_base_lithium_ion_model.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_base_lithium_ion_model.py index 315896b29f..fbc916d4a5 100644 --- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_base_lithium_ion_model.py +++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_base_lithium_ion_model.py @@ -29,6 +29,25 @@ def test_default_parameters(self): ) os.chdir(cwd) + def test_insert_reference_electrode(self): + model = pybamm.lithium_ion.SPM() + model.insert_reference_electrode() + self.assertIn("Negative electrode 3E potential [V]", model.variables) + self.assertIn("Positive electrode 3E potential [V]", model.variables) + self.assertIn("Reference electrode potential [V]", model.variables) + + model = pybamm.lithium_ion.SPM({"working electrode": "positive"}) + model.insert_reference_electrode() + self.assertNotIn("Negative electrode potential [V]", model.variables) + self.assertIn("Positive electrode 3E potential [V]", model.variables) + self.assertIn("Reference electrode potential [V]", model.variables) + + model = pybamm.lithium_ion.SPM({"dimensionality": 2}) + with self.assertRaisesRegex( + NotImplementedError, "Reference electrode can only be" + ): + model.insert_reference_electrode() + if __name__ == "__main__": print("Add -v for more debug output") diff --git a/tests/unit/test_parameters/test_parameter_values.py b/tests/unit/test_parameters/test_parameter_values.py index 3610b53424..7b5fd38bf0 100644 --- a/tests/unit/test_parameters/test_parameter_values.py +++ b/tests/unit/test_parameters/test_parameter_values.py @@ -267,6 +267,15 @@ def test_process_symbol(self): self.assertEqual(processed_a.value, 4) self.assertEqual(processed_x, x) + # process EvaluateAt + evaluate_at = pybamm.EvaluateAt(x, aa) + processed_evaluate_at = parameter_values.process_symbol(evaluate_at) + self.assertIsInstance(processed_evaluate_at, pybamm.EvaluateAt) + self.assertEqual(processed_evaluate_at.children[0], x) + self.assertEqual(processed_evaluate_at.position, 4) + with self.assertRaisesRegex(ValueError, "'position' in 'EvaluateAt'"): + parameter_values.process_symbol(pybamm.EvaluateAt(x, x)) + # process broadcast whole_cell = ["negative electrode", "separator", "positive electrode"] broad = pybamm.PrimaryBroadcast(a, whole_cell) diff --git a/tests/unit/test_spatial_methods/test_base_spatial_method.py b/tests/unit/test_spatial_methods/test_base_spatial_method.py index 37b4eb6a0b..d48ea69a7b 100644 --- a/tests/unit/test_spatial_methods/test_base_spatial_method.py +++ b/tests/unit/test_spatial_methods/test_base_spatial_method.py @@ -36,6 +36,8 @@ def test_basics(self): spatial_method.delta_function(None, None) with self.assertRaises(NotImplementedError): spatial_method.internal_neumann_condition(None, None, None, None) + with self.assertRaises(NotImplementedError): + spatial_method.evaluate_at(None, None, None) def test_get_auxiliary_domain_repeats(self): # Test the method to read number of repeats from auxiliary domains diff --git a/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py b/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py index 91a5b70044..16a3bbde2c 100644 --- a/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py +++ b/tests/unit/test_spatial_methods/test_finite_volume/test_finite_volume.py @@ -551,6 +551,30 @@ def test_full_broadcast_domains(self): disc = pybamm.Discretisation(mesh, spatial_methods) disc.process_model(model) + def test_evaluate_at(self): + mesh = get_p2d_mesh_for_testing() + spatial_methods = { + "macroscale": pybamm.FiniteVolume(), + "negative particle": pybamm.FiniteVolume(), + "positive particle": pybamm.FiniteVolume(), + } + disc = pybamm.Discretisation(mesh, spatial_methods) + + n = mesh["negative electrode"].npts + var = pybamm.StateVector(slice(0, n), domain="negative electrode") + + idx = 3 + position = pybamm.Scalar(mesh["negative electrode"].nodes[idx]) + evaluate_at = pybamm.EvaluateAt(var, position) + evaluate_at_disc = disc.process_symbol(evaluate_at) + + self.assertIsInstance(evaluate_at_disc, pybamm.MatrixMultiplication) + self.assertIsInstance(evaluate_at_disc.left, pybamm.Matrix) + self.assertIsInstance(evaluate_at_disc.right, pybamm.StateVector) + + y = np.arange(n)[:, np.newaxis] + self.assertEqual(evaluate_at_disc.evaluate(y=y), y[idx]) + if __name__ == "__main__": print("Add -v for more debug output")