Generalized Multi-Step Charge Rate Optimization

This example demonstrates how to optimize the charge rate in battery experiments using PyBaMM’s new InputParameter functionality while maximizing the negative electrode potential constraint. This approach allows dynamic adjustment of charging currents during optimization.

Key features:

  • Uses PyBaMM InputParameters for dynamic charge current control

  • Optimizes charge rate to maximize charging speed

  • Ensures negative electrode potential stays within safe operating limits

  • Demonstrates advanced experiment design with parameterized protocols

import ionworkspipeline as iwp
from ionworkspipeline.data_fits.models.metrics import Minimum, Time
from ionworkspipeline.data_fits.models.metrics.actions import GreaterThan, Minimize
import matplotlib.pyplot as plt
import numpy as np
import pybamm

Setup Model and Base Parameters

We’ll use a DFN model with the Chen2020 parameter set as our baseline configuration.

# Setup base model and parameters
model = pybamm.lithium_ion.DFN()
params_baseline = pybamm.ParameterValues("Chen2020")

print("Base model and parameter set configured")
print(f"Model: {model.name}")
print("Parameter set: Chen2020")
Base model and parameter set configured
Model: Doyle-Fuller-Newman model
Parameter set: Chen2020

Define Experiment with InputParameters

We’ll create an experiment protocol that uses PyBaMM’s InputParameter functionality to make the charge current optimizable. This allows us to dynamically adjust the charging rate during optimization.

# Number of charge steps with corresponding potential transitions
num_charge_steps = 6
potentials = [3.6, 3.7, 3.9, 4.0, 4.1, 4.2]

# Define optimization parameters for charge current
parameters = {}
for i in range(num_charge_steps):
    parameters[f"I_charge_{i}"] = iwp.Parameter(
        f"I_charge_{i}", initial_value=-4.0, bounds=(-10.0, -0.1)
    )

# Define experiment steps with parameterized charge current
charge_steps = []
i = 0
for _, value in parameters.items():
    charge_steps.append(pybamm.step.current(value, termination=f"> {potentials[i]} V"))
    i += 1

experiment_steps = [
    pybamm.step.current(5.0, termination="< 3.0 V"),  # Discharge at 5A (fixed)
    pybamm.step.rest(duration="10 minutes"),  # Rest period
    *charge_steps,
    pybamm.step.voltage(4.20, termination="< C/20"),  # CV hold
]

# Create experiment with the parameterized steps
experiment = pybamm.Experiment(experiment_steps, period="10 seconds")

print("Experiment protocol with InputParameter defined:")
print(experiment_steps)
Experiment protocol with InputParameter defined:
[Step(5.0, duration=86400, termination=< 3.0 V), Step(0, duration=10 minutes), Step(I_charge_0, duration=86400, termination=> 3.6 V), Step(I_charge_1, duration=86400, termination=> 3.7 V), Step(I_charge_2, duration=86400, termination=> 3.9 V), Step(I_charge_3, duration=86400, termination=> 4.0 V), Step(I_charge_4, duration=86400, termination=> 4.1 V), Step(I_charge_5, duration=86400, termination=> 4.2 V), Step(4.2, duration=86400, termination=< C/20)]

Define Optimization Metrics

We’ll create metrics that focus on:

  1. Minimizing charge time (faster charging)

  2. Maximizing negative electrode potential (safety constraint)

  3. Monitoring overall charging performance

# Create optimization actions
actions = {
    # Minimize charge time (maximize charge rate)
    "Final time [s]": Minimize(Time(variable="Time [s]", value=-1), weight=0.01),
}

# Constraints define bounds that must be satisfied (using GreaterThan/LessThan actions)
constraints = {
    # Penalize negative electrode potential below 0 V (safety constraint)
    "Anode potential constraint violation [V]": GreaterThan(
        Minimum(
            "Negative electrode surface potential difference at separator interface [V]"
        ),
        value=0.0,
    ),
}

print(
    f"Defined {len(actions)} optimization metrics and {len(constraints)} constraints:"
)
print("  1. Charge Time (minimize) - Faster charging")
print("  2. Min Negative Electrode Potential constraint - Safety constraint")
Defined 1 optimization metrics and 1 constraints:
  1. Charge Time (minimize) - Faster charging
  2. Min Negative Electrode Potential constraint - Safety constraint

Set up Custom Cost Function

We’ll create a weighted cost function that balances fast charging with safety constraints.

print("Custom cost function created with variable weights")
Custom cost function created with variable weights

Run Charge Rate Optimization

Now we’ll run the optimization to find the optimal charge current that maximizes charging speed while maintaining safe electrode potentials.

# Set up design objective with InputParameter support
objective_options = {
    "model": model,
    "simulation_kwargs": {
        "experiment": experiment,
    },
}

objective = iwp.objectives.DesignObjective(
    options=objective_options,
    actions=actions,
    constraints=constraints,
)

# Create and run optimization
optim = iwp.optimizers.Pints(
    max_iterations=10,  # Increase for better results (~200)
    max_unchanged_iterations=80,
    max_unchanged_iterations_threshold=1e-4,
    log_to_screen=True,
    sigma0=0.2,
    # population_size=30, # Enable during parallel optimisation
)
design_datafit = iwp.DataFit(
    objective,
    parameters=parameters,
    options={"seed": 42},
    optimizer=optim,
    # parallel=True,
)

print("Running charge rate optimization...")
print("This may take a few minutes as we optimize for safe fast charging...")

results = design_datafit.run(params_baseline)

print("\nCharge rate optimization completed!")
print(f"Final cost: {results.costs:.6f}")
print("\nOptimal parameters:")
for name in parameters:
    value = results.parameter_values[name]
    print(f"  Optimal charge current: {value:.3f}A")

optimal_charge_current = results.parameter_values
Running charge rate optimization...
This may take a few minutes as we optimize for safe fast charging...
Minimising error measure
Using CMAES
Population size: 9
--------------------------------------------------------------------------------
 Iter.  Eval.            Best         Current     Time m:s
--------------------------------------------------------------------------------
     1      9    1.000186e+06    1.000186e+06     0:02
     2     18    9.955188e+01    9.955188e+01     0:04
     3     27    9.955188e+01    9.955188e+01     0:05
     4     36    9.955188e+01    9.955188e+01     0:07
     5     45    9.955188e+01    9.955188e+01     0:09
    10     90    9.870133e+01    9.870133e+01     0:18
--------------------------------------------------------------------------------
Maximum number of iterations (10) reached.
Optimization complete

Charge rate optimization completed!
Final cost: 98.701328

Optimal parameters:
  Optimal charge current: -4.502A
  Optimal charge current: -4.302A
  Optimal charge current: -5.421A
  Optimal charge current: -5.519A
  Optimal charge current: -4.570A
  Optimal charge current: -3.040A

Validate and Compare Results

Let’s compare the optimized charge rate with a baseline charge rate to see the improvements.

# Helper function to extract charging metrics
def extract_charge_metrics(solution, charge_step_index=2):
    """Extract key charging performance metrics from simulation solution."""
    charge_sol = solution.sub_solutions[charge_step_index:]

    # Extract charging metrics
    charge_time = (charge_sol[-1].t[-1] - charge_sol[0].t[0]) / 3600  # Convert to hours

    # Negative electrode potential analysis
    anode_potential = []
    for sol in charge_sol:
        anode_potential.extend(
            sol[
                "Negative electrode surface potential difference at separator interface [V]"
            ]()
        )
    min_anode_potential = min(anode_potential)
    mean_anode_potential = np.asarray(anode_potential).mean()

    # Voltage analysis
    voltage = charge_sol[-1]["Voltage [V]"]()
    final_voltage = voltage[-1]

    return {
        "charge_time": charge_time,
        "min_anode_potential": min_anode_potential,
        "mean_anode_potential": mean_anode_potential,
        "final_voltage": final_voltage,
        "solution": solution,
        "charge_sol": charge_sol,
    }


# Run simulations for comparison
input_names = list(parameters.keys())
baseline_current = -4.0  # Baseline charge current
baseline_inputs = {name: baseline_current for name in input_names}
optimal_inputs = {}
for name in input_names:
    optimal_inputs[name] = results.parameter_values[name]

print("Running comparison simulations...")
print(f"Baseline charge current: {baseline_current}A")

# Baseline simulation
sim_baseline = iwp.Simulation(
    model=model, parameter_values=params_baseline, experiment=experiment
)
sol_baseline = sim_baseline.solve(
    initial_soc=1.0,
    inputs=baseline_inputs,
)

# Optimal simulation
sim_optimal = iwp.Simulation(
    model=model, parameter_values=params_baseline, experiment=experiment
)
sol_optimal = sim_optimal.solve(
    initial_soc=1.0,
    inputs=optimal_inputs,
)

# Extract metrics
baseline_metrics = extract_charge_metrics(sol_baseline)
optimal_metrics = extract_charge_metrics(sol_optimal)

print("\nCharging Performance Comparison:")
print("=" * 75)
print(f"{'Metric':<35} {'Baseline':<15} {'Optimal':<15} {'Improvement'}")
print("-" * 75)

# Calculate improvements
time_improvement = (
    (baseline_metrics["charge_time"] - optimal_metrics["charge_time"])
    / baseline_metrics["charge_time"]
    * 100
)
min_anode_change = (
    optimal_metrics["min_anode_potential"] - baseline_metrics["min_anode_potential"]
)

print(
    f"{'Charge Time [h]':<35} {baseline_metrics['charge_time']:<15.4f} {optimal_metrics['charge_time']:<15.4f} {time_improvement:+.1f}%"
)
print(
    f"{'Min Anode Potential [V]':<35} {baseline_metrics['min_anode_potential']:<15.4f} {optimal_metrics['min_anode_potential']:<15.4f} {min_anode_change:+.4f}"
)
print(
    f"{'Mean Anode Potential [V]':<35} {baseline_metrics['mean_anode_potential']:<15.4f} {optimal_metrics['mean_anode_potential']:<15.4f} {(optimal_metrics['mean_anode_potential'] - baseline_metrics['mean_anode_potential']):+.4f}"
)
print("=" * 75)
Running comparison simulations...
Baseline charge current: -4.0A
Charging Performance Comparison:
===========================================================================
Metric                              Baseline        Optimal         Improvement
---------------------------------------------------------------------------
Charge Time [h]                     1.7207          1.6586          +3.6%
Min Anode Potential [V]             -0.0075         0.0024          +0.0099
Mean Anode Potential [V]            0.0474          0.0424          -0.0050
===========================================================================

Parameter Analysis and Visualization

Let’s analyze the optimized parameters and create comprehensive visualizations of the results.

# Extract parameter comparison data
print("Parameter Optimization Results:")
print("=" * 60)
for name, param in parameters.items():
    baseline_val = param.initial_value
    optimal_val = results.parameter_values[name]
    change = optimal_val - baseline_val
    percent_change = (change / baseline_val) * 100

    short_name = name.replace("I_charge_", "Charge Current ")
    print(f"\n{short_name}:")
    print(f"  Baseline:       {baseline_val:.4f}A")
    print(f"  Optimal:        {optimal_val:.4f}A")
    print(f"  Change:         {change:+.4f}A ({percent_change:+.2f}%)")

print("\n" + "=" * 60)

# Create comprehensive visualization
fig, axes = plt.subplots(2, 2, figsize=(12.5, 10))

# Plot 1: Voltage comparison
t_baseline = sol_baseline.t / 3600
t_optimal = sol_optimal.t / 3600

axes[0, 0].plot(
    t_baseline,
    sol_baseline["Voltage [V]"](sol_baseline.t),
    label="Baseline",
    linewidth=2,
    alpha=0.8,
)
axes[0, 0].plot(
    t_optimal,
    sol_optimal["Voltage [V]"](sol_optimal.t),
    label="Optimal",
    linewidth=2,
    alpha=0.8,
    linestyle="--",
)
axes[0, 0].set_xlabel("Time [h]")
axes[0, 0].set_ylabel("Voltage [V]")
axes[0, 0].set_title("Cell Voltage: Baseline vs Optimal")
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3)

# Plot 2: Current comparison
axes[0, 1].plot(
    t_baseline,
    sol_baseline["Current [A]"](sol_baseline.t),
    label="Baseline",
    linewidth=2,
    alpha=0.8,
)
axes[0, 1].plot(
    t_optimal,
    sol_optimal["Current [A]"](sol_optimal.t),
    label="Optimal",
    linewidth=2,
    alpha=0.8,
    linestyle="--",
)
axes[0, 1].set_xlabel("Time [h]")
axes[0, 1].set_ylabel("Current [A]")
axes[0, 1].set_title("Current Profile: Baseline vs Optimal")
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3)

# Plot 3: Anode potential during charge
# Combine charge steps for both solutions
charge_steps_baseline = sol_baseline.sub_solutions[2:-1]
charge_steps_optimal = sol_optimal.sub_solutions[2:-1]

# For baseline
charge_sol_baseline = charge_steps_baseline[0]
for step in charge_steps_baseline[1:]:
    charge_sol_baseline = charge_sol_baseline + step

# For optimal
charge_sol_optimal = charge_steps_optimal[0]
for step in charge_steps_optimal[1:]:
    charge_sol_optimal = charge_sol_optimal + step

t_charge_baseline = charge_sol_baseline.t / 3600 - charge_sol_baseline.t[0] / 3600
t_charge_optimal = charge_sol_optimal.t / 3600 - charge_sol_optimal.t[0] / 3600

anode_baseline = charge_sol_baseline[
    "Negative electrode surface potential difference at separator interface [V]"
](charge_sol_baseline.t)
anode_optimal = charge_sol_optimal[
    "Negative electrode surface potential difference at separator interface [V]"
](charge_sol_optimal.t)

axes[1, 0].plot(
    t_charge_baseline, anode_baseline, label="Baseline", linewidth=2, alpha=0.8
)
axes[1, 0].plot(
    t_charge_optimal,
    anode_optimal,
    label="Optimal",
    linewidth=2,
    alpha=0.8,
    linestyle="--",
)
axes[1, 0].axhline(y=0.0, color="r", linestyle=":", alpha=0.7, label="Safety limit")
axes[1, 0].set_xlabel("Charge Time [h]")
axes[1, 0].set_ylabel("Electrode - Electrolyte Potential [V]")
axes[1, 0].set_title("Negative Electrode - Electrolyte Potential During Charge")
axes[1, 0].legend()
axes[1, 0].grid(True, alpha=0.3)

# Plot 4: Performance metrics comparison
metrics_names = [
    "Charge Time [h]",
    "Min Anode\nPotential [V]",
    "Mean Anode\nPotential [V]",
]
baseline_values = [
    baseline_metrics["charge_time"],
    baseline_metrics["min_anode_potential"],
    baseline_metrics["mean_anode_potential"],
]
optimal_values = [
    optimal_metrics["charge_time"],
    optimal_metrics["min_anode_potential"],
    optimal_metrics["mean_anode_potential"],
]

x = np.arange(len(metrics_names))
width = 0.35

bars1 = axes[1, 1].bar(
    x - width / 2,
    baseline_values,
    width,
    label="Baseline",
    alpha=0.8,
    color="lightblue",
)
bars2 = axes[1, 1].bar(
    x + width / 2, optimal_values, width, label="Optimal", alpha=0.8, color="orange"
)

# Add value labels on bars
for i, (baseline, optimal) in enumerate(
    zip(baseline_values, optimal_values, strict=False)
):
    y_offset = abs(baseline) * 0.02 if baseline != 0 else 0.005
    axes[1, 1].text(
        i - width / 2,
        baseline + y_offset,
        f"{baseline:.3f}",
        ha="center",
        va="bottom",
        fontsize=9,
    )
    axes[1, 1].text(
        i + width / 2,
        optimal + y_offset,
        f"{optimal:.3f}",
        ha="center",
        va="bottom",
        fontsize=9,
    )

axes[1, 1].set_xlabel("Performance Metrics")
axes[1, 1].set_ylabel("Value")
axes[1, 1].set_ylim(top=3)
axes[1, 1].set_title("Performance Metrics: Baseline vs Optimal")
axes[1, 1].set_xticks(x)
axes[1, 1].set_xticklabels(metrics_names)
axes[1, 1].legend()
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

print("\nCharge rate optimization analysis completed!")
Parameter Optimization Results:
============================================================

Charge Current 0:
  Baseline:       -4.0000A
  Optimal:        -4.5017A
  Change:         -0.5017A (+12.54%)

Charge Current 1:
  Baseline:       -4.0000A
  Optimal:        -4.3020A
  Change:         -0.3020A (+7.55%)

Charge Current 2:
  Baseline:       -4.0000A
  Optimal:        -5.4206A
  Change:         -1.4206A (+35.52%)

Charge Current 3:
  Baseline:       -4.0000A
  Optimal:        -5.5191A
  Change:         -1.5191A (+37.98%)

Charge Current 4:
  Baseline:       -4.0000A
  Optimal:        -4.5696A
  Change:         -0.5696A (+14.24%)

Charge Current 5:
  Baseline:       -4.0000A
  Optimal:        -3.0401A
  Change:         +0.9599A (-24.00%)

============================================================
../../../_images/e6ce1d7d71a975c005f7f4f40f64ce2284134065946e01e7f909d0a2e6c3ec3a.png
Charge rate optimization analysis completed!