Skip to content

simulation/

API reference for the simulation package.

simulation.simulation

Core physics model and energy simulation.

simulation.simulation

Simulation for MSXVI.

Date: 2025-11-07 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

VehicleParams dataclass

A data class to store parameters for the model.

Source code in simulation/simulation.py 
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
@dataclass
class VehicleParams:
    """A data class to store parameters for the model."""

    # fmt: off
    mass: float = 450.0                             # mass
    drag_coeff: float = 0.18                        # drag coefficient
    front_area: float = 1.357                       # frontal area
    C_RR: float = 0.004                             # rolling resistance
    solar_area: float = 4.0                         # Solar array area
    panel_eff: float = 0.243                        # electrical panel efficiency (fraction)
    bat_max_energy: float = 40 * 3.63 * 36 * 3600   # battery energy capacity, Joules (40 - Ah, 3.63 - nominal cell voltage, 36 - number of cells, 3600 - watt hours to joules)
    air_density: float = 1.293                      # air density
    gravity_const: float = 9.81                     # gravity
    drive_eff: float = 0.94                         # drivetrain efficiency on traction losses in the set [0, 1]

SimResult

Bases: NamedTuple

The expected output values.

Source code in simulation/simulation.py 
33
34
35
36
37
38
class SimResult(NamedTuple):
    """The expected output values."""

    final_distance_m: float
    final_soc_J: float
    traces: Dict[str, np.ndarray]

rolling_power(v, M, g, C_RR)

Calculates rolling resistance power loss.

Source code in simulation/simulation.py 
44
45
46
def rolling_power(v: np.ndarray, M: float, g: float, C_RR: float) -> np.ndarray:
    """Calculates rolling resistance power loss."""
    return (M * g * C_RR) * v

drag_power(v, rho, Cd, A)

Calculates aerodynamic drag power loss.

Source code in simulation/simulation.py 
49
50
51
def drag_power(v: np.ndarray, rho: float, Cd: float, A: float) -> np.ndarray:
    """Calculates aerodynamic drag power loss."""
    return 0.5 * rho * Cd * A * v**3

grade_power(v, theta_rad, M, g)

Calculates gravitational power (positive uphill, negative downhill).

Source code in simulation/simulation.py 
54
55
56
def grade_power(v: np.ndarray, theta_rad: np.ndarray, M: float, g: float) -> np.ndarray:
    """Calculates gravitational power (positive uphill, negative downhill)."""
    return M * g * np.sin(theta_rad) * v

solar_power(G_wm2, A_solar, panel_eff)

Calculates solar power generation.

Source code in simulation/simulation.py 
59
60
61
def solar_power(G_wm2: np.ndarray, A_solar: float, panel_eff: float) -> np.ndarray:
    """Calculates solar power generation."""
    return A_solar * panel_eff * G_wm2

simulate(v, dt, d0, theta_fn, ghi_fn, params=VehicleParams())

Vectorized simulation over a fixed horizon.

Parameters:

Name Type Description Default
v ndarray

Speed profile (m/s).

 required
dt float

Timestep (s).

 required
d0 float

Starting distance (m).

 required
theta_fn Callable[[ndarray], ndarray]

Callable that returns road grade (deg) at given distance(s).

 required
ghi_fn Callable[[ndarray], ndarray]

Callable that returns GHI (W/m^2) at given distance(s).

 required
params VehicleParams

Vehicle parameters (see dataclass).

 VehicleParams()

Returns:

Type Description
SimResult

SimResult(final_distance_m, final_soc_J, traces).
Source code in simulation/simulation.py 
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
def simulate(
    v: np.ndarray,
    dt: float,
    d0: float,
    theta_fn: Callable[[np.ndarray], np.ndarray],
    ghi_fn: Callable[[np.ndarray], np.ndarray],
    params: VehicleParams = VehicleParams(),
) -> SimResult:
    """Vectorized simulation over a fixed horizon.

    Args:
        v: Speed profile (m/s).
        dt: Timestep (s).
        d0: Starting distance (m).
        theta_fn: Callable that returns road grade (deg) at given distance(s).
        ghi_fn: Callable that returns GHI (W/m^2) at given distance(s).
        params: Vehicle parameters (see dataclass).

    Returns:
        SimResult(final_distance_m, final_soc_J, traces).
    """
    v = np.asarray(v, dtype=float)
    N = v.size
    p = params

    # Distance timeline (Euler)
    d = np.empty(N + 1, dtype=float)
    d[0] = d0
    d[1:] = d0 + np.cumsum(v * dt)

    # Query terrain and irradiance at start-of-step positions
    theta_deg = theta_fn(d[:-1])
    ghi = ghi_fn(d[:-1])

    theta_rad = np.deg2rad(theta_deg)

    # Drivetrain power
    P_rr = rolling_power(v, p.mass, p.gravity_const, p.C_RR)
    P_drag = drag_power(v, p.air_density, p.drag_coeff, p.front_area)
    P_grade_raw = grade_power(v, theta_rad, p.mass, p.gravity_const)
    # Apply drivetrain efficiency
    # TODO Account for Regen
    P_grade_loss_pos = np.maximum(P_grade_raw, 0.0) / max(p.drive_eff, 1e-9)
    P_rr_drive = P_rr / max(p.drive_eff, 1e-9)
    P_drag_drive = P_drag / max(p.drive_eff, 1e-9)

    # Solar Power
    P_solar = solar_power(ghi, p.solar_area, p.panel_eff)

    # Net battery power draw
    # (- Means drawing from battery)
    # (Losses (positive) minus solar)
    P_net = (P_rr_drive + P_drag_drive + P_grade_loss_pos) - P_solar

    E_rr = P_rr * dt
    E_drag = P_drag * dt
    E_grade = P_grade_raw * dt
    E_solar = P_solar * dt
    E_net = P_net * dt

    Ebat_raw = np.empty(N + 1, dtype=float)
    Ebat_raw[0] = p.bat_max_energy
    Ebat_raw[1:] = p.bat_max_energy - np.cumsum(E_net)

    # For display only
    Ebat = np.clip(Ebat_raw, 0.0, p.bat_max_energy)

    # debugging purposes in case numbers look funky
    traces = {
        "distance_m": d,
        "Ebat_J": Ebat,
        "Ebat_raw_J": Ebat_raw,
        "P_rr_W": P_rr,
        "P_drag_W": P_drag,
        "P_grade_W": P_grade_raw,
        "P_solar_W": P_solar,
        "P_net_W": P_net,
        "E_rr_J": E_rr,
        "E_drag_J": E_drag,
        "E_grade_J": E_grade,
        "E_solar_J": E_solar,
        "E_net_J": E_net,
        "theta_deg": theta_deg,
        "ghi": ghi,
        "v": v,
    }
    return SimResult(
        final_distance_m=float(d[-1]), final_soc_J=float(Ebat_raw[-1]), traces=traces
    )

wh_from_joules(E_J)

Convert energy from Joules to Watt-hours.

Parameters:

Name Type Description Default
E_J ndarray | float

Energy in Joules.

 required

Returns:

Type Description
ndarray | float

Energy in Watt-hours.
Source code in simulation/simulation.py 
161
162
163
164
165
166
167
168
169
170
def wh_from_joules(E_J: np.ndarray | float) -> np.ndarray | float:
    """Convert energy from Joules to Watt-hours.

    Args:
        E_J: Energy in Joules.

    Returns:
        Energy in Watt-hours.
    """
    return np.asarray(E_J) / 3600.0

simulation.optimizer

Velocity profile optimization.

simulation.optimizer

Optimizer

Date: 2025-11-09 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

OptimizeConfig dataclass

A configuration dataclass for the optimizer.

Source code in simulation/optimizer.py 
17
18
19
20
21
22
23
24
25
26
27
28
29
@dataclass
class OptimizeConfig:
    """A configuration dataclass for the optimizer."""

    # fmt: off
    dt: float = 10.0            # Time step (s)
    horizon: int = 28800        # How many seconds we want to simulate (s)
    d0: float = 0.0             # Starting distance (m)
    vmin: float = 8.9           # Minimum allowed speed (m/s), default value as per ASC 2026 regs
    vmax: float = 29.0          # Maximum allowed speed (m/s), default value as per ASC 2026 regs
    method: str = "SLSQP"       # Optimization method from scipy
    max_iter: int = 1000        # Maximum itterations (regardless of if we reach convergence or not)
    min_soc: float = 0.2        # Minimum SOC threshold at the end of the simulation

objective(vs, dt, d0, theta_fn, ghi_fn, params, cfg)

Objective for SciPy minimize.

Parameters:

Name Type Description Default
vs ndarray

The speed profile (m/s).

 required
theta_fn Callable

Callable returning road grade (deg) at given distance(s).

 required
ghi_fn Callable

Callable returning GHI (W/m^2) at given distance(s).

 required
params VehicleParams

For simulation.

 required

Returns:

Type Description
float

Negative distance (so minimizing -> maximize distance).
Source code in simulation/optimizer.py 
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
def objective(
    vs: np.ndarray,
    dt: float,
    d0: float,
    theta_fn: Callable,
    ghi_fn: Callable,
    params: VehicleParams,
    cfg: OptimizeConfig,
) -> float:
    """Objective for SciPy minimize.

    Args:
        vs: The speed profile (m/s).
        theta_fn: Callable returning road grade (deg) at given distance(s).
        ghi_fn: Callable returning GHI (W/m^2) at given distance(s).
        params: For simulation.

    Returns:
        Negative distance (so minimizing -> maximize distance).
    """
    res = simulate(vs, dt, d0, theta_fn, ghi_fn, params)
    min_reserve = cfg.min_soc * params.bat_max_energy
    margin = np.min(res.traces["Ebat_raw_J"][1:] - min_reserve)

    if margin < 0:
        return 1e12 + 1e6 * (-margin)

    return -res.final_distance_m

SLSQP_velocity(cfg, theta_fn, ghi_fn, params)

Runs an SLSQP optimization on given slope and irradiance callables.

Parameters:

Name Type Description Default
cfg OptimizeConfig

Optimizer configuration data (see OptimizeConfig dataclass).

 required
theta_fn Callable

Callable returning road grade (deg) at given distance(s).

 required
ghi_fn Callable

Callable returning GHI (W/m^2) at given distance(s).

 required

Returns:

Type Description
Tuple[ndarray, float]

Best velocity array and objective value.
Source code in simulation/optimizer.py 
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
def SLSQP_velocity(
    cfg: OptimizeConfig, theta_fn: Callable, ghi_fn: Callable, params: VehicleParams
) -> Tuple[np.ndarray, float]:
    """Runs an SLSQP optimization on given slope and irradiance callables.

    Args:
        cfg: Optimizer configuration data (see OptimizeConfig dataclass).
        theta_fn: Callable returning road grade (deg) at given distance(s).
        ghi_fn: Callable returning GHI (W/m^2) at given distance(s).

    Returns:
        Best velocity array and objective value.
    """
    N = int(cfg.horizon / cfg.dt)

    vs0 = np.full(N, cfg.vmin)

    bounds = [(cfg.vmin, cfg.vmax)] * N

    # Per timestep SOC >= cfg.min_soc constraints
    min_reserve = cfg.min_soc * params.bat_max_energy  # In Joules

    def soc_constraints(vs):
        """Return SOC - reserve for every timestep; each must be >= 0."""
        res = simulate(vs, cfg.dt, cfg.d0, theta_fn, ghi_fn, params)
        return res.traces["Ebat_raw_J"][1:] - min_reserve

    constraints = [{"type": "ineq", "fun": soc_constraints}]

    iteration_log = []

    def callback(xk):
        """Called by scipy's optimize function. Shows traces of each itteration"""
        res = simulate(xk, cfg.dt, cfg.d0, theta_fn, ghi_fn, params)
        obj = -res.final_distance_m
        soc_pct = res.final_soc_J / params.bat_max_energy * 100
        iteration_log.append(
            {
                "iter": len(iteration_log),
                "objective": obj,
                "distance_km": res.final_distance_m / 1000,
                "final_soc_pct": soc_pct,
                "mean_speed": float(np.mean(xk)),
            }
        )
        entry = iteration_log[-1]
        print(
            f"  Iter {entry['iter']:4d} | "
            f"dist={entry['distance_km']:.2f} km | "
            f"SOC={entry['final_soc_pct']:.1f}% | "
            f"avg_v={entry['mean_speed']:.2f} m/s"
        )

    # Pass cfg through to the objective so it can access tunable weights.
    result = minimize(
        objective,
        vs0,
        args=(cfg.dt, cfg.d0, theta_fn, ghi_fn, params, cfg),
        method=cfg.method,
        bounds=bounds,
        constraints=constraints,
        callback=callback,
        options={"maxiter": cfg.max_iter, "disp": True},
    )
    if not result.success:
        print(f"Warning: Optimization did not converge. Message: {result.message}")

    best_vs = result.x
    best_obj = result.fun
    return best_vs, best_obj

exhaustive_search_velocity(cfg, theta_fn, ghi_fn, params)

Perform an exhaustive search (try every single velocity vector) to find optimized velocity

Mainly used for benchmarking, hence small N

Parameters:

Name Type Description Default
cfg OptimizeConfig

Optimizer configuration data (see OptimizeConfig dataclass).

 required
theta_fn Callable

Callable that takes an array of indices or distances and returns road grade (deg) array.

 required
ghi_fn Callable

Callable that takes an array of indices or distances and returns GHI (W/m^2) array.

 required
params VehicleParams

Vehicle parameters.

 required

Returns:

Type Description
Tuple[ndarray, float]

Best velocity array and objective value.
Source code in simulation/optimizer.py 
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
def exhaustive_search_velocity(
    cfg: OptimizeConfig, theta_fn: Callable, ghi_fn: Callable, params: VehicleParams
) -> Tuple[np.ndarray, float]:
    """Perform an exhaustive search (try every single velocity vector) to find optimized velocity

    Mainly used for benchmarking, hence small N

    Args:
        cfg: Optimizer configuration data (see OptimizeConfig dataclass).
        theta_fn: Callable that takes an array of indices or distances and returns road grade (deg) array.
        ghi_fn: Callable that takes an array of indices or distances and returns GHI (W/m^2) array.
        params: Vehicle parameters.

    Returns:
        Best velocity array and objective value.
    """

    # N = int(cfg.horizon / cfg.dt)
    N = 10  # For testing, set to 10. Change back to above for full search (warning: very slow!)

    v_grid = np.linspace(cfg.vmin, cfg.vmax, num=4)

    best_vs = None
    best_obj = float("inf")
    min_reserve = cfg.min_soc * params.bat_max_energy

    print(
        f"Starting exhaustive search for N={N}. Total permutations: {len(v_grid) ** N}"
    )

    for vs_tuple in itertools.product(v_grid, repeat=N):
        vs = np.array(vs_tuple)
        res = simulate(vs, cfg.dt, cfg.d0, theta_fn, ghi_fn, params)
        margin = np.min(res.traces["Ebat_raw_J"][1:] - min_reserve)
        if margin < 0:
            continue  # Battery below SOC, skip this profile
        obj = -res.final_distance_m
        if obj < best_obj:
            best_obj = obj
            best_vs = vs.copy()
            formatted_vs = "[" + ", ".join(f"{float(v):.2f}" for v in vs) + "]"
            print(
                f"New best: dist={res.final_distance_m / 1000:.2f} km | "
                f"vs={formatted_vs}"
            )

    if best_vs is None:
        best_vs = np.full(N, cfg.vmin)
        best_obj = 1e12

    return best_vs, best_obj

simulation.config

Runtime configuration management.

simulation.config

Configuration management for strategy simulations.

Date: 2025-11-25 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

SimConfig dataclass

Configuration parameters for simulation and optimization.

Source code in simulation/config.py 
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
@dataclass
class SimConfig:
    """Configuration parameters for simulation and optimization."""

    # fmt: off
    # Optimization parameters
    dt: float = 1800.0          # Timestep
    horizon: float = 28800.0    # Driving time per day (s), default 8 hours
    vmin: float = 8.9           # Minimum speed
    vmax: float = 30.0          # Maximum speed
    method: str = "SLSQP"       # Optimization method
    max_iter: int = 2000        # Maximum iterations
    min_soc: float = 0.2        # Energy penalty weight

    # Data sources
    use_solcast: bool = False              # Use Solcast API for GHI data
    solcast_api_key: Optional[str] = None  # Solcast API key
    gpx_file: str = "0_FullBaseRoute.gpx"  # GPX filename
    # fmt: on

    def display(self) -> None:
        """Display all configuration parameters with indices."""
        params = [
            ("dt", self.dt, "s", "Timestep"),
            ("horizon", self.horizon, "s", "Driving time per day"),
            ("vmin", self.vmin, "m/s", "Minimum speed"),
            ("vmax", self.vmax, "m/s", "Maximum speed"),
            ("method", self.method, "", "Optimization method"),
            ("max_iter", self.max_iter, "", "Max iterations"),
            ("min_soc", self.min_soc, "", "Energy penalty weight"),
            ("use_solcast", self.use_solcast, "", "Use Solcast API"),
            ("gpx_file", self.gpx_file, "", "GPX filename"),
        ]

        print(f"\n{'-' * 50}")
        print("Current Configuration")
        print(f"{'-' * 50}")
        for idx, (name, value, unit, description) in enumerate(params, 1):
            unit_str = f" {unit}" if unit else ""
            print(f"{idx:2d}. {description:25s} = {value}{unit_str}")
        print(f"{'-' * 50}\n")

    def update_param(self, param_name: str, value: any) -> bool:
        """Update a configuration parameter.

        Args:
            param_name: Name of the parameter to update.
            value: New value for the parameter.

        Returns:
            True if successful, False otherwise.
        """
        if hasattr(self, param_name):
            # Type conversion based on current type
            current_val = getattr(self, param_name)
            try:
                if isinstance(current_val, bool):
                    new_val = (
                        value.lower() in ("true", "1", "yes", "y")
                        if isinstance(value, str)
                        else bool(value)
                    )
                elif isinstance(current_val, int):
                    new_val = int(value)
                elif isinstance(current_val, float):
                    new_val = float(value)
                else:
                    new_val = value
                setattr(self, param_name, new_val)
                return True
            except (ValueError, AttributeError):
                return False
        return False

display()

Display all configuration parameters with indices.

Source code in simulation/config.py 
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
def display(self) -> None:
    """Display all configuration parameters with indices."""
    params = [
        ("dt", self.dt, "s", "Timestep"),
        ("horizon", self.horizon, "s", "Driving time per day"),
        ("vmin", self.vmin, "m/s", "Minimum speed"),
        ("vmax", self.vmax, "m/s", "Maximum speed"),
        ("method", self.method, "", "Optimization method"),
        ("max_iter", self.max_iter, "", "Max iterations"),
        ("min_soc", self.min_soc, "", "Energy penalty weight"),
        ("use_solcast", self.use_solcast, "", "Use Solcast API"),
        ("gpx_file", self.gpx_file, "", "GPX filename"),
    ]

    print(f"\n{'-' * 50}")
    print("Current Configuration")
    print(f"{'-' * 50}")
    for idx, (name, value, unit, description) in enumerate(params, 1):
        unit_str = f" {unit}" if unit else ""
        print(f"{idx:2d}. {description:25s} = {value}{unit_str}")
    print(f"{'-' * 50}\n")

update_param(param_name, value)

Update a configuration parameter.

Parameters:

Name Type Description Default
param_name str

Name of the parameter to update.

 required
value any

New value for the parameter.

 required

Returns:

Type Description
bool

True if successful, False otherwise.
Source code in simulation/config.py 
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
def update_param(self, param_name: str, value: any) -> bool:
    """Update a configuration parameter.

    Args:
        param_name: Name of the parameter to update.
        value: New value for the parameter.

    Returns:
        True if successful, False otherwise.
    """
    if hasattr(self, param_name):
        # Type conversion based on current type
        current_val = getattr(self, param_name)
        try:
            if isinstance(current_val, bool):
                new_val = (
                    value.lower() in ("true", "1", "yes", "y")
                    if isinstance(value, str)
                    else bool(value)
                )
            elif isinstance(current_val, int):
                new_val = int(value)
            elif isinstance(current_val, float):
                new_val = float(value)
            else:
                new_val = value
            setattr(self, param_name, new_val)
            return True
        except (ValueError, AttributeError):
            return False
    return False

simulation.scenarios

Test and race day scenario runners.

simulation.scenarios

Scenario runners for test and race day simulations.

Date: 2025-11-25 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

run_test_scenario(test_path, config)

Run a test scenario using mock YAML or CSV data.

Parameters:

Name Type Description Default
test_path str

Path to test file (YAML or CSV).

 required
config SimConfig

Simulation configuration.

 required

Returns:

Type Description
SimResult

Simulation result.
Source code in simulation/scenarios.py 
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
def run_test_scenario(test_path: str, config: SimConfig) -> SimResult:
    """Run a test scenario using mock YAML or CSV data.

    Args:
        test_path: Path to test file (YAML or CSV).
        config: Simulation configuration.

    Returns:
        Simulation result.
    """
    print(f"\nLoading test scenario: {test_path}")
    if test_path.lower().endswith(".csv"):
        theta_deg_arr, ghi_arr, meta = load_mock_csv(test_path)
    else:
        theta_deg_arr, ghi_arr, meta = load_mock_yaml(test_path)

    dt = meta.get("dt", config.dt)
    d0 = meta.get("d0", 0.0)
    N_steps = theta_deg_arr.size
    horizon = N_steps * dt

    # Build position-based lookup functions from test arrays
    avg_v = (config.vmin + config.vmax) / 2
    total_dist = N_steps * dt * avg_v
    dist_points = np.linspace(0, total_dist, N_steps)

    def theta_fn(d):
        return np.interp(d, dist_points, theta_deg_arr)

    def ghi_fn(d):
        return np.interp(d, dist_points, ghi_arr)

    params = VehicleParams()

    print("\nOptimization Methods:")
    print("  1. SLSQP (fast, gradient-based)")
    print("  2. Exhaustive Search (slow, brute-force)")
    choice = input("\nSelect optimization method: ").strip()

    if choice == "2":
        method = "exhaustive"
    else:
        method = "SLSQP"

    opt_cfg = OptimizeConfig(
        dt=dt,
        horizon=horizon,
        d0=d0,
        vmin=config.vmin,
        vmax=config.vmax,
        method=method,
        max_iter=config.max_iter,
        min_soc=config.min_soc,
    )

    if method == "SLSQP":
        print("Optimizing velocity profile using SLSQP...")
        best_vs, best_obj = SLSQP_velocity(opt_cfg, theta_fn, ghi_fn, params)
    elif method == "exhaustive":
        print("Optimizing velocity profile using Exhaustive Search...")
        best_vs, best_obj = exhaustive_search_velocity(
            opt_cfg, theta_fn, ghi_fn, params
        )

    print(f"Optimization complete. Objective: {best_obj:.2f}")
    print("Simulating...")
    res = simulate(best_vs, dt, d0, theta_fn, ghi_fn, params)

    _print_results(res, best_vs, params)

    # Ask if user wants to see plots
    plot_choice = input("\nGenerate plots? (y/n): ").strip().lower()
    if plot_choice in ("y", "yes"):
        dist_km = res.traces["distance_m"] / 1000
        soc_wh = wh_from_joules(res.traces["Ebat_J"])
        bat_max_wh = wh_from_joules(params.bat_max_energy)
        generate_plots(dist_km, res, soc_wh, config.min_soc, bat_max_wh)

    return res

run_raceday_scenario(config)

Run a race day scenario using GPX data and optionally Solcast.

Parameters:

Name Type Description Default
config SimConfig

Simulation configuration.

 required

Returns:

Type Description
SimResult

Simulation result.
Source code in simulation/scenarios.py 
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
def run_raceday_scenario(config: SimConfig) -> SimResult:
    """Run a race day scenario using GPX data and optionally Solcast.

    Args:
        config: Simulation configuration.

    Returns:
        Simulation result.
    """
    print("\nLoading race day scenario...")

    # Load GPX data
    data_dir = Path(__file__).parent.parent / "data" / "asc_24_(temp)"
    gpx_path = data_dir / config.gpx_file

    if not gpx_path.exists():
        raise FileNotFoundError(f"GPX file not found: {gpx_path}")

    print(f"Loading GPX: {gpx_path.name}")
    pts = load_gpx_points(str(gpx_path), AscTrack.NashvilleToPaducah)
    dist_m, theta_deg_arr = compute_segments(pts)

    # Get GHI data
    if config.use_solcast and SOLCAST_AVAILABLE and config.solcast_api_key:
        print("Fetching GHI data from Solcast...")
        # TODO
        ghi_arr = np.linspace(700, 900, len(pts))
    else:
        print("Using mock GHI data (Solcast disabled or unavailable)")
        ghi_arr = np.linspace(700, 900, len(pts))

    # Build position-based lookup functions from GPS data
    def theta_fn(d):
        return np.interp(d, dist_m, theta_deg_arr)

    def ghi_fn(d):
        return np.interp(d, dist_m, ghi_arr)

    # Setup simulation
    dt = config.dt
    horizon = config.horizon

    # Use default vehicle params from VehicleParams dataclass
    params = VehicleParams()

    print("\nOptimization Methods:")
    print("  1. SLSQP (fast, gradient-based)")
    print("  2. Exhaustive Search (slow, brute-force)")
    choice = input("\nSelect optimization method [1/2]: ").strip()

    if choice == "2":
        method = "exhaustive"
    else:
        method = "SLSQP"

    opt_cfg = OptimizeConfig(
        dt=dt,
        horizon=horizon,
        d0=0.0,
        vmin=config.vmin,
        vmax=config.vmax,
        method=method,
        max_iter=config.max_iter,
        min_soc=config.min_soc,
    )

    if method == "SLSQP":
        print("Optimizing velocity profile using SLSQP (may take a while)...")
        best_vs, best_obj = SLSQP_velocity(opt_cfg, theta_fn, ghi_fn, params)
    else:
        print(
            "Optimizing velocity profile using Exhaustive Search (may take a while)..."
        )
        best_vs, best_obj = exhaustive_search_velocity(
            opt_cfg, theta_fn, ghi_fn, params
        )

    print(f"Optimization complete. Objective: {best_obj:.2f}")
    print("Simulating...")
    res = simulate(best_vs, dt, 0.0, theta_fn, ghi_fn, params)

    _print_results(res, best_vs, params)

    # Ask if user wants to see plots
    plot_choice = input("\nGenerate plots? (y/n): ").strip().lower()
    if plot_choice in ("y", "yes"):
        dist_km = res.traces["distance_m"] / 1000
        soc_wh = wh_from_joules(res.traces["Ebat_J"])
        bat_max_wh = wh_from_joules(params.bat_max_energy)
        generate_plots(dist_km, res, soc_wh, config.min_soc, bat_max_wh)

    return res

simulation.plots

Result visualization: speed and battery plots.

simulation.plots

Generate plots for results.

Date: 2025-11-25 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

generate_plots(dist_km, res, soc_wh, min_soc, bat_max_wh=None)

Generates speed and battery plots from simulation results.

Parameters:

Name Type Description Default
dist_km ndarray

Distance array in kilometers.

 required
res SimResult

SimResult object containing simulation traces.

 required
soc_wh ndarray

Battery state of charge array in Watt-hours.

 required
min_soc float

Minimum SOC Line

 required
bat_max_wh float

Maximum battery capacity in Wh

 None
Source code in simulation/plots.py 
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
def generate_plots(
    dist_km: np.ndarray,
    res: SimResult,
    soc_wh: np.ndarray,
    min_soc: float,
    bat_max_wh: float = None,
) -> None:
    """Generates speed and battery plots from simulation results.

    Args:
        dist_km: Distance array in kilometers.
        res: SimResult object containing simulation traces.
        soc_wh: Battery state of charge array in Watt-hours.
        min_soc: Minimum SOC Line
        bat_max_wh: Maximum battery capacity in Wh
    """
    fig, ax = plt.subplots(2, 1, figsize=(9, 8), sharex=True)

    # Align arrays: traces may differ in length by 1
    v = res.traces["v"]
    n = min(len(dist_km), len(v), len(soc_wh))
    dist_km = dist_km[:n]
    v = v[:n]
    soc_wh = soc_wh[:n]

    ax[0].step(dist_km, v, where="post", linewidth=1.5)
    ax[0].plot(dist_km, v, "o", markersize=4)
    ax[0].set_ylabel("Speed (m/s)")
    ax[0].grid(True, alpha=0.3)

    ax[1].plot(dist_km, soc_wh, marker="o", markersize=4, label="Battery SOC")
    if bat_max_wh is not None:
        threshold_wh = min_soc * bat_max_wh
        ax[1].axhline(
            y=threshold_wh,
            color="r",
            linestyle="--",
            linewidth=1.5,
            label=f"{min_soc * 100}% minimum",
        )
        ax[1].legend()
    ax[1].set_ylabel("Battery (Wh)")
    ax[1].set_xlabel("Distance (km)")
    ax[1].grid(True, alpha=0.3)

    plt.tight_layout()
    plt.show()

simulation.mock_data

YAML-based test data loader for mock scenarios.

simulation.mock_data

Mock data loader for YAML-based test scenarios.

Date: 2025-11-23 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

load_mock_csv(path)

Load a CSV mock profile for testing.

Parameters:

Name Type Description Default
path str

Path to CSV file containing mock scenario data.

 required

Returns:

Type Description
Tuple[ndarray, ndarray, Dict[str, Any]]

Tuple of (theta_deg, ghi, meta) where meta contains dt and d0.
Source code in simulation/mock_data.py 
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
def load_mock_csv(path: str) -> Tuple[np.ndarray, np.ndarray, Dict[str, Any]]:
    """Load a CSV mock profile for testing.

    Args:
        path: Path to CSV file containing mock scenario data.

    Returns:
        Tuple of (theta_deg, ghi, meta) where meta contains dt and d0.
    """
    rows = []
    with open(path, "r", encoding="utf-8") as f:
        reader = csv.DictReader(f)
        for row in reader:
            rows.append(row)
    if not rows:
        raise ValueError("CSV file is empty or invalid format.")
    # Assume all rows have the same dt and d0
    dt = float(rows[0].get("dt", 1800))
    d0 = float(rows[0].get("d0", 0.0))
    ghi = np.array([float(r["ghi"]) for r in rows], dtype=float)
    theta_deg = np.array([float(r["theta_deg"]) for r in rows], dtype=float)
    meta = {"dt": dt, "d0": d0}
    return theta_deg, ghi, meta

load_mock_yaml(path)

Load a YAML mock profile for testing (legacy).

Source code in simulation/mock_data.py 
45
46
47
48
49
50
51
52
53
54
55
56
57
def load_mock_yaml(path: str) -> Tuple[np.ndarray, np.ndarray, Dict[str, Any]]:
    """Load a YAML mock profile for testing (legacy)."""
    if yaml is None:
        raise ImportError("PyYAML is not installed.")
    with open(path, "r", encoding="utf-8") as f:
        data = yaml.safe_load(f) or {}
    ghi = np.asarray(data.get("ghi", []), dtype=float)
    theta_deg = np.asarray(data.get("theta_deg", []), dtype=float)
    meta = {
        "dt": float(data.get("dt", 1800)),
        "d0": float(data.get("d0", 0.0)),
    }
    return theta_deg, ghi, meta

simulation.map_visualization

GPX route parsing, segment computation, and Grafana map export.

simulation.map_visualization

GPX utilities + mapping for strategy XVI.

Date: 2025-11-10 Author: Midnight Sun Team #24 - MSXVI Group: Strategy_XVI

Note

This is mainly for ASC. For FSGP elevation data isn't AS important since it's a grand prix. However, would be cool to have lap by lap coloured map gradients for best speeds for FSGP.

AscTrack

Bases: IntEnum

To store track indices.

Note

This is not fully listed out since this repo still uses asc2024 track data. Will be updated when we get asc2026 track data.

Source code in simulation/map_visualization.py 
29
30
31
32
33
34
35
36
37
class AscTrack(IntEnum):
    """To store track indices.

    Note:
        This is not fully listed out since this repo still uses asc2024 track data.
        Will be updated when we get asc2026 track data.
    """

    NashvilleToPaducah = 0

haversine(lat1, lon1, lat2, lon2)

Distance in meters between two points using GPS data.

Source code in simulation/map_visualization.py 
41
42
43
44
45
46
47
def haversine(lat1: float, lon1: float, lat2: float, lon2: float) -> float:
    """Distance in meters between two points using GPS data."""
    R = 6371000.0
    lat1, lon1, lat2, lon2 = np.radians([lat1, lon1, lat2, lon2])
    dlat, dlon = lat2 - lat1, lon2 - lon1
    a = np.sin(dlat / 2.0) ** 2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon / 2.0) ** 2
    return float(2.0 * R * np.arcsin(np.sqrt(a)))

load_gpx_points(path, track)

Loads GPX points from file.

Parameters:

Name Type Description Default
path str

File path to GPX file.

 required
track IntEnum

Specific track along the route (e.g., NashvilleToPaducah).

 required

Returns:

Type Description
ndarray

np.ndarray of shape [N, 3] as (lat, lon, elevation).
Source code in simulation/map_visualization.py 
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
def load_gpx_points(path: str, track: IntEnum) -> np.ndarray:
    """Loads GPX points from file.

    Args:
        path: File path to GPX file.
        track: Specific track along the route (e.g., NashvilleToPaducah).

    Returns:
        np.ndarray of shape [N, 3] as (lat, lon, elevation).
    """
    with open(path, "r", encoding="utf-8") as f:
        gpx = gpxpy.parse(f)

    pts = []
    track_idx = int(track)

    if track_idx >= len(gpx.tracks):
        raise ValueError(f"Track index {track_idx} not found in GPX file")

    trk = gpx.tracks[track_idx]

    for seg in trk.segments:
        for p in seg.points:
            pts.append((p.latitude, p.longitude, p.elevation))

    if not pts:
        raise ValueError("No GPX points found in the specified track.")

    return np.asarray(pts, dtype=float)

compute_segments(pts)

Computes cumulative distances and grade angles between GPS points.

Parameters:

Name Type Description Default
pts ndarray

[N,3] array of (lat, lon, elevation).

 required

Returns:

Type Description
Tuple[ndarray, ndarray]

(distance, grade_deg) - cumulative distance and grade in degrees.
Source code in simulation/map_visualization.py 
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
def compute_segments(pts: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
    """Computes cumulative distances and grade angles between GPS points.

    Args:
        pts: [N,3] array of (lat, lon, elevation).

    Returns:
        (distance, grade_deg) - cumulative distance and grade in degrees.
    """
    lat, lon, ele = pts[:, 0], pts[:, 1], pts[:, 2]
    N = len(pts)
    dist = np.zeros(N, dtype=float)
    for i in range(1, N):
        dist[i] = dist[i - 1] + haversine(lat[i - 1], lon[i - 1], lat[i], lon[i])

    # Compute grade
    d_ele = np.diff(ele, prepend=ele[0])
    d_dist = np.diff(dist, prepend=0.0)

    grade_deg = np.zeros(N, dtype=float)
    mask = d_dist > 1e-6  # Only compute where there's meaningful distance
    grade_deg[mask] = np.degrees(np.arctan(d_ele[mask] / d_dist[mask]))

    return dist, grade_deg

interpolate_to_time_grid(theta_deg, ghi, dist_m, dt, avg_speed=15.0)

Interpolates GPS-based data to match simulation timesteps.

Parameters:

Name Type Description Default
theta_deg ndarray

Grade angles at GPS points.

 required
ghi ndarray

Global horizontal irradiance at GPS points.

 required
dist_m ndarray

Cumulative distance at GPS points.

 required
dt float

Time step for simulation (seconds).

 required
avg_speed float

Estimated average speed for calculating time grid (m/s).

 15.0

Returns:

Type Description
Tuple[ndarray, ndarray]

(theta_interp, ghi_interp) - interpolated arrays matching simulation steps.
Source code in simulation/map_visualization.py 
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
def interpolate_to_time_grid(
    theta_deg: np.ndarray,
    ghi: np.ndarray,
    dist_m: np.ndarray,
    dt: float,
    avg_speed: float = 15.0,
) -> Tuple[np.ndarray, np.ndarray]:
    """Interpolates GPS-based data to match simulation timesteps.

    Args:
        theta_deg: Grade angles at GPS points.
        ghi: Global horizontal irradiance at GPS points.
        dist_m: Cumulative distance at GPS points.
        dt: Time step for simulation (seconds).
        avg_speed: Estimated average speed for calculating time grid (m/s).

    Returns:
        (theta_interp, ghi_interp) - interpolated arrays matching simulation steps.
    """
    # TODO
    return theta_deg, ghi

color_for_speed(v, vmin=10.0, vmax=20.0)

Returns color gradient hex code based on speed.

Parameters:

Name Type Description Default
v float

Current speed (m/s).

 required
vmin float

Minimum speed in range.

 10.0
vmax float

Maximum speed in range.

 20.0

Returns:

Type Description
str

Hex color code (e.g., '#ff0000').
Source code in simulation/map_visualization.py 
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
def color_for_speed(v: float, vmin: float = 10.0, vmax: float = 20.0) -> str:
    """Returns color gradient hex code based on speed.

    Args:
        v: Current speed (m/s).
        vmin: Minimum speed in range.
        vmax: Maximum speed in range.

    Returns:
        Hex color code (e.g., '#ff0000').
    """
    t = float((v - vmin) / max(vmax - vmin, 1e-9))
    t = min(max(t, 0.0), 1.0)

    # Color gradient: blue -> cyan -> green -> yellow -> red
    stops = [
        (0.00, (0, 0, 255)),
        (0.25, (0, 255, 255)),
        (0.50, (0, 255, 0)),
        (0.75, (255, 255, 0)),
        (1.00, (255, 0, 0)),
    ]

    for i in range(1, len(stops)):
        if t <= stops[i][0]:
            t0, c0 = stops[i - 1]
            t1, c1 = stops[i]
            u = (t - t0) / max(t1 - t0, 1e-9)
            r = int(c0[0] + u * (c1[0] - c0[0]))
            g = int(c0[1] + u * (c1[1] - c0[1]))
            b = int(c0[2] + u * (c1[2] - c0[2]))
            return f"#{r:02x}{g:02x}{b:02x}"

    return "#ff0000"

export_grafana_json(gpx_path, track, speed, theta_deg=None, ghi=None, vmin=10.0, vmax=20.0, outfile=None)

Builds a JSON object suitable for Grafana map visualization.

Each entry includes lat/lon, segment info, and color-coded speed.

Parameters:

Name Type Description Default
gpx_path str

Path to GPX file.

 required
track AscTrack

Track index to load.

 required
speed ndarray

Numpy array of target speeds (m/s).

 required
theta_deg Optional[ndarray]

Grade angles.

 None
ghi Optional[ndarray]

Global horizontal irradiance.

 None
vmin float

Minimum speed for color scaling.

 10.0
vmax float

Maximum speed for color scaling.

 20.0
outfile Optional[str]

Optional output file path for JSON.

 None

Returns:

Type Description
Dict

Dictionary containing segment data for Grafana.
Source code in simulation/map_visualization.py 
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
def export_grafana_json(
    gpx_path: str,
    track: AscTrack,
    speed: np.ndarray,
    theta_deg: Optional[np.ndarray] = None,
    ghi: Optional[np.ndarray] = None,
    vmin: float = 10.0,
    vmax: float = 20.0,
    outfile: Optional[str] = None,
) -> Dict:
    """Builds a JSON object suitable for Grafana map visualization.

    Each entry includes lat/lon, segment info, and color-coded speed.

    Args:
        gpx_path: Path to GPX file.
        track: Track index to load.
        speed: Numpy array of target speeds (m/s).
        theta_deg: Grade angles.
        ghi: Global horizontal irradiance.
        vmin: Minimum speed for color scaling.
        vmax: Maximum speed for color scaling.
        outfile: Optional output file path for JSON.

    Returns:
        Dictionary containing segment data for Grafana.
    """
    # Load GPS points
    pts = load_gpx_points(gpx_path, track)

    # Align array lengths
    n_segments = min(len(pts) - 1, len(speed), len(theta_deg), len(ghi))

    # Build segment data
    data: List[Dict] = []
    for i in range(n_segments):
        entry = {
            "segment": int(i),
            "lat_start": float(pts[i, 0]),
            "lon_start": float(pts[i, 1]),
            "lat_end": float(pts[i + 1, 0]),
            "lon_end": float(pts[i + 1, 1]),
            "grade_deg": float(theta_deg[i]),
            "ghi_wm2": float(ghi[i]),
            "speed_mps": float(speed[i]),
            "color": color_for_speed(speed[i], vmin, vmax),
        }
        data.append(entry)

    json_obj = {
        "segments": data,
        "meta": {"vmin": vmin, "vmax": vmax, "total_segments": n_segments},
    }

    if outfile:
        with open(outfile, "w", encoding="utf-8") as f:
            json.dump(json_obj, f, indent=2)
        print(f"Saved Grafana JSON to {outfile}")

    return json_obj