From b1a9184232bce26cecf97a6ceaa094c982d580d0 Mon Sep 17 00:00:00 2001 From: Ali Hussain Umar Bhatti <65511923+alibh95@users.noreply.github.com> Date: Thu, 29 Apr 2021 01:45:54 +0500 Subject: [PATCH 1/2] Simulate Drive Cycle In Experiment This commit simulates drive cycle ('A' or 'V' or 'W') in experiment only if drive cycle of particular type is initial instruction in experiment and there are no further drive cycle related instructions. --- examples/scripts/experiment_drive_cycle.py | 87 +++++++++++++++ pybamm/experiments/experiment.py | 2 +- pybamm/simulation.py | 124 ++++++++++++++------- pybamm/solvers/base_solver.py | 9 ++ 4 files changed, 180 insertions(+), 42 deletions(-) create mode 100644 examples/scripts/experiment_drive_cycle.py diff --git a/examples/scripts/experiment_drive_cycle.py b/examples/scripts/experiment_drive_cycle.py new file mode 100644 index 0000000000..5ee7fbab2c --- /dev/null +++ b/examples/scripts/experiment_drive_cycle.py @@ -0,0 +1,87 @@ +# +# Constant-current constant-voltage charge with US06 Drive Cycle using Experiment Class. +# +import pybamm +import pandas as pd +import os + +os.chdir(pybamm.__path__[0] + "/..") + +pybamm.set_logging_level("INFO") + +# import drive cycle from file +drive_cycle_current = pd.read_csv("pybamm/input/drive_cycles/US06.csv", + comment="#", + header=None).to_numpy() + + +# Map Drive Cycle +def Map_Drive_Cycle(x, min_ip_value, max_ip_value, min_op_value, max_op_value): + return (x - min_ip_value) * (max_op_value - min_op_value) / ( + max_ip_value - min_ip_value) + min_op_value + + +# Map Current to Voltage +temp_volts = np.array([]) +min_ip_value = drive_cycle_current[:, 1].min() +max_ip_value = drive_cycle_current[:, 1].max() +min_Voltage = 3.5 +max_Voltage = 4.1 +for I in drive_cycle_current[:, 1]: + temp_volts = np.append( + temp_volts, + Map_Drive_Cycle(I, min_ip_value, max_ip_value, min_Voltage, + max_Voltage)) + +drive_cycle_voltage = drive_cycle_current +drive_cycle_voltage = np.column_stack((np.delete(drive_cycle_voltage, 1, + 1), np.array(temp_volts))) + +# Map Current to Power +temp_volts = np.array([]) +min_ip_value = drive_cycle_current[:, 1].min() +max_ip_value = drive_cycle_current[:, 1].max() +min_Power = 2.5 +max_Power = 5.5 +for I in drive_cycle_current[:, 1]: + temp_volts = np.append( + temp_volts, + Map_Drive_Cycle(I, min_ip_value, max_ip_value, min_Power, max_Power)) + +drive_cycle_power = drive_cycle_current +drive_cycle_power = np.column_stack((np.delete(drive_cycle_power, 1, + 1), np.array(temp_volts))) +experiment = pybamm.Experiment([ +# "Charge at 1 A until 4.1 V", +# "Hold at 4.1 V until 50 mA", +# "Rest for 1 hour", + "Run US06_A (A), +# "Rest for 1 hour", +] +# + [ +# "Charge at 1 A until 4.1 V", +# "Hold at 4.1 V until 50 mA", +# "Rest for 1 hour", +# "Run US06_V (V)", +# "Rest for 1 hour", +# ] + [ +# "Charge at 1 A until 4.1 V", +# "Hold at 4.1 V until 50 mA", +# "Rest for 1 hour", +# "Run US06_W (W)", +# "Rest for 1 hour", +# ] + ,drive_cycles={ + "US06_A": drive_cycle_current, + "US06_V": drive_cycle_voltage, + "US06_W": drive_cycle_power, + }) + +model = pybamm.lithium_ion.DFN() +sim = pybamm.Simulation(model, + experiment=experiment, + solver=pybamm.CasadiSolver()) +sim.solve() + +# Show all plots +sim.plot() diff --git a/pybamm/experiments/experiment.py b/pybamm/experiments/experiment.py index 724d913565..2ef16d8878 100644 --- a/pybamm/experiments/experiment.py +++ b/pybamm/experiments/experiment.py @@ -233,7 +233,7 @@ def read_string(self, cond, drive_cycles): examples ) ) - return electric + (time,) + (period,), events + return electric + (time,) + (period,) + (cond,), events def extend_drive_cycle(self, drive_cycle, end_time): "Extends the drive cycle to enable for event" diff --git a/pybamm/simulation.py b/pybamm/simulation.py index e5c4f0e3a3..ea637dc77f 100644 --- a/pybamm/simulation.py +++ b/pybamm/simulation.py @@ -170,44 +170,82 @@ def set_up_experiment(self, model, experiment): self._experiment_inputs = [] self._experiment_times = [] for op, events in zip(experiment.operating_conditions, experiment.events): - if op[1] in ["A", "C"]: - # Update inputs for constant current + if isinstance(op[0], np.ndarray): + # self.operating_mode = "drive cycle" ########!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + # If ndarray is recived from, create interpolant + # create interpolant + timescale = self._parameter_values.evaluate(model.timescale) + drive_cycle_interpolant = pybamm.Interpolant( + op[0][:, 0], op[0][:, 1], timescale * pybamm.t + ) + # print(type(drive_cycle_interpolant)) #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! if op[1] == "A": - I = op[0] - else: - # Scale C-rate with capacity to obtain current - capacity = self._parameter_values["Nominal cell capacity [A.h]"] - I = op[0] * capacity - operating_inputs = { - "Current switch": 1, - "Voltage switch": 0, - "Power switch": 0, - "Current input [A]": I, - "Voltage input [V]": 0, # doesn't matter - "Power input [W]": 0, # doesn't matter - } - elif op[1] == "V": - # Update inputs for constant voltage - V = op[0] - operating_inputs = { - "Current switch": 0, - "Voltage switch": 1, - "Power switch": 0, - "Current input [A]": 0, # doesn't matter - "Voltage input [V]": V, - "Power input [W]": 0, # doesn't matter - } - elif op[1] == "W": - # Update inputs for constant power - P = op[0] - operating_inputs = { - "Current switch": 0, - "Voltage switch": 0, - "Power switch": 1, - "Current input [A]": 0, # doesn't matter - "Voltage input [V]": 0, # doesn't matter - "Power input [W]": P, - } + operating_inputs = { + "Current switch": 1, + "Voltage switch": 0, + "Power switch": 0, + "Current input [A]": drive_cycle_interpolant, + "Voltage input [V]": 0, # doesn't matter + "Power input [W]": 0, # doesn't matter + } + if op[1] == "V": + operating_inputs = { + "Current switch": 0, + "Voltage switch": 1, + "Power switch": 0, + "Current input [A]": 0, + "Voltage input [V]": drive_cycle_interpolant, # doesn't matter + "Power input [W]": 0, # doesn't matter + } + if op[1] == "W": + operating_inputs = { + "Current switch": 0, + "Voltage switch": 0, + "Power switch": 1, + "Current input [A]": 0, + "Voltage input [V]": 0, # doesn't matter + "Power input [W]": drive_cycle_interpolant, # doesn't matter + } + else: + # self.operating_mode = "experiment" ########!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + if op[1] in ["A", "C"]: + # Update inputs for constant current + if op[1] == "A": + I = op[0] + else: + # Scale C-rate with capacity to obtain current + capacity = self._parameter_values["Nominal cell capacity [A.h]"] + I = op[0] * capacity + operating_inputs = { + "Current switch": 1, + "Voltage switch": 0, + "Power switch": 0, + "Current input [A]": I, + "Voltage input [V]": 0, # doesn't matter + "Power input [W]": 0, # doesn't matter + } + elif op[1] == "V": + # Update inputs for constant voltage + V = op[0] + operating_inputs = { + "Current switch": 0, + "Voltage switch": 1, + "Power switch": 0, + "Current input [A]": 0, # doesn't matter + "Voltage input [V]": V, + "Power input [W]": 0, # doesn't matter + } + elif op[1] == "W": + # Update inputs for constant power + P = op[0] + operating_inputs = { + "Current switch": 0, + "Voltage switch": 0, + "Power switch": 1, + "Current input [A]": 0, # doesn't matter + "Voltage input [V]": 0, # doesn't matter + "Power input [W]": P, + } # Update period operating_inputs["period"] = op[3] # Update events @@ -319,7 +357,7 @@ def set_up_model_for_experiment_old(self, model): self.model = new_model self.op_conds_to_model_and_param = { - op_cond[:2]: (new_model, self.parameter_values) + op_cond[-1]: (new_model, self.parameter_values) for op_cond in set(self.experiment.operating_conditions) } self.op_conds_to_built_models = None @@ -339,7 +377,7 @@ def set_up_model_for_experiment_new(self, model): ): # Create model for this operating condition if it has not already been seen # before - if op_cond[:2] not in self.op_conds_to_model_and_param: + if op_cond[-1] not in self.op_conds_to_model_and_param: if op_inputs["Current switch"] == 1: # Current control # Make a new copy of the model (we will update events later)) @@ -439,8 +477,10 @@ def set_up_model_for_experiment_new(self, model): {"Power function [W]": op_inputs["Power input [W]"]}, check_already_exists=False, ) + + # print(new_parameter_values) #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!) - self.op_conds_to_model_and_param[op_cond[:2]] = ( + self.op_conds_to_model_and_param[op_cond[-1]] = ( new_model, new_parameter_values, ) @@ -691,7 +731,9 @@ def solve( dt = self._experiment_times[idx] op_conds_str = self.experiment.operating_conditions_strings[idx] op_conds_elec = self.experiment.operating_conditions[idx][:2] - model = self.op_conds_to_built_models[op_conds_elec] + # print('Operating Conditions', op_conds_elec) + # model = self.op_conds_to_built_models[op_conds_elec] + model = self.op_conds_to_built_models[op_conds_str] # Use 1-indexing for printing cycle number as it is more # human-intuitive pybamm.logger.notice( diff --git a/pybamm/solvers/base_solver.py b/pybamm/solvers/base_solver.py index 77069c78ce..736d63a53c 100644 --- a/pybamm/solvers/base_solver.py +++ b/pybamm/solvers/base_solver.py @@ -907,7 +907,15 @@ def step( # Set up external variables and inputs external_variables = external_variables or {} inputs = inputs or {} + if isinstance(inputs['Current input [A]'], pybamm.Interpolant): + del inputs['Current input [A]'] + elif isinstance(inputs['Voltage input [V]'], pybamm.Interpolant): + del inputs['Voltage input [V]'] + elif isinstance(inputs['Power input [W]'], pybamm.Interpolant): + del inputs['Power input [W]'] + print("First Inputs",inputs) ext_and_inputs = {**external_variables, **inputs} + print("First EXT and Inputs",ext_and_inputs) # Check that any inputs that may affect the scaling have not changed # Set model timescale @@ -1199,6 +1207,7 @@ def __init__(self, function, model): def __call__(self, inputs): if self.form == "casadi": if isinstance(inputs, dict): + print(inputs) #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! inputs = casadi.vertcat(*[x for x in inputs.values()]) return self._function(0, self.y_dummy, inputs) else: From 351183b2e8bd6fd11da998625d0fc029d8f99c1b Mon Sep 17 00:00:00 2001 From: Ali Hussain Umar Bhatti <65511923+alibh95@users.noreply.github.com> Date: Fri, 25 Jun 2021 15:35:02 +0500 Subject: [PATCH 2/2] Format Code --- pybamm/simulation.py | 18 ++++++------------ pybamm/solvers/base_solver.py | 3 --- 2 files changed, 6 insertions(+), 15 deletions(-) diff --git a/pybamm/simulation.py b/pybamm/simulation.py index ea637dc77f..08a9773049 100644 --- a/pybamm/simulation.py +++ b/pybamm/simulation.py @@ -171,14 +171,12 @@ def set_up_experiment(self, model, experiment): self._experiment_times = [] for op, events in zip(experiment.operating_conditions, experiment.events): if isinstance(op[0], np.ndarray): - # self.operating_mode = "drive cycle" ########!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! # If ndarray is recived from, create interpolant # create interpolant timescale = self._parameter_values.evaluate(model.timescale) drive_cycle_interpolant = pybamm.Interpolant( op[0][:, 0], op[0][:, 1], timescale * pybamm.t ) - # print(type(drive_cycle_interpolant)) #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! if op[1] == "A": operating_inputs = { "Current switch": 1, @@ -193,8 +191,8 @@ def set_up_experiment(self, model, experiment): "Current switch": 0, "Voltage switch": 1, "Power switch": 0, - "Current input [A]": 0, - "Voltage input [V]": drive_cycle_interpolant, # doesn't matter + "Current input [A]": 0, # doesn't matter + "Voltage input [V]": drive_cycle_interpolant, "Power input [W]": 0, # doesn't matter } if op[1] == "W": @@ -202,12 +200,11 @@ def set_up_experiment(self, model, experiment): "Current switch": 0, "Voltage switch": 0, "Power switch": 1, - "Current input [A]": 0, + "Current input [A]": 0, # doesn't matter "Voltage input [V]": 0, # doesn't matter - "Power input [W]": drive_cycle_interpolant, # doesn't matter + "Power input [W]": drive_cycle_interpolant, } else: - # self.operating_mode = "experiment" ########!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! if op[1] in ["A", "C"]: # Update inputs for constant current if op[1] == "A": @@ -460,7 +457,8 @@ def set_up_model_for_experiment_new(self, model): for event in new_model.events if event.name not in ["Minimum voltage", "Maximum voltage"] ] - +# Make Interpolant Here + print("OK GEE") # Update parameter values new_parameter_values = self.parameter_values.copy() if op_inputs["Current switch"] == 1: @@ -731,8 +729,6 @@ def solve( dt = self._experiment_times[idx] op_conds_str = self.experiment.operating_conditions_strings[idx] op_conds_elec = self.experiment.operating_conditions[idx][:2] - # print('Operating Conditions', op_conds_elec) - # model = self.op_conds_to_built_models[op_conds_elec] model = self.op_conds_to_built_models[op_conds_str] # Use 1-indexing for printing cycle number as it is more # human-intuitive @@ -754,9 +750,7 @@ def solve( ) steps.append(step_solution) current_solution = step_solution - cycle_solution = cycle_solution + step_solution - # Only allow events specified by experiment if not ( cycle_solution is None diff --git a/pybamm/solvers/base_solver.py b/pybamm/solvers/base_solver.py index 736d63a53c..34f9f9a2ec 100644 --- a/pybamm/solvers/base_solver.py +++ b/pybamm/solvers/base_solver.py @@ -913,9 +913,7 @@ def step( del inputs['Voltage input [V]'] elif isinstance(inputs['Power input [W]'], pybamm.Interpolant): del inputs['Power input [W]'] - print("First Inputs",inputs) ext_and_inputs = {**external_variables, **inputs} - print("First EXT and Inputs",ext_and_inputs) # Check that any inputs that may affect the scaling have not changed # Set model timescale @@ -1207,7 +1205,6 @@ def __init__(self, function, model): def __call__(self, inputs): if self.form == "casadi": if isinstance(inputs, dict): - print(inputs) #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! inputs = casadi.vertcat(*[x for x in inputs.values()]) return self._function(0, self.y_dummy, inputs) else: