#492 get direct voltage model working

Author: valentinsulzer

examples/scripts/DFN.py

@@ -8,13 +8,14 @@
 pybamm.set_logging_level("INFO")
 
 # load model
-model = pybamm.lithium_ion.DFN()
+model = pybamm.lithium_ion.DFN({"operating mode": "voltage"})
 
 # create geometry
 geometry = model.default_geometry
 
 # load parameter values and process model and geometry
 param = model.default_parameter_values
+param["Voltage function"] = 4.1
 param.process_model(model)
 param.process_geometry(geometry)
 

pybamm/models/full_battery_models/base_battery_model.py

@@ -534,13 +534,13 @@ def set_external_circuit_submodel(self):
                 self.param
             )
         elif self.options["operating mode"] == "voltage":
-            self.submodels["external circuit"] = pybamm.external_circuit.VoltageControl(
-                self.param
-            )
+            self.submodels[
+                "external circuit"
+            ] = pybamm.external_circuit.VoltageFunctionControl(self.param)
         elif self.options["operating mode"] == "power":
-            self.submodels["external circuit"] = pybamm.external_circuit.PowerControl(
-                self.param
-            )
+            self.submodels[
+                "external circuit"
+            ] = pybamm.external_circuit.PowerFunctionControl(self.param)
         elif callable(self.options["operating mode"]):
             self.submodels[
                 "external circuit"

pybamm/models/full_battery_models/lead_acid/full.py

@@ -69,7 +69,9 @@ def set_interfacial_submodel(self):
     def set_solid_submodel(self):
         if self.options["surface form"] is False:
             submod_n = pybamm.electrode.ohm.Full(self.param, "Negative", self.reactions)
-            submod_p = pybamm.electrode.ohm.Full(self.param, "Positive", self.reactions)
+            submod_p = pybamm.electrode.ohm.Full(
+                self.param, "Positive", self.reactions, self.options["operating mode"]
+            )
         else:
             submod_n = pybamm.electrode.ohm.SurfaceForm(self.param, "Negative")
             submod_p = pybamm.electrode.ohm.SurfaceForm(self.param, "Positive")
@@ -124,3 +126,12 @@ def set_side_reaction_submodels(self):
             self.submodels[
                 "negative oxygen interface"
             ] = pybamm.interface.lead_acid_oxygen.NoReaction(self.param, "Negative")
+
+    def set_external_circuit_submodel(self):
+        """ See :meth:`BaseBatteryModel.set_external_circuit_submodel` """
+        if self.options["operating mode"] == "voltage":
+            self.submodels["external circuit"] = pybamm.external_circuit.VoltageControl(
+                self.param
+            )
+        else:
+            super().set_external_circuit_submodel()

pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py

@@ -19,19 +19,6 @@ def __init__(self, options=None, name="Unnamed lithium-ion model"):
 
     def set_standard_output_variables(self):
         super().set_standard_output_variables()
-        # Current
-        i_cell = pybamm.standard_parameters_lithium_ion.current_with_time
-        i_cell_dim = (
-            pybamm.standard_parameters_lithium_ion.dimensional_current_density_with_time
-        )
-        I = pybamm.standard_parameters_lithium_ion.dimensional_current_with_time
-        self.variables.update(
-            {
-                "Total current density": i_cell,
-                "Total current density [A.m-2]": i_cell_dim,
-                "Current [A]": I,
-            }
-        )
 
         # Time
         time_scale = pybamm.standard_parameters_lithium_ion.tau_discharge
@@ -40,7 +27,6 @@ def set_standard_output_variables(self):
                 "Time [s]": pybamm.t * time_scale,
                 "Time [min]": pybamm.t * time_scale / 60,
                 "Time [h]": pybamm.t * time_scale / 3600,
-                "Discharge capacity [A.h]": I * pybamm.t * time_scale / 3600,
             }
         )
 

pybamm/models/full_battery_models/lithium_ion/dfn.py

@@ -88,7 +88,7 @@ def set_solid_submodel(self):
             self.param, "Negative", self.reactions
         )
         self.submodels["positive electrode"] = pybamm.electrode.ohm.Full(
-            self.param, "Positive", self.reactions
+            self.param, "Positive", self.reactions, self.options["operating mode"]
         )
 
     def set_electrolyte_submodel(self):
@@ -102,6 +102,15 @@ def set_electrolyte_submodel(self):
             self.param, self.reactions
         )
 
+    def set_external_circuit_submodel(self):
+        """ See :meth:`BaseBatteryModel.set_external_circuit_submodel` """
+        if self.options["operating mode"] == "voltage":
+            self.submodels["external circuit"] = pybamm.external_circuit.VoltageControl(
+                self.param
+            )
+        else:
+            super().set_external_circuit_submodel()
+
     @property
     def default_geometry(self):
         dimensionality = self.options["dimensionality"]

pybamm/models/submodels/electrode/base_electrode.py

@@ -134,8 +134,6 @@ def _get_standard_current_collector_potential_variables(self, phi_s_cn, phi_s_cp
         # Voltage is local current collector potential difference at the tabs, in 1D
         # this will be equal to the local current collector potential difference
         phi_s_cp_dim = U_ref + phi_s_cp * pot_scale
-        V = pybamm.boundary_value(phi_s_cp, "positive tab")
-        V_dim = pybamm.boundary_value(phi_s_cp_dim, "positive tab")
 
         variables = {
             "Negative current collector potential": phi_s_cn,
@@ -144,8 +142,6 @@ def _get_standard_current_collector_potential_variables(self, phi_s_cn, phi_s_cp
             "Positive current collector potential [V]": phi_s_cp_dim,
             "Local voltage": V_cc,
             "Local voltage [V]": U_ref + V_cc * pot_scale,
-            "Terminal voltage": V,
-            "Terminal voltage [V]": V_dim,
         }
 
         return variables
@@ -191,5 +187,6 @@ def _get_standard_whole_cell_variables(self, variables):
                     phi_s_cn, phi_s_cp
                 )
             )
+
         return variables
 

pybamm/models/submodels/electrode/ohm/full_ohm.py

@@ -19,8 +19,9 @@ class Full(BaseModel):
     **Extends:** :class:`pybamm.electrode.ohm.BaseModel`
     """
 
-    def __init__(self, param, domain, reactions):
+    def __init__(self, param, domain, reactions, operating_mode="current"):
         super().__init__(param, domain, reactions)
+        self.operating_mode = operating_mode
 
     def get_fundamental_variables(self):
 
@@ -71,19 +72,23 @@ def set_boundary_conditions(self, variables):
         phi_s = variables[self.domain + " electrode potential"]
         phi_s_cn = variables["Negative current collector potential"]
         tor = variables[self.domain + " electrode tortuosity"]
-        i_boundary_cc = variables["Current collector current density"]
 
         if self.domain == "Negative":
             lbc = (phi_s_cn, "Dirichlet")
             rbc = (pybamm.Scalar(0), "Neumann")
 
         elif self.domain == "Positive":
             lbc = (pybamm.Scalar(0), "Neumann")
-            sigma_eff = self.param.sigma_p * tor
-            rbc = (
-                i_boundary_cc / pybamm.boundary_value(-sigma_eff, "right"),
-                "Neumann",
-            )
+            if self.operating_mode == "voltage":
+                V = variables["Terminal voltage"]
+                rbc = (V, "Dirichlet")
+            else:
+                sigma_eff = self.param.sigma_p * tor
+                i_boundary_cc = variables["Current collector current density"]
+                rbc = (
+                    i_boundary_cc / pybamm.boundary_value(-sigma_eff, "right"),
+                    "Neumann",
+                )
 
         self.boundary_conditions[phi_s] = {"left": lbc, "right": rbc}
 

pybamm/models/submodels/external_circuit/__init__.py

@@ -1,6 +1,9 @@
 from .base_external_circuit import BaseModel
 from .current_control_external_circuit import CurrentControl
-from .function_control_external_circuit import FunctionControl
+from .function_control_external_circuit import (
+    FunctionControl,
+    VoltageFunctionControl,
+    PowerFunctionControl,
+)
 from .voltage_control_external_circuit import VoltageControl
-from .power_control_external_circuit import PowerControl
 

pybamm/models/submodels/external_circuit/current_control_external_circuit.py

@@ -25,3 +25,15 @@ def get_fundamental_variables(self):
         }
 
         return variables
+
+    def get_coupled_variables(self, variables):
+        # Update terminal voltage
+        phi_s_cp_dim = variables["Positive current collector potential [V]"]
+        phi_s_cp = variables["Positive current collector potential"]
+
+        V = pybamm.boundary_value(phi_s_cp, "positive tab")
+        V_dim = pybamm.boundary_value(phi_s_cp_dim, "positive tab")
+        variables["Terminal voltage"] = V
+        variables["Terminal voltage [V]"] = V_dim
+
+        return variables

pybamm/models/submodels/external_circuit/function_control_external_circuit.py

@@ -35,6 +35,16 @@ def get_fundamental_variables(self):
 
         return variables
 
+    def get_coupled_variables(self, variables):
+        # Update terminal voltage
+        phi_s_cp_dim = variables["Positive current collector potential [V]"]
+        phi_s_cp = variables["Positive current collector potential"]
+
+        V = pybamm.boundary_value(phi_s_cp, "positive tab")
+        V_dim = pybamm.boundary_value(phi_s_cp_dim, "positive tab")
+        variables["Terminal voltage"] = V
+        variables["Terminal voltage [V]"] = V_dim
+
     def set_initial_conditions(self, variables):
         # Initial condition as a guess for consistent initial conditions
         i_cell = variables["Total current density"]
@@ -46,3 +56,37 @@ def set_algebraic(self, variables):
         # or a combination (e.g. I*V for power control)
         i_cell = variables["Total current density"]
         self.algebraic[i_cell] = self.external_circuit_class(variables)
+
+
+class VoltageFunctionControl(FunctionControl):
+    """
+    External circuit with voltage control, implemented as an extra algebraic equation.
+    """
+
+    def __init__(self, param):
+        super().__init__(param, ConstantVoltage())
+
+
+class ConstantVoltage:
+    num_switches = 0
+
+    def __call__(self, variables):
+        V = variables["Terminal voltage [V]"]
+        return V - pybamm.FunctionParameter("Voltage function", pybamm.t)
+
+
+class PowerFunctionControl(FunctionControl):
+    """External circuit with power control. """
+
+    def __init__(self, param):
+        super().__init__(param, ConstantPower())
+
+
+class ConstantPower:
+    num_switches = 0
+
+    def __call__(self, variables):
+        I = variables["Current [A]"]
+        V = variables["Terminal voltage [V]"]
+        return I * V - pybamm.FunctionParameter("Power function", pybamm.t)
+

pybamm/models/submodels/external_circuit/power_control_external_circuit.py

@@ -2,20 +2,47 @@
 # External circuit with voltage control
 #
 import pybamm
-from .function_control_external_circuit import FunctionControl
+from .base_external_circuit import BaseModel
 
 
-class VoltageControl(FunctionControl):
-    """External circuit with voltage control. """
+class VoltageControl(BaseModel):
+    """
+    External circuit with voltage control, implemented directly.
+    """
 
     def __init__(self, param):
-        super().__init__(param, ConstantVoltage())
+        super().__init__(param)
 
+    def get_fundamental_variables(self):
+        # Voltage is given as a function of time
+        V_dim = pybamm.FunctionParameter("Voltage function", pybamm.t)
 
-class ConstantVoltage:
-    num_switches = 0
+        param = self.param
+        V = (V_dim - (param.U_p_ref - param.U_n_ref)) / param.potential_scale
 
-    def __call__(self, variables):
-        V = variables["Terminal voltage [V]"]
-        return V - pybamm.FunctionParameter("Voltage function", pybamm.t)
+        variables = {"Terminal voltage [V]": V_dim, "Terminal voltage": V}
 
+        return variables
+
+    def get_coupled_variables(self, variables):
+        # Update currrent
+        phi_s_p = variables["Positive electrode potential"]
+        tor = variables["Positive electrode tortuosity"]
+
+        param = self.param
+        i_boundary_cc = (
+            -param.sigma_p
+            * pybamm.BoundaryValue(tor, "right")
+            * pybamm.BoundaryGradient(phi_s_p, "right")
+        )
+        i_cell = pybamm.BoundaryValue(i_boundary_cc, "right")
+        I = i_cell * abs(param.I_typ)
+        i_cell_dim = I / (param.n_electrodes_parallel * param.A_cc)
+
+        variables = {
+            "Total current density": i_cell,
+            "Total current density [A.m-2]": i_cell_dim,
+            "Current [A]": I,
+        }
+
+        return variables