#709 allow current to be input parameter

Author: valentinsulzer

examples/notebooks/change-input-current.ipynb

@@ -46,10 +46,7 @@
     "param = model.default_parameter_values\n",
     "\n",
     "# change the current function to be an input parameter\n",
-    "def current_function(t):\n",
-    "    return pybamm.InputParameter(\"current\")\n",
-    "\n",
-    "param[\"Current function [A]\"] = current_function"
+    "param[\"Current function [A]\"] = \"[input]\""
    ]
   },
   {
@@ -69,7 +66,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "839ef921550749ecaae748a2068ee9cc",
+       "model_id": "770ba7f159054c5bad11038ac7be3f72",
        "version_major": 2,
        "version_minor": 0
       },
@@ -97,7 +94,7 @@
     "solver = model.default_solver\n",
     "n = 300\n",
     "t_eval = np.linspace(0, 0.02, n)\n",
-    "solution = solver.solve(model, t_eval, inputs={\"current\": 16})\n",
+    "solution = solver.solve(model, t_eval, inputs={\"Current function [A]\": 16})\n",
     "\n",
     "# plot\n",
     "quick_plot = pybamm.QuickPlot(solution)\n",
@@ -121,7 +118,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "cfb285a7877c4853ae7144891b38e4cd",
+       "model_id": "4625a91d401b456d8e43aae201f4a380",
        "version_major": 2,
        "version_minor": 0
       },
@@ -135,7 +132,7 @@
    ],
    "source": [
     "# Solve the model at the given time points\n",
-    "solution = solver.solve(model, t_eval, inputs={\"current\": 0})\n",
+    "solution = solver.solve(model, t_eval, inputs={\"Current function [A]\": 0})\n",
     "\n",
     "# plot\n",
     "quick_plot = pybamm.QuickPlot(solution)\n",
@@ -161,7 +158,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "d31307d03a2b4082a5ad7b4d91f85823",
+       "model_id": "88c9c957b666420e9cf6f615ff43817c",
        "version_major": 2,
        "version_minor": 0
       },
@@ -285,7 +282,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "bb1a66bff842443cb81abb44a4b4a6a1",
+       "model_id": "a8b90dfe860644c68af15c7c07e3adeb",
        "version_major": 2,
        "version_minor": 0
       },

examples/notebooks/change-settings.ipynb

@@ -104,7 +104,7 @@
     {
      "data": {
       "text/plain": [
-       "<pybamm.models.full_battery_models.lithium_ion.dfn.DFN at 0x10597e518>"
+       "<pybamm.models.full_battery_models.lithium_ion.dfn.DFN at 0x109dcd710>"
       ]
      },
      "execution_count": 4,
@@ -153,11 +153,7 @@
    "outputs": [],
    "source": [
     "param = dfn.default_parameter_values\n",
-    "\n",
-    "def current_function(t):\n",
-    "    return pybamm.InputParameter(\"current\")\n",
-    "\n",
-    "param[\"Current function [A]\"] = current_function"
+    "param[\"Current function [A]\"] = \"[input]\""
    ]
   },
   {
@@ -253,9 +249,9 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "Solved the Doyle-Fuller-Newman model in 2.111 seconds\n",
-      "Solved the Single Particle Model in 0.386 seconds\n",
-      "Solved the Single Particle Model with electrolyte in 0.399 seconds\n"
+      "Solved the Doyle-Fuller-Newman model in 2.609 seconds\n",
+      "Solved the Single Particle Model in 0.462 seconds\n",
+      "Solved the Single Particle Model with electrolyte in 0.491 seconds\n"
      ]
     }
    ],
@@ -266,7 +262,7 @@
     "solver = pybamm.CasadiSolver()\n",
     "for model_name, model in models.items():\n",
     "    start = timer.time()\n",
-    "    solution = solver.solve(model, t_eval, inputs={\"current\": 1})\n",
+    "    solution = solver.solve(model, t_eval, inputs={\"Current function [A]\": 1})\n",
     "    end = timer.time()\n",
     "    print(\"Solved the {} in {:.3f} seconds\".format(model.name, end-start))\n",
     "    solutions[model_name] = solution"
@@ -349,7 +345,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "af31d03838e4456aac24954eb1b5d0df",
+       "model_id": "36cafdf1324c44b19fa3bf54a8e579cb",
        "version_major": 2,
        "version_minor": 0
       },
@@ -371,7 +367,7 @@
    "cell_type": "markdown",
    "metadata": {},
    "source": [
-    "# A note on changing parameters:"
+    "# Changing parameters"
    ]
   },
   {
@@ -389,7 +385,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "7226178875b24a0eacb8c38ab1c4fafc",
+       "model_id": "fc6e5fae66414fdda8e177864e49cb02",
        "version_major": 2,
        "version_minor": 0
       },
@@ -404,7 +400,7 @@
    "source": [
     "# update parameter values and solve again\n",
     "for model_name, model in models.items():\n",
-    "    solutions[model_name] = model.default_solver.solve(model, t_eval, inputs={\"current\": 5})\n",
+    "    solutions[model_name] = model.default_solver.solve(model, t_eval, inputs={\"Current function [A]\": 5})\n",
     "\n",
     "# Plot\n",
     "list_of_models = list(models.values())\n",

examples/notebooks/compare-comsol-discharge-curve.ipynb

@@ -219,10 +219,8 @@
     "models = [loqs, composite, full]\n",
     "\n",
     "# process parameters\n",
-    "def current_function(t):\n",
-    "    return pybamm.InputParameter(\"current\")\n",
     "param = models[0].default_parameter_values\n",
-    "param[\"Current function [A]\"] = current_function\n",
+    "param[\"Current function [A]\"] = \"[input]\"\n",
     "for model in models:\n",
     "    param.process_model(model)"
    ]
@@ -267,9 +265,9 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
-      "Solved the LOQS model in 0.133 seconds\n",
-      "Solved the Composite model in 1.034 seconds\n",
-      "Solved the Full model in 1.187 seconds\n"
+      "Solved the LOQS model in 0.127 seconds\n",
+      "Solved the Composite model in 1.075 seconds\n",
+      "Solved the Full model in 1.143 seconds\n"
      ]
     }
    ],
@@ -280,7 +278,7 @@
     "solver = pybamm.CasadiSolver()\n",
     "for model in models:\n",
     "    start = timer.time()\n",
-    "    solution = solver.solve(model, t_eval, inputs={\"current\": 1})\n",
+    "    solution = solver.solve(model, t_eval, inputs={\"Current function [A]\": 1})\n",
     "    end = timer.time()\n",
     "    print(\"Solved the {} in {:.3f} seconds\".format(model.name, end-start))\n",
     "    solutions[model] = solution"
@@ -345,7 +343,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "6145d7740fa04796aa4dae9bb09c125e",
+       "model_id": "a08ab3900db7448f81089148ae3ee749",
        "version_major": 2,
        "version_minor": 0
       },
@@ -381,7 +379,7 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "3578ee915f79432c931ffecf4b543c1b",
+       "model_id": "71e4761492ee47d6836cb51fd0882812",
        "version_major": 2,
        "version_minor": 0
       },
@@ -395,7 +393,7 @@
    ],
    "source": [
     "for model in models:\n",
-    "    solutions[model] = solver.solve(model, t_eval, inputs={\"current\": 20})\n",
+    "    solutions[model] = solver.solve(model, t_eval, inputs={\"Current function [A]\": 20})\n",
     "\n",
     "# Plot\n",
     "solution_values = [solutions[model] for model in models]\n",

examples/notebooks/models/compare-lithium-ion.ipynb

@@ -16,14 +16,10 @@
 
 
 # load parameters and process model and geometry
-def current_function(t):
-    return pybamm.InputParameter("current")
-
-
 param = model.default_parameter_values
 param["Electrode width [m]"] = 1
 param["Electrode height [m]"] = 1
-param["Current function [A]"] = current_function
+param["Current function [A]"] = "[input]"
 param.process_model(model)
 param.process_geometry(geometry)
 
@@ -72,7 +68,7 @@ def current_function(t):
     # solve model at comsol times
     t = comsol_time / tau
     solution = pybamm.CasadiSolver(mode="fast").solve(
-        model, t, inputs={"current": current}
+        model, t, inputs={"Current function [A]": current}
     )
 
     # discharge capacity

examples/notebooks/models/lead-acid.ipynb

@@ -278,6 +278,8 @@ def check_and_update_parameter_values(self, values):
                     data = values["C-rate"][1]
                     data[:, 1] = data[:, 1] * capacity
                     value = (values["C-rate"][0] + "_to_Crate", data)
+                elif values["C-rate"] == "[input]":
+                    value = CrateToCurrent(values["C-rate"], capacity, typ="input")
                 else:
                     value = values["C-rate"] * capacity
                 self._dict_items["Current function [A]"] = value
@@ -288,6 +290,10 @@ def check_and_update_parameter_values(self, values):
                     data = values["Current function [A]"][1]
                     data[:, 1] = data[:, 1] / capacity
                     value = (values["Current function [A]"][0] + "_to_current", data)
+                elif values["Current function [A]"] == "[input]":
+                    value = CurrentToCrate(
+                        values["Current function [A]"], capacity, typ="input"
+                    )
                 else:
                     value = values["Current function [A]"] / capacity
                 self._dict_items["C-rate"] = value
@@ -472,6 +478,9 @@ def _process_symbol(self, symbol):
                 function = pybamm.Scalar(
                     function_name, name=symbol.name
                 ) * pybamm.ones_like(*new_children)
+            elif isinstance(function_name, pybamm.InputParameter):
+                # Replace the function with an input parameter
+                function = function_name
             else:
                 # otherwise evaluate the function to create a new PyBaMM object
                 function = function_name(*new_children)
@@ -547,20 +556,28 @@ def evaluate(self, symbol):
 class CurrentToCrate:
     "Convert a current function to a C-rate function"
 
-    def __init__(self, function, capacity):
-        self.function = function
+    def __init__(self, current, capacity, typ="function"):
+        self.current = current
         self.capacity = capacity
+        self.type = typ
 
     def __call__(self, t):
-        return self.function(t) / self.capacity
+        if self.type == "function":
+            return self.current(t) / self.capacity
+        elif self.type == "input":
+            return pybamm.InputParameter("Current function [A]") / self.capacity
 
 
 class CrateToCurrent:
     "Convert a C-rate function to a current function"
 
-    def __init__(self, function, capacity):
-        self.function = function
+    def __init__(self, Crate, capacity, typ="function"):
+        self.Crate = Crate
         self.capacity = capacity
+        self.type = typ
 
     def __call__(self, t):
-        return self.function(t) * self.capacity
+        if self.type == "function":
+            return self.Crate(t) * self.capacity
+        elif self.type == "input":
+            return pybamm.InputParameter("C-rate") * self.capacity

examples/scripts/compare_comsol/discharge_curve.py

@@ -18,12 +18,9 @@ def test_leading_order_convergence(self):
         composite_model = pybamm.lead_acid.Composite()
         full_model = pybamm.lead_acid.Full()
 
-        def current_function(t):
-            return pybamm.InputParameter("Current")
-
         # Same parameters, same geometry
         parameter_values = full_model.default_parameter_values
-        parameter_values["Current function [A]"] = current_function
+        parameter_values["Current function [A]"] = "[input]"
         parameter_values.process_model(leading_order_model)
         parameter_values.process_model(composite_model)
         parameter_values.process_model(full_model)
@@ -51,23 +48,17 @@ def get_max_error(current):
             pybamm.logger.info("current = {}".format(current))
             # Solve, make sure times are the same and use tight tolerances
             t_eval = np.linspace(0, 0.6)
-            solver_loqs = leading_order_model.default_solver
-            solver_loqs.rtol = 1e-8
-            solver_loqs.atol = 1e-8
-            solution_loqs = solver_loqs.solve(
-                leading_order_model, t_eval, inputs={"Current": current}
+            solver = pybamm.CasadiSolver()
+            solver.rtol = 1e-8
+            solver.atol = 1e-8
+            solution_loqs = solver.solve(
+                leading_order_model, t_eval, inputs={"Current function [A]": current}
             )
-            solver_comp = composite_model.default_solver
-            solver_comp.rtol = 1e-8
-            solver_comp.atol = 1e-8
-            solution_comp = solver_comp.solve(
-                composite_model, t_eval, inputs={"Current": current}
+            solution_comp = solver.solve(
+                composite_model, t_eval, inputs={"Current function [A]": current}
             )
-            solver_full = full_model.default_solver
-            solver_full.rtol = 1e-8
-            solver_full.atol = 1e-8
-            solution_full = solver_full.solve(
-                full_model, t_eval, inputs={"Current": current}
+            solution_full = solver.solve(
+                full_model, t_eval, inputs={"Current function [A]": current}
             )
 
             # Post-process variables
@@ -105,5 +96,4 @@ def get_max_error(current):
 
     if "-v" in sys.argv:
         debug = True
-    pybamm.set_logging_level("DEBUG")
     unittest.main()

pybamm/parameters/parameter_values.py

@@ -115,6 +115,18 @@ def test_check_and_update_parameter_values(self):
             param["C-rate"][1], np.hstack([x, 0.2 * x])
         )
 
+        # With input parameters
+        # if only C-rate and capacity provided, update current
+        values = {"C-rate": "[input]", "Cell capacity [A.h]": 10}
+        param = pybamm.ParameterValues(values)
+        self.assertEqual(param["Current function [A]"](2).evaluate(u={"C-rate": 1}), 10)
+        # if only current and capacity provided, update C-rate
+        values = {"Current function [A]": "[input]", "Cell capacity [A.h]": 10}
+        param = pybamm.ParameterValues(values)
+        self.assertEqual(
+            param["C-rate"](5).evaluate(u={"Current function [A]": 5}), 0.5
+        )
+
     def test_process_symbol(self):
         parameter_values = pybamm.ParameterValues({"a": 1, "b": 2, "c": 3})
         # process parameter
@@ -297,6 +309,16 @@ def test_process_function_parameter(self):
         processed_diff_func = parameter_values.process_symbol(diff_func)
         self.assertEqual(processed_diff_func.evaluate(u={"a": 3}), 123)
 
+        # function itself as input (different to the variable being an input)
+        parameter_values = pybamm.ParameterValues({"func": "[input]"})
+        a = pybamm.Scalar(3)
+        func = pybamm.FunctionParameter("func", a)
+        processed_func = parameter_values.process_symbol(func)
+        self.assertEqual(processed_func.evaluate(u={"func": 13}), 13)
+
+        
+
+
     def test_process_inline_function_parameters(self):
         def D(c):
             return c ** 2