diff --git a/utils/replay_sim.py b/utils/replay_sim.py index 0670b2a..566d53f 100644 --- a/utils/replay_sim.py +++ b/utils/replay_sim.py @@ -1,20 +1,30 @@ #!/usr/bin/env python3 -"""Replay a BrewPi log through the temp controller with candidate PID params. +"""Closed-loop replay simulation for the BrewPi temp controller. -Feeds real temp_ist readings from a log back into a fresh TempController -(Normal mode) tick by tick and compares the resulting power output to the -original logged values. Because the real sensor measurements are injected at -each tick, no plant model is needed: the replay shows exactly how a different -set of PID gains would have responded to the same temperature trajectory. +Loads a real run log and drives a simulated Pot plant with a fresh +TempController (Normal mode) using the same setpoints (temp_soll, +heatrate_soll_set) as the original run. The controller's power output +feeds back into the plant model each tick, generating a simulated +temperature trajectory (sim_temp_ist) that can be compared against the +real one from the log (real_temp_ist). + +Use this to verify PID params: tweak gains and see whether the simulated +plant reaches setpoint faster, overshoots less, or oscillates less than +the original logged run, all without touching real hardware. Usage: python utils/replay_sim.py python utils/replay_sim.py \\ --heat-kp 0.1 --heat-ki 0.03 --heat-kt 1.0 - python utils/replay_sim.py --config alt_config.json + python utils/replay_sim.py \\ + --plant-M 25 --plant-L 0.15 --ambient 18 -By default reads PID params from config.json in the repo root. Any --
- -flag overrides just that one value; the rest come from the config. +Plant params default to water (C=4190 J/kg·K), M=20 kg, L=0.2 W/kg·K, +Td=30 s — pass --plant-* to match the actual batch. Max heater power is +inferred from the highest power_eff step in the log. + +By default reads PID params from config.json in the repo root. Any +--
- flag overrides just that one value. """ import argparse @@ -31,6 +41,7 @@ from matplotlib.ticker import AutoMinorLocator sys.path.insert(0, str(Path(__file__).parent.parent)) from components.pid.temp_controller import TempController +from components.plant.pot import Pot def _apply_gain_overrides(params, args): @@ -44,32 +55,50 @@ def _apply_gain_overrides(params, args): def _infer_max_power(samples): - """Largest discrete heater step seen in the log — used to normalise power_set.""" + """Largest discrete heater step seen in the log — heater's effective max.""" return max((s['power_eff'] for s in samples), default=1.0) or 1.0 def _infer_heatrate_soll_set(samples): - """Heat-rate setpoint from the log: max rate_soll seen while power > 0. + """Estimate heatrate_soll_set from the log. rate_soll = heatrate_soll_set * pid_hold.get_y(); when the hold PID is - saturated at 1.0 (full-power heating) the two are equal, so the maximum - over all heating samples is a good estimate of heatrate_soll_set.""" + saturated at 1.0 during active heating the two are equal, so the max + rate_soll seen while the heater is on is a good upper bound.""" candidates = [s['rate_soll'] for s in samples if s['power_set'] > 0] return max(candidates, default=1.0) or 1.0 -def run_replay(samples, pid_params, heatrate_soll_set, dt): +def run_replay(samples, pid_params, heatrate_soll_set, plant_params, ambient, max_power, dt): + """Closed-loop simulation: controller drives plant, plant feeds temp back. + + The plant is seeded at the real starting temperature from the log so the + simulation begins at a realistic operating point. Setpoints follow the + real run's temp_soll schedule sample by sample.""" tc = TempController(dt) tc.set_params(pid_params) tc.set_enabled(True) + plant = Pot(dt) + plant.set_plant_params(plant_params) + plant.set_ambient_temperature(ambient) + plant.initial(samples[0]['temp_ist']) + out = [] for s in samples: - tc.set_theta_ist(s['temp_ist']) + plant.process() + sim_temp = plant.get_temperature() + + tc.set_theta_ist(sim_temp) tc.set_theta_soll(s['temp_soll']) tc.set_heatrate_soll(heatrate_soll_set) tc.process() + + power_watts = max(0.0, tc.get_power() * max_power) + plant.set_power(power_watts) + out.append({ + 'sim_temp_ist': sim_temp, 'power': tc.get_power(), 'state': tc.state.name, 'heatrate_ist': tc.get_heatrate_ist(), @@ -94,26 +123,27 @@ def plot_replay(samples, replay, name, max_power, params_label): t = [s['t'] for s in samples] fig, (ax_temp, ax_rate, ax_power) = plt.subplots(3, 1, sharex=True, facecolor='white') - fig.suptitle('{} — replay simulation'.format(name), fontsize='small') + fig.suptitle('{} — replay simulation ({})'.format(name, params_label), fontsize='small') - # Temperature: identical in both runs (real sensor readings are injected) - ax_temp.plot(t, [s['temp_ist'] for s in samples], '-b', linewidth=1, label='theta_ist') - ax_temp.plot(t, [s['temp_soll'] for s in samples], '-r', linewidth=1, label='theta_soll') + # Temperature: real log vs simulated plant trajectory + ax_temp.plot(t, [s['temp_ist'] for s in samples], '-b', linewidth=1, label='real_temp_ist') + ax_temp.plot(t, [r['sim_temp_ist'] for r in replay], '-', color='C1', linewidth=1, label='sim_temp_ist') + ax_temp.plot(t, [s['temp_soll'] for s in samples], '-r', linewidth=1, alpha=0.5, label='temp_soll') ax_temp.set_ylabel('Temp [°C]') ax_temp.legend(loc='upper left', fontsize='small') # Heat rate - ax_rate.plot(t, [s['rate_ist'] for s in samples], '-b', linewidth=1, alpha=0.6, label='rate_ist (log)') - ax_rate.plot(t, [r['heatrate_ist'] for r in replay], '--', color='C1', linewidth=1, label='rate_ist (replay)') - ax_rate.plot(t, [s['rate_soll'] for s in samples], '-r', linewidth=1, alpha=0.4, label='rate_soll (log)') + ax_rate.plot(t, [s['rate_ist'] for s in samples], '-b', linewidth=1, alpha=0.6, label='rate_ist (real)') + ax_rate.plot(t, [r['heatrate_ist'] for r in replay], '-', color='C1', linewidth=1, label='rate_ist (sim)') + ax_rate.plot(t, [s['rate_soll'] for s in samples], '-r', linewidth=1, alpha=0.4, label='rate_soll (real)') ax_rate.set_ylabel('Rate [K/min]') ax_rate.legend(loc='upper left', fontsize='small') - # Power (both normalised to 0–1) + # Power: log (normalised) vs sim controller output log_norm = [s['power_set'] / max_power for s in samples] - replay_pwr = [r['power'] for r in replay] - ax_power.plot(t, log_norm, '-b', linewidth=1, label='power (log, norm.)') - ax_power.plot(t, replay_pwr, '-', color='C1', linewidth=1, label='power ({})'.format(params_label)) + replay_pwr = [r['power'] for r in replay] + ax_power.plot(t, log_norm, '-b', linewidth=1, label='power (real, norm.)') + ax_power.plot(t, replay_pwr, '-', color='C1', linewidth=1, label='power (sim)') ax_power.set_ylabel('Power [0–1]') ax_power.set_xlabel('t [s]') ax_power.legend(loc='upper left', fontsize='small') @@ -143,8 +173,19 @@ def main(): parser.add_argument('--config', default=str(default_config), help='Config file for base PID params (default: %(default)s)') parser.add_argument('--rate-soll', type=float, default=None, - help='Heat-rate setpoint [K/min] fed to the replay controller ' - '(default: inferred from log as max rate_soll while heating)') + help='Heat-rate setpoint [K/min] (default: inferred from log)') + parser.add_argument('--ambient', type=float, default=None, + help='Ambient temperature [°C] (default: from config ambient_temperature)') + + # Plant model params + parser.add_argument('--plant-M', type=float, default=20.0, metavar='kg', + help='Plant mass [kg] (default: %(default)s)') + parser.add_argument('--plant-C', type=float, default=4190.0, metavar='J/kgK', + help='Specific heat capacity [J/(kg·K)] (default: %(default)s — water)') + parser.add_argument('--plant-L', type=float, default=0.2, metavar='W/kgK', + help='Heat loss coefficient [W/(kg·K)] (default: %(default)s)') + parser.add_argument('--plant-Td', type=float, default=30.0, metavar='s', + help='Transport delay [s] (default: %(default)s)') for section, prefix in (('Hold', 'hold'), ('Heat', 'heat'), ('Cool', 'cool')): for gain in ('kp', 'ki', 'kd', 'kt'): @@ -173,12 +214,15 @@ def main(): with open(config_path) as f: config = json.load(f) - samples = log['Samples'] - name = log.get('Name', '{}_{}'.format(args.date_time, args.sud_name)) - pid_params = _apply_gain_overrides(config['TempCtrl'], args) - dt = (samples[1]['t'] - samples[0]['t']) if len(samples) > 1 else 1.0 - max_power = _infer_max_power(samples) - rate_soll = args.rate_soll if args.rate_soll is not None else _infer_heatrate_soll_set(samples) + samples = log['Samples'] + name = log.get('Name', '{}_{}'.format(args.date_time, args.sud_name)) + pid_params = _apply_gain_overrides(config['TempCtrl'], args) + dt = (samples[1]['t'] - samples[0]['t']) if len(samples) > 1 else 1.0 + max_power = _infer_max_power(samples) + rate_soll = args.rate_soll if args.rate_soll is not None else _infer_heatrate_soll_set(samples) + ambient = args.ambient if args.ambient is not None else config.get('ambient_temperature', 20.0) + + plant_params = {'M': args.plant_M, 'C': args.plant_C, 'L': args.plant_L, 'Td': args.plant_Td} any_override = any( getattr(args, '{}_{}'.format(p, g)) is not None @@ -190,13 +234,15 @@ def main(): print('Log: {} ({} samples, dt={:.1f} s)'.format(name, len(samples), dt)) print('Power: max={:.0f} W'.format(max_power)) print('Rate: heatrate_soll_set={:.2f} K/min'.format(rate_soll)) + print('Plant: M={M} kg C={C} J/kgK L={L} W/kgK Td={Td} s ambient={amb} °C'.format( + amb=ambient, **plant_params)) print('Params: {}'.format('overridden' if any_override else 'from config')) print() for section in ('Hold', 'Heat', 'Cool'): p = pid_params[section] print(' {}: kp={kp} ki={ki} kd={kd} kt={kt}'.format(section, **p)) - replay = run_replay(samples, pid_params, rate_soll, dt) + replay = run_replay(samples, pid_params, rate_soll, plant_params, ambient, max_power, dt) plot_replay(samples, replay, name, max_power, params_label) plt.show()