diff --git a/README.md b/README.md index ee8b75d..a9b8e61 100644 --- a/README.md +++ b/README.md @@ -80,7 +80,7 @@ commands back on (e.g. `{"TempCtrl": {"Soll": {"Temp": 65}}}`). - Python 3.8+ - Server (`brewpi/requirements.txt`): `numpy`, `scipy`, `websockets`, `dpath`, and `pyserial`/`spidev` if using real hardware backends. -- Client (`client/requirements.txt`): `PyQt5` (and `PyQt5.Qwt` for plotting). +- Client (`client/requirements.txt`): `PyQt5`. Install with: diff --git a/client/brewpi.ui b/client/brewpi.ui index 90641b6..5462c6c 100644 --- a/client/brewpi.ui +++ b/client/brewpi.ui @@ -517,7 +517,7 @@ Controller - + 20 @@ -526,20 +526,14 @@ 281 - - 0.000000000000000 + + 0 - - 100.000000000000000 + + 100 - - 7 - - - 3 - - - 0.000000000000000 + + Qt::Vertical @@ -600,7 +594,7 @@ Heater - + 10 @@ -609,20 +603,14 @@ 281 - - 0.000000000000000 + + 0 - - 100.000000000000000 + + 100 - - 7 - - - 3 - - - 0.000000000000000 + + Qt::Vertical @@ -674,7 +662,7 @@ Activate - + 40 @@ -683,20 +671,14 @@ 281 - - 0.000000000000000 + + 0 - - 100.000000000000000 + + 100 - - 7 - - - 3 - - - 0.000000000000000 + + Qt::Vertical @@ -3745,13 +3727,6 @@ - - - QwtSlider - QWidget -
qwt_slider.h
- - diff --git a/client/brewpi_gui.py b/client/brewpi_gui.py index 56dbb27..c0b16ad 100755 --- a/client/brewpi_gui.py +++ b/client/brewpi_gui.py @@ -1,5 +1,5 @@ #!/usr/bin/env python3 -from brewpi_win import Ui_MainWindow +from main_window import Ui_MainWindow from PyQt5 import QtWidgets, QtGui import sys from ws.client.ws_client import WsClient @@ -105,7 +105,7 @@ class Window(QtWidgets.QMainWindow, Ui_MainWindow): if 'Temp' in submsg: self.lcdNumber_temp_soll.display(str(submsg['Temp'])) if self.slider_temp_soll_initial_update: - self.Slider_temp_soll.setValue(submsg['Temp']) + self.Slider_temp_soll.setValue(int(round(submsg['Temp']))) self.slider_temp_soll_initial_update = False if 'Rate' in submsg: subsubmsg = submsg['Rate'] @@ -133,14 +133,14 @@ class Window(QtWidgets.QMainWindow, Ui_MainWindow): elif "Power" in key: self.lcdNumber_power_heater.display(msg['Power']) if self.slider_pwr_initial_update: - self.Slider_pwr_soll.setValue(msg['Power']) + self.Slider_pwr_soll.setValue(int(round(msg['Power']))) self.slider_pwr_initial_update = False elif "Capabilities" in key: submsg = msg['Capabilities'] if "Power" in submsg: submsg = submsg['Power'] - self.Slider_pwr_soll.setLowerBound(submsg['Min']) - self.Slider_pwr_soll.setUpperBound(submsg['Max']) + self.Slider_pwr_soll.setMinimum(int(submsg['Min'])) + self.Slider_pwr_soll.setMaximum(int(submsg['Max'])) def on_stirrer_changed(self, msg): print("on_stirrer_changed {}".format(msg)) @@ -151,14 +151,14 @@ class Window(QtWidgets.QMainWindow, Ui_MainWindow): self.checkbox_stirrer_activate_initial_update = False elif "Speed" in key: if self.slider_speed_initial_update: - self.Slider_speed_soll.setValue(msg['Speed']) + self.Slider_speed_soll.setValue(int(round(msg['Speed']))) self.slider_speed_initial_update = False elif "Capabilities" in key: submsg = msg['Capabilities'] if "Power" in submsg: submsg = submsg['Power'] - self.Slider_speed_soll.setLowerBound(submsg['Min']) - self.Slider_speed_soll.setUpperBound(submsg['Max']) + self.Slider_speed_soll.setMinimum(int(submsg['Min'])) + self.Slider_speed_soll.setMaximum(int(submsg['Max'])) def closeEvent(self, event): print("Exitting gracefully!") diff --git a/client/brewpi_win.py b/client/main_window.py similarity index 98% rename from client/brewpi_win.py rename to client/main_window.py index 93f3f51..0c811d6 100644 --- a/client/brewpi_win.py +++ b/client/main_window.py @@ -2,9 +2,10 @@ # Form implementation generated from reading ui file 'brewpi.ui' # -# Created by: PyQt5 UI code generator 5.14.1 +# Created by: PyQt5 UI code generator 5.15.10 # -from PyQt5.Qwt import * +# WARNING: Any manual changes made to this file will be lost when pyuic5 is +# run again. Do not edit this file unless you know what you are doing. from PyQt5 import QtCore, QtGui, QtWidgets @@ -189,13 +190,11 @@ class Ui_MainWindow(object): self.groupBox_2.setSizePolicy(sizePolicy) self.groupBox_2.setMinimumSize(QtCore.QSize(0, 0)) self.groupBox_2.setObjectName("groupBox_2") - self.Slider_temp_soll = QwtSlider(self.groupBox_2) + self.Slider_temp_soll = QtWidgets.QSlider(self.groupBox_2) self.Slider_temp_soll.setGeometry(QtCore.QRect(20, 70, 63, 281)) - self.Slider_temp_soll.setLowerBound(0.0) - self.Slider_temp_soll.setUpperBound(100.0) - self.Slider_temp_soll.setScaleMaxMajor(7) - self.Slider_temp_soll.setScaleMaxMinor(3) - self.Slider_temp_soll.setScaleStepSize(0.0) + self.Slider_temp_soll.setMinimum(0) + self.Slider_temp_soll.setMaximum(100) + self.Slider_temp_soll.setOrientation(QtCore.Qt.Vertical) self.Slider_temp_soll.setObjectName("Slider_temp_soll") self.label_6 = QtWidgets.QLabel(self.groupBox_2) self.label_6.setGeometry(QtCore.QRect(10, 370, 67, 31)) @@ -214,13 +213,11 @@ class Ui_MainWindow(object): self.horizontalLayout.addWidget(self.groupBox_2) self.groupBox = QtWidgets.QGroupBox(self.horizontalTabWidgetPage1) self.groupBox.setObjectName("groupBox") - self.Slider_pwr_soll = QwtSlider(self.groupBox) + self.Slider_pwr_soll = QtWidgets.QSlider(self.groupBox) self.Slider_pwr_soll.setGeometry(QtCore.QRect(10, 70, 81, 281)) - self.Slider_pwr_soll.setLowerBound(0.0) - self.Slider_pwr_soll.setUpperBound(100.0) - self.Slider_pwr_soll.setScaleMaxMajor(7) - self.Slider_pwr_soll.setScaleMaxMinor(3) - self.Slider_pwr_soll.setScaleStepSize(0.0) + self.Slider_pwr_soll.setMinimum(0) + self.Slider_pwr_soll.setMaximum(100) + self.Slider_pwr_soll.setOrientation(QtCore.Qt.Vertical) self.Slider_pwr_soll.setObjectName("Slider_pwr_soll") self.checkbox_heater_activate = QtWidgets.QCheckBox(self.groupBox) self.checkbox_heater_activate.setGeometry(QtCore.QRect(30, 370, 81, 31)) @@ -235,13 +232,11 @@ class Ui_MainWindow(object): self.checkbox_stirrer_activate = QtWidgets.QCheckBox(self.groupBox_3) self.checkbox_stirrer_activate.setGeometry(QtCore.QRect(40, 370, 81, 23)) self.checkbox_stirrer_activate.setObjectName("checkbox_stirrer_activate") - self.Slider_speed_soll = QwtSlider(self.groupBox_3) + self.Slider_speed_soll = QtWidgets.QSlider(self.groupBox_3) self.Slider_speed_soll.setGeometry(QtCore.QRect(40, 70, 63, 281)) - self.Slider_speed_soll.setLowerBound(0.0) - self.Slider_speed_soll.setUpperBound(100.0) - self.Slider_speed_soll.setScaleMaxMajor(7) - self.Slider_speed_soll.setScaleMaxMinor(3) - self.Slider_speed_soll.setScaleStepSize(0.0) + self.Slider_speed_soll.setMinimum(0) + self.Slider_speed_soll.setMaximum(100) + self.Slider_speed_soll.setOrientation(QtCore.Qt.Vertical) self.Slider_speed_soll.setObjectName("Slider_speed_soll") self.label_11 = QtWidgets.QLabel(self.groupBox_3) self.label_11.setGeometry(QtCore.QRect(0, 30, 141, 31)) @@ -1277,7 +1272,7 @@ class Ui_MainWindow(object): self.label_5.setText(_translate("MainWindow", "Power [W]")) self.label_7.setText(_translate("MainWindow", "State")) self.label_state.setText(_translate("MainWindow", "State")) - self.plainTextUri.setPlainText(_translate("MainWindow", "ws://localhost:8765")) + self.plainTextUri.setPlainText(_translate("MainWindow", "ws://brewpi:8765")) self.label_8.setText(_translate("MainWindow", "Host")) self.menu_File.setTitle(_translate("MainWindow", "&File")) self.menu_Help.setTitle(_translate("MainWindow", "&Help")) diff --git a/components/pid/TODO.md b/components/pid/TODO.md new file mode 100644 index 0000000..460adc1 --- /dev/null +++ b/components/pid/TODO.md @@ -0,0 +1,47 @@ +# components/pid design backlog + +Findings from a design review, not yet acted on. See git history on +`claude_refactor` for items already fixed (Smith predictor's delay +correction enabled, the two TempController classes deduplicated into +`temp_controller_base.py`). + +- [x] **FSM thresholds aren't configurable.** Moved into an optional + `TempCtrl.Thresholds` config section (`HoldIdle`, `HoldHeat`, + `IdleHeat`, `IdleHold`, `HeatHold`, `HeatIdle`), merged over + `DEFAULT_THRESHOLDS` in `tc_constants.py` so existing configs without + the section keep working unchanged. + +- [ ] **`Pid.scale()` is a gain-scheduling hack, not anti-windup.** + `pid.py`'s `scale(k)` multiplies every gain by `1/heatrate_soll_set` + (called from `temp_controller*.py`'s `process()`), but it's framed + as overshoot compensation, conflated with the real anti-windup logic + a few lines below. `self.k` also never resets in `reset()`, so a + stale scale factor can survive a state-transition reset. Consider + moving this scheduling logic out of the generic `Pid` class and into + the temp controller, and resetting `k` in `reset()`. + +- [ ] **Three Kalman filters share one tuning.** In + `temp_controller_smith.py`, `kalman_model`, `kalman_model_delay`, + and `kalman_plant` are all constructed with the same + `params['Kalman']` despite tracking signals with different noise + characteristics. Allow independent `var_P`/`var_Q`/`var_R` per + filter in config. + +- [ ] **matplotlib imported at module level in production code.** + `temp_controller.py`, `temp_controller_smith.py`, and `kalman.py` + all `from matplotlib.pyplot import ...` just to support their + `__main__` self-test plots, which couples the live server's import + graph (and its dependency list) to a GUI plotting library it never + uses at runtime. Move the plotting demos into separate scripts. + +- [ ] **No automated tests.** The `__main__` blocks (and + `kalman_eval.py`) are manual, eyeballed-plot harnesses. Add + regression tests for FSM transitions, anti-windup clamping, and + Kalman convergence — this code drives a physical heater. + +- [ ] **Config access is unchecked `dict[key]` everywhere.** + `params['Kalman']`, `model_params['gain']`, etc., throughout, with + no schema validation at load time. We already hit this bug class + once (`config.json.sim`'s `gain`/`Model` nesting mismatch). A small + schema/dataclass validation layer at config-load would surface + errors immediately instead of mid-`__init__`. diff --git a/components/pid/kalman.py b/components/pid/kalman.py index cb90568..c31358f 100644 --- a/components/pid/kalman.py +++ b/components/pid/kalman.py @@ -1,6 +1,5 @@ import numpy as np from numpy.linalg import inv -from matplotlib.pyplot import plot, figure, subplot, grid, show class Kalman: @@ -76,57 +75,3 @@ class Kalman: I = np.eye(self.N) self.P = (I - K * self.H) * P return self.X - - -# Main -if __name__ == '__main__': - dt = 1.0 - - params = { - 'var_P' : 1, - 'var_Q' : 0.0001, - 'var_R' : 0.5 - } - k = Kalman(dt, params) - - _x1 = np.empty(0) - _y1 = np.empty(0) - _x2 = np.empty(0) - _y2 = np.empty(0) - - N = int(1000/dt) - seqn = range(0, N) - - _seqn = range(0, 2*N) - X = np.array([dt, 0]) - for n in seqn: - Z = k.process_measurement(X, 0.5) - # Kalman.print("Z:", Z) - Xp = k.process(Z) - - _x1 = np.append(_x1, Z[0]) - _x2 = np.append(_x2, Z[1]) - _y1 = np.append(_y1, Xp[0]) - _y2 = np.append(_y2, Xp[1]) - - for n in seqn: - X[0] = X[0] + dt - Z = k.process_measurement(X, 0.5) - # Kalman.print("Z:", Z) - Xp = k.process(Z) - - _x1 = np.append(_x1, Z[0]) - _x2 = np.append(_x2, Z[1]) - _y1 = np.append(_y1, Xp[0]) - _y2 = np.append(_y2, Xp[1]) - - figure(1) - subplot(2, 1, 1) - plot(_seqn, _x1, 'bx', _seqn, _y1, '-r', linewidth=1) - grid(True) - subplot(2, 1, 2) - plot(_seqn, _x2, 'bx', _seqn, _y2, '-r', linewidth=1) - grid(True) - show() - - print("End of program") diff --git a/components/pid/tc_constants.py b/components/pid/tc_constants.py index 9703336..0809cd5 100644 --- a/components/pid/tc_constants.py +++ b/components/pid/tc_constants.py @@ -8,12 +8,14 @@ class States(Enum): HOLD = 2 -THRESH_HOLD_IDLE = 0.1 -THRESH_HOLD_HEAT = 1.0 -THRESH_IDLE_HEAT = 1.0 -THRESH_IDLE_HOLD = 0.1 -THRESH_HEAT_HOLD = 1.0 -THRESH_HEAT_IDLE = 1.0 +DEFAULT_THRESHOLDS = { + "HoldIdle": 0.1, + "HoldHeat": 1.0, + "IdleHeat": 1.0, + "IdleHold": 0.1, + "HeatHold": 1.0, + "HeatIdle": 1.0 +} class Test: diff --git a/components/pid/temp_controller.py b/components/pid/temp_controller.py index 0b34212..3f50315 100644 --- a/components/pid/temp_controller.py +++ b/components/pid/temp_controller.py @@ -1,75 +1,18 @@ -from components.plant.pot import Pot -from matplotlib.pyplot import plot, figure, subplot, grid, show, legend -from components import APid -from components.pid import Pid, Kalman -from components.pid.tc_constants import * -import numpy as np +from components.pid.temp_controller_base import TempControllerBase -class TempController(APid): +class TempController(TempControllerBase): def __init__(self, dt, params): - APid.__init__(self) - self.pid_hold = Pid(dt) - self.pid_rate = Pid(dt) - self.theta_ist_set = 0 - self.theta_soll_set = 0 - self.heatrate_ist_set = 0 - self.heatrate_soll_set = 1.0 - self.heatrate_soll = 1.0 - self.theta_ist = 0 + TempControllerBase.__init__(self, dt, params) + self.dt = dt + self.last_theta_ist = 20 self.heatrate_ist = 0 - self.params = params - self.kalman = Kalman(dt, params['Kalman']) - self.y = -1 - self.state = States.INIT - self.use_kalman = True - self.pid_hold.set_params(params['Hold']) - self.pid_rate.set_params(params['Heat']) - self.is_startup = True - - def set_theta_ist(self, value): - self.theta_ist_set = value - if self.is_startup: - self.is_startup = False - self.kalman.initial((value, 0)) - - def get_theta_ist(self): - return self.theta_ist - - def set_heatrate_ist(self, value): - self.heatrate_ist_set = value - - def get_heatrate_ist(self): - return self.heatrate_ist - - def set_theta_soll(self, value): - self.theta_soll_set = value - - def get_theta_soll(self): - return self.theta_soll - - def get_theta_soll_set(self): - return self.theta_soll_set - - def set_heatrate_soll(self, value): - self.heatrate_soll_set = value - - def get_heatrate_soll(self): - return self.heatrate_soll - - def get_heatrate_soll_set(self): - return self.heatrate_soll_set - def process(self): # Process Kalman - if self.use_kalman: - Z = self.kalman.process_measurement((self.theta_ist_set, 0), 0.0) - xp = self.kalman.process(Z) - self.theta_ist = xp[0, 0] - self.heatrate_ist = xp[1, 0] * 60 - else: - self.theta_ist = self.theta_ist_set - self.heatrate_ist = self.heatrate_ist_set + self.theta_ist = self.theta_ist_set + heatrate = (self.theta_ist - self.last_theta_ist)/self.dt*60 + self.heatrate_ist = heatrate + self.last_theta_ist = self.theta_ist # Compensate for max heat rate to reduce overshoot if self.heatrate_soll_set > 0: @@ -82,118 +25,3 @@ class TempController(APid): diff = self.theta_soll_set - self.theta_ist self.process_fsm(diff) self.process_pid(theta_err, heatrate_err) - - def process_fsm(self, diff): - # Process state - state_next = self.state - if self.state == States.INIT: - if not self.is_startup: - state_next = States.IDLE - elif self.state == States.IDLE: - if diff >= THRESH_IDLE_HEAT: - state_next = States.HEAT - self.pid_rate.reset() - elif diff >= -THRESH_IDLE_HOLD: - state_next = States.HOLD - self.pid_rate.reset() - elif self.state == States.HOLD: - if diff >= THRESH_HOLD_HEAT: - state_next = States.HEAT - self.pid_rate.reset() - elif diff <= -THRESH_HOLD_IDLE: - state_next = States.IDLE - elif self.state == States.HEAT: - if diff <= -THRESH_HEAT_IDLE: - state_next = States.IDLE - elif diff <= THRESH_HEAT_HOLD: - state_next = States.HOLD - self.pid_hold.reset() - - if state_next != self.state: - self.state = state_next - print("New state = {}".format(state_next)) - - def process_pid(self, theta_err, heatrate_err): - self.pid_hold.process(theta_err, -self.theta_ist) - self.pid_rate.process(heatrate_err, -self.heatrate_ist) - - if self.state == States.IDLE: - self.y = 0 - else: - self.y = self.pid_rate.get_y() - - def get_power(self): - return self.y - - -if __name__ == '__main__': - dt = 1.0 - - temp_ist = 0 - temp_soll = 20 - ctrl = TempController(dt, Test.tc_ctrl_params) - plant = Pot(dt, Test.tc_pot_params) - _temp_ist = np.empty(0) - _temp_soll = np.empty(0) - _y = np.empty(0) - _fb = np.empty(0) - _t = np.empty(0) - _heatrate_ist_kalman = np.empty(0) - _temp_ist_kalman = np.empty(0) - a = 0.5 - fb = 0 - rho = 0.02 - temps = [{'Temp': 20, 'Duration': 1000}, {'Temp': 40, 'Duration': 1000}, {'Temp': 50, 'Duration': 1000}, {'Temp': 60, 'Duration': 1000}, {'Temp': 70, 'Duration': 1000}, {'Temp': 80, 'Duration': 1000}, {'Temp': 78, 'Duration': 1000}] - - t = 0 - for temp in temps: - temp_soll = temp['Temp'] - hold_counter = temp['Duration'] - hold = False - ctrl.set_theta_soll(temp_soll) - ctrl.set_heatrate_soll(1.0) - - while True: - if hold: - if hold_counter == 0: - break - hold_counter -= 1 - plant.process() - temp_ist = plant.get_temperature() + 0.0 * np.random.randn() - ctrl.set_theta_ist(temp_ist) - ctrl.process() - - y = 3500*ctrl.get_power() - power = max(0, y) - plant.set_power(power) - fb = plant.get_power() - if abs(temp_ist - temp_soll) < 0.1: - hold = True - - temp_ist -= rho - _temp_ist = np.append(_temp_ist, temp_ist) - _temp_soll = np.append(_temp_soll, temp_soll) - _y = np.append(_y, y) - _fb = np.append(_fb, fb) - _t = np.append(_t, t) - _heatrate_ist_kalman = np.append(_heatrate_ist_kalman, max(-1, min(3, ctrl.heatrate_ist))) - _temp_ist_kalman = np.append(_temp_ist, ctrl.theta_ist) - t += 1 - - figure(1) - subplot(3, 1, 1) - plot(_t, _temp_ist, _t, _temp_soll, 'r-', linewidth=1) - legend(["ist", "soll"]) - grid(True) - subplot(3, 1, 2) - plot(_t, _y, '-b', _t, _fb, '-r', linewidth=1) - legend(["y", "pot"]) - grid(True) - subplot(3, 1, 3) - plot(_t, _heatrate_ist_kalman, '-b', linewidth=1) - legend(["heatrate"]) - grid(True) - show() - - print("End of program") - diff --git a/components/pid/temp_controller_base.py b/components/pid/temp_controller_base.py new file mode 100644 index 0000000..36c6069 --- /dev/null +++ b/components/pid/temp_controller_base.py @@ -0,0 +1,106 @@ +from components import APid +from components.pid.pid import Pid +from components.pid.tc_constants import States, DEFAULT_THRESHOLDS + + +class TempControllerBase(APid): + def __init__(self, dt, params): + APid.__init__(self) + self.pid_hold = Pid(dt) + self.pid_rate = Pid(dt) + self.theta_ist_set = 0 + self.theta_soll_set = 0 + self.heatrate_ist_set = 0 + self.heatrate_soll_set = 1.0 + self.heatrate_soll = 1.0 + self.theta_ist = 0 + self.heatrate_ist = 0 + self.params = params + self.thresholds = {**DEFAULT_THRESHOLDS, **params.get('Thresholds', {})} + self.y = -1 + self.state = States.INIT + self.pid_hold.set_params(params['Hold']) + self.pid_rate.set_params(params['Heat']) + self.is_startup = True + + def on_state_entered(self, state): + pass + + def post_pid(self): + pass + + def set_theta_ist(self, value): + self.theta_ist_set = value + if self.is_startup: + self.is_startup = False + + def get_theta_ist(self): + return self.theta_ist + + def set_heatrate_ist(self, value): + self.heatrate_ist_set = value + + def get_heatrate_ist(self): + return self.heatrate_ist + + def set_theta_soll(self, value): + self.theta_soll_set = value + + def get_theta_soll(self): + return self.theta_soll + + def get_theta_soll_set(self): + return self.theta_soll_set + + def set_heatrate_soll(self, value): + self.heatrate_soll_set = value + + def get_heatrate_soll(self): + return self.heatrate_soll + + def get_heatrate_soll_set(self): + return self.heatrate_soll_set + + def process_fsm(self, diff): + state_next = self.state + if self.state == States.INIT: + if not self.is_startup: + state_next = States.IDLE + elif self.state == States.IDLE: + if diff >= self.thresholds['IdleHeat']: + state_next = States.HEAT + self.pid_rate.reset() + elif diff >= -self.thresholds['IdleHold']: + state_next = States.HOLD + self.pid_rate.reset() + elif self.state == States.HOLD: + if diff >= self.thresholds['HoldHeat']: + state_next = States.HEAT + self.pid_rate.reset() + elif diff <= -self.thresholds['HoldIdle']: + state_next = States.IDLE + elif self.state == States.HEAT: + if diff <= -self.thresholds['HeatIdle']: + state_next = States.IDLE + elif diff <= self.thresholds['HeatHold']: + state_next = States.HOLD + self.pid_hold.reset() + + if state_next != self.state: + self.state = state_next + print("New state = {}".format(state_next)) + self.on_state_entered(state_next) + + def process_pid(self, theta_err, heatrate_err): + self.pid_hold.process(theta_err, -self.theta_ist) + self.pid_rate.process(heatrate_err, -self.heatrate_ist) + + if self.state == States.IDLE: + self.y = 0 + else: + self.y = self.pid_rate.get_y() + + self.post_pid() + + def get_power(self): + return self.y diff --git a/components/pid/temp_controller_smith.py b/components/pid/temp_controller_smith.py index 6d33516..4576c35 100755 --- a/components/pid/temp_controller_smith.py +++ b/components/pid/temp_controller_smith.py @@ -1,249 +1,61 @@ from components.plant.pot import Pot -from matplotlib.pyplot import plot, figure, subplot, grid, show, legend -from components import APid -from components.pid import Pid, Kalman -from components.pid.tc_constants import * -import numpy as np +from components.pid.temp_controller_base import TempControllerBase +from components.pid.tc_constants import States -class TempController(APid): +class TempController(TempControllerBase): def __init__(self, dt, params, model_params): - APid.__init__(self) - self.pid_hold = Pid(dt) - self.pid_rate = Pid(dt) - self.theta_ist_set = 0 - self.theta_soll_set = 0 - self.heatrate_ist_set = 0 - self.heatrate_soll_set = 1.0 - self.heatrate_soll = 1.0 - self.theta_ist = 0 + TempControllerBase.__init__(self, dt, params) + self.dt = dt + self.last_theta_ist = 20 self.heatrate_ist = 0 - self.params = params - self.model_params = model_params - self.kalman_model = Kalman(dt, params['Kalman']) - self.kalman_model_delay = Kalman(dt, params['Kalman']) - self.kalman_plant = Kalman(dt, params['Kalman']) - self.y = -1 - self.state = States.INIT - self.use_kalman = True - self.pid_hold.set_params(params['Hold']) - self.pid_rate.set_params(params['Heat']) - self.model = Pot(dt, model_params) - self.is_startup = True + + # Fast model: same plant model but with zero transport delay, used + # to predict the current temperature without the dead time. + self.model = Pot(dt, {**model_params, 'Td': 0}) + # Delayed model: keeps the plant's assumed transport delay, so it + # can be compared like-for-like against the real (delayed) measurement. + self.model_delay = Pot(dt, model_params) + + self.theta_ist_plant = 0 + self.theta_ist_model = 0 + self.theta_ist_model_delay = 0 def set_model_power(self, power): self.model.set_power(power) + self.model_delay.set_power(power) - def set_theta_ist(self, value): - self.theta_ist_set = value - if self.is_startup: - self.is_startup = False - self.kalman_model.initial((value, 0)) - self.kalman_model_delay.initial((value, 0)) - self.kalman_plant.initial((value, 0)) + def on_state_entered(self, state): + if state == States.HEAT: + self.model.initial(self.theta_ist) + self.model_delay.initial(self.theta_ist) - def get_theta_ist(self): - return self.theta_ist - - def set_heatrate_ist(self, value): - self.heatrate_ist_set = value - - def get_heatrate_ist(self): - return self.heatrate_ist - - def set_theta_soll(self, value): - self.theta_soll_set = value - - def get_theta_soll(self): - return self.theta_soll - - def get_theta_soll_set(self): - return self.theta_soll_set - - def set_heatrate_soll(self, value): - self.heatrate_soll_set = value - - def get_heatrate_soll(self): - return self.heatrate_soll - - def get_heatrate_soll_set(self): - return self.heatrate_soll_set + def post_pid(self): + self.model.process() + self.model_delay.process() def process(self): - # Process Kalman of Plant - Z_plant = self.kalman_plant.process_measurement((self.theta_ist_set, 0), 0.0) - xp_plant = self.kalman_plant.process(Z_plant) - theta_ist_plant = xp_plant[0, 0] - heatrate_ist_plant = xp_plant[1, 0] * 60 + self.theta_ist_plant = self.theta_ist_set + self.theta_ist_model = self.model.get_temperature() + self.theta_ist_model_delay = self.model_delay.get_temperature() - # Process Kalman of Model - k_model = self.kalman_model.process_measurement((self.model.get_temperature_intermediate(), 0), 0.0) - xp_model = self.kalman_model.process(k_model) - theta_ist_model = xp_model[0, 0] - heatrate_ist_model = xp_model[1, 0] * 60 + # Smith predictor: use the fast model's prediction, corrected by the + # mismatch between the real (delayed) plant and the delayed model, + # so the dead time drops out of the feedback path. + self.theta_ist = self.theta_ist_model + (self.theta_ist_plant - self.theta_ist_model_delay) - # Process Kalman of delayed Model - k_model_delay = self.kalman_model_delay.process_measurement((self.model.get_temperature(), 0), 0.0) - xp_model_delay = self.kalman_model_delay.process(k_model_delay) - theta_ist_model_delay = xp_model_delay[0, 0] - dtheta_ist_model_delay = xp_model_delay[1, 0] * 60 - - self.theta_ist_plant = theta_ist_plant - self.dtheta_ist_plant = heatrate_ist_plant - self.theta_ist_model = theta_ist_model - self.dtheta_ist_model = heatrate_ist_model - self.theta_ist_model_delay = theta_ist_model_delay - self.dtheta_ist_model_delay = dtheta_ist_model_delay - self.theta_ist = theta_ist_plant - self.heatrate_ist = heatrate_ist_plant + heatrate = (self.theta_ist - self.last_theta_ist)/self.dt*60 + self.heatrate_ist = heatrate + self.last_theta_ist = self.theta_ist # Compensate for max heat rate to reduce overshoot if self.heatrate_soll_set > 0: self.pid_hold.scale(1.0/self.heatrate_soll_set) self.heatrate_soll = self.heatrate_soll_set * self.pid_hold.get_y() + theta_err = self.theta_soll_set - self.theta_ist + heatrate_err = self.heatrate_soll - self.heatrate_ist - if 0: - theta_err = self.theta_soll_set - (theta_ist_plant - theta_ist_model_delay + theta_ist_model) - else: - theta_err = self.theta_soll_set - theta_ist_plant - - heatrate_err = self.heatrate_soll - (heatrate_ist_plant - dtheta_ist_model_delay + heatrate_ist_model) diff = self.theta_soll_set - self.theta_ist - self.process_fsm(diff) self.process_pid(theta_err, heatrate_err) - - def process_fsm(self, diff): - # Process state - state_next = self.state - if self.state == States.INIT: - if not self.is_startup: - state_next = States.IDLE - elif self.state == States.IDLE: - if diff >= THRESH_IDLE_HEAT: - state_next = States.HEAT - self.pid_rate.reset() - elif diff >= -THRESH_IDLE_HOLD: - state_next = States.HOLD - self.pid_rate.reset() - elif self.state == States.HOLD: - if diff >= THRESH_HOLD_HEAT: - state_next = States.HEAT - self.pid_rate.reset() - elif diff <= -THRESH_HOLD_IDLE: - state_next = States.IDLE - elif self.state == States.HEAT: - if diff <= -THRESH_HEAT_IDLE: - state_next = States.IDLE - elif diff <= THRESH_HEAT_HOLD: - state_next = States.HOLD - self.pid_hold.reset() - - if state_next != self.state: - self.state = state_next - print("New state = {}".format(state_next)) - if state_next == States.HEAT: - self.model.initial(self.theta_ist) - self.kalman_model.initial((self.theta_ist, 0)) - self.kalman_model_delay.initial((self.theta_ist, 0)) - - def process_pid(self, theta_err, heatrate_err): - self.pid_hold.process(theta_err, -self.theta_ist) - self.pid_rate.process(heatrate_err, -self.heatrate_ist) - - if self.state == States.IDLE: - self.y = 0 - else: - self.y = self.pid_rate.get_y() - - self.model.process() - - def get_power(self): - return self.y - - -if __name__ == '__main__': - dt = 1.0 - - temp_soll = 20 - ctrl = TempController(dt, Test.tc_ctrl_params, Test.tc_model_params) - plant = Pot(dt, Test.tc_pot_params) - _y = np.empty(0) - _fb = np.empty(0) - _t = np.empty(0) - _temp_soll = np.empty(0) - _temp_ist_kalman = np.empty(0) - _heatrate_ist_kalman = np.empty(0) - _temp_ist_kalman_plant = np.empty(0) - _heatrate_ist_kalman_plant = np.empty(0) - _temp_ist_kalman_model = np.empty(0) - _heatrate_ist_kalman_model = np.empty(0) - a = 0.5 - fb = 0 - rho = 0.02 - temps = [{'Temp': 20, 'Duration': 1000}, {'Temp': 40, 'Duration': 1000}, {'Temp': 50, 'Duration': 1000}, {'Temp': 60, 'Duration': 1000}, {'Temp': 70, 'Duration': 1000}, {'Temp': 80, 'Duration': 1000}, {'Temp': 78, 'Duration': 1000}] - - t = 0 - - for temp in temps: - temp_soll = temp['Temp'] - hold_counter = temp['Duration'] - hold = False - ctrl.set_theta_soll(temp_soll) - ctrl.set_heatrate_soll(1.0) - - while True: - if hold: - if hold_counter == 0: - break - hold_counter -= 1 - plant.process() - temp_ist = plant.get_temperature() - ctrl.set_theta_ist(temp_ist) - ctrl.process() - - y = 3500*ctrl.get_power() - power = max(0, y) - plant.set_power(power) - fb = plant.get_power() - if abs(temp_ist - temp_soll) < 0.1: - hold = True - - _y = np.append(_y, y) - _fb = np.append(_fb, fb) - _t = np.append(_t, t) - _temp_soll = np.append(_temp_soll, temp_soll) - _temp_ist_kalman_plant = np.append(_temp_ist_kalman_plant, ctrl.theta_ist_plant) - _heatrate_ist_kalman_plant = np.append(_heatrate_ist_kalman_plant, max(-1, min(3, ctrl.dtheta_ist_plant))) - _temp_ist_kalman_model = np.append(_temp_ist_kalman_model, ctrl.theta_ist_model_delay) - _heatrate_ist_kalman_model = np.append(_heatrate_ist_kalman_model, max(-1, min(3, ctrl.dtheta_ist_model_delay))) - t += 1 - - figure(1) - subplot(3, 1, 1) - plot(_t, _temp_ist_kalman_plant, _t, _temp_soll, 'r-', linewidth=1) - legend(["ist", "soll"]) - grid(True) - subplot(3, 1, 2) - plot(_t, _y, '-b', _t, _fb, '-r', linewidth=1) - legend(["y", "pot"]) - grid(True) - subplot(3, 1, 3) - plot(_t, _heatrate_ist_kalman_plant, '-b', linewidth=1) - legend(["heatrate"]) - grid(True) - - figure(2) - subplot(2, 1, 1) - plot(_t, _temp_ist_kalman_plant, '-b', _t, _temp_ist_kalman_model, 'r-', linewidth=1) - legend(["plant", "model"]) - grid(True) - subplot(2, 1, 2) - plot(_t, _heatrate_ist_kalman_plant, '-b', _t, _heatrate_ist_kalman_model, '-r', linewidth=1) - legend(["plant", "model"]) - grid(True) - show() - - print("End of program") - diff --git a/components/plant/TODO.md b/components/plant/TODO.md new file mode 100644 index 0000000..19a8d4e --- /dev/null +++ b/components/plant/TODO.md @@ -0,0 +1,51 @@ +# components/plant design backlog + +Findings from reviewing `pot.py`'s water-bath model against real-plant +behavior (sim/real control mismatch reported by user). Not yet acted on. + +- [ ] **Heat-loss term has wrong units.** `pot.py:37`: + `leak = 1/self.M * self.L * (self.theta_amb - self.temp)/self.theta_amb`. + The power term right above it is `power/(self.M * self.C)` (W / + (J/K) -> K/s). `leak` should follow the same pattern — + `self.L * (self.theta_amb - self.temp) / (self.M * self.C)` — but + instead it's missing `/C` and divides by `theta_amb` (~20) instead, + which isn't a physically meaningful normalization. With + `C=4190` this suppresses ambient heat loss by ~3-4 orders of + magnitude versus what the power term implies, so the simulated pot + barely cools at all. The Smith predictor's internal model is also a + `Pot` instance with the identical bug, so this never showed up in + sim-vs-sim testing — only against the real plant. Likely the + primary cause of the sim/real control mismatch. + +- [ ] **`kn` is loaded but never applied.** `pot.py:18` reads + `params['kn']` but nothing in `process()`/`get_temperature()` uses + it. Git history shows it used to add Gaussian sensor noise + (`self.kn * np.random.normal(...)`) to the output, since commented + out and then dropped during the heat-diffusion refactor. Sim is + currently fully deterministic/noise-free, making the Kalman filter + and PID derivative term look better-behaved in sim than they will + against a real noisy sensor. Re-enable as measurement noise (and + decide if process noise is also wanted). + +- [ ] **`theta_amb` is hardcoded, not configured.** `brewpi.py:40` + passes a literal `20` instead of reading it from config or a live + ambient sensor. Real ambient temperature drifts and is never + modeled as a disturbance, so the controller's robustness to it is + untested. + +- [ ] **Thermal lag modeled as single-pole lag, not transport delay.** + `HeatDiffusion` (`heat_diffusion.py`) implements `alpha=dt/Td` + exponential smoothing, not true dead time. There's already a + `Delay` class (`plant/delay.py`) for pure transport delay, dropped + in favor of `HeatDiffusion` in commit `16cb3d3`. If the real pot's + heater-to-sensor coupling behaves more like dead time (heat must + propagate through the water before reaching the sensor) than a + single lag, `Td` won't represent the same dynamic in sim vs. reality, + and the Smith predictor's delay compensation (which assumes `Td` + matches the real plant) will be off. + +- [ ] **No actuator/process realism.** No heater power saturation + enforced inside `Pot` itself (handled ad hoc in test harnesses), no + evaporation/lid heat loss, no varying thermal mass `M` as volume + changes (e.g. boil-off, adding water/grain) — `M` is a fixed + constant for the whole simulation. \ No newline at end of file diff --git a/components/plant/pot.py b/components/plant/pot.py index 26d6c17..9aa9341 100644 --- a/components/plant/pot.py +++ b/components/plant/pot.py @@ -1,6 +1,6 @@ from components.aplant import APlant from components.plant.heat_diffusion import HeatDiffusion - +from components.plant import delay class Pot(APlant): def __init__(self, dt, params, theta_amb=20): @@ -10,35 +10,57 @@ class Pot(APlant): self.e = 0 self.x = 0 + # Plant power gain (efficiency) self.gain = params['gain'] - self.C = params['C'] - self.M = params['M'] - self.L = params['L'] - self.Td = params['Td'] - self.kn = params['kn'] - self.temp = params['theta'] - self.temp_intermediate = params['theta'] + # Plant specific thermal capacity [W*s/(kg*K)] + self.C = params['C'] + + # Plant mass [kg] + self.M = params['M'] + + # Plant energy loss coefficient [W/(kg*K)] + # Negative input power as a function of plant mass and temperature difference T_plant and T_ambient + # P_loss = L * Mass * (T_plant - T_ambient) + self.L = params['L'] + + # Energy transport propagation delay + self.Td = params['Td'] + + # TBD + self.kn = params['kn'] + + # Plant temperature [°C] + self.temp = params['theta'] + + # Ambient temperature [°C] self.theta_amb = theta_amb + # Set power [W] self.power_set = 0 - self.delay = HeatDiffusion(self.dt, self.Td, params['theta']) + + # Plant delay [s] + self.delay = delay.Delay(dt, self.Td, 0) + self.p_in = 0 + self.p_pot = 0 def initial(self, temp): - self.temp_intermediate = temp self.temp = temp - self.delay.initial(temp) def activate(self, enable): pass def process(self): - power = self.gain * self.power_set - leak = 1/self.M * self.L * (self.theta_amb - self.temp)/self.theta_amb - self.temp_intermediate += (power / (self.M * self.C) + leak) * self.dt + self.p_in = self.gain * self.power_set + self.delay.put(self.p_in) - # Delay - self.temp = self.delay.process(self.temp_intermediate) + p_loss = self.L * self.M * (self.temp - self.theta_amb) + self.p_pot = self.delay.get() - p_loss + + print(f"p_loss: {p_loss}") + print(f"theta_amb: {self.theta_amb}") + print(f"temp: {self.temp}") + self.temp = min(100, self.temp + self.p_pot/(self.M * self.C) * self.dt) def is_activated(self): return True @@ -52,5 +74,6 @@ class Pot(APlant): def get_temperature(self): return self.temp - def get_temperature_intermediate(self): - return self.temp_intermediate + def get_p_pot(self): + return self.p_pot + diff --git a/config.json.templ b/config.json.templ index 889b62c..c0a323c 100644 --- a/config.json.templ +++ b/config.json.templ @@ -25,6 +25,14 @@ "ki": 0.01, "kd": 0.0, "kt": 1.5 + }, + "Thresholds": { + "HoldIdle": 0.1, + "HoldHeat": 1.0, + "IdleHeat": 1.0, + "IdleHold": 0.1, + "HeatHold": 1.0, + "HeatIdle": 1.0 } }, "Model": { diff --git a/scripts/__init__.py b/scripts/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/scripts/demos/__init__.py b/scripts/demos/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/scripts/demos/pid/__init__.py b/scripts/demos/pid/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/scripts/demos/pid/demo_kalman.py b/scripts/demos/pid/demo_kalman.py new file mode 100644 index 0000000..da11e36 --- /dev/null +++ b/scripts/demos/pid/demo_kalman.py @@ -0,0 +1,56 @@ +import numpy as np +from matplotlib.pyplot import plot, figure, subplot, grid, show +from components.pid.kalman import Kalman + + +if __name__ == '__main__': + dt = 1.0 + + params = { + 'var_P' : 1, + 'var_Q' : 0.0001, + 'var_R' : 0.5 + } + k = Kalman(dt, params) + + _x1 = np.empty(0) + _y1 = np.empty(0) + _x2 = np.empty(0) + _y2 = np.empty(0) + + N = int(1000/dt) + seqn = range(0, N) + + _seqn = range(0, 2*N) + X = np.array([dt, 0]) + for n in seqn: + Z = k.process_measurement(X, 0.5) + # Kalman.print("Z:", Z) + Xp = k.process(Z) + + _x1 = np.append(_x1, Z[0]) + _x2 = np.append(_x2, Z[1]) + _y1 = np.append(_y1, Xp[0]) + _y2 = np.append(_y2, Xp[1]) + + for n in seqn: + X[0] = X[0] + dt + Z = k.process_measurement(X, 0.5) + # Kalman.print("Z:", Z) + Xp = k.process(Z) + + _x1 = np.append(_x1, Z[0]) + _x2 = np.append(_x2, Z[1]) + _y1 = np.append(_y1, Xp[0]) + _y2 = np.append(_y2, Xp[1]) + + figure(1) + subplot(2, 1, 1) + plot(_seqn, _x1, 'bx', _seqn, _y1, '-r', linewidth=1) + grid(True) + subplot(2, 1, 2) + plot(_seqn, _x2, 'bx', _seqn, _y2, '-r', linewidth=1) + grid(True) + show() + + print("End of program") diff --git a/components/pid/kalman_eval.py b/scripts/demos/pid/demo_kalman_eval.py similarity index 100% rename from components/pid/kalman_eval.py rename to scripts/demos/pid/demo_kalman_eval.py diff --git a/scripts/demos/pid/demo_temp_controller.py b/scripts/demos/pid/demo_temp_controller.py new file mode 100644 index 0000000..9e0b8a1 --- /dev/null +++ b/scripts/demos/pid/demo_temp_controller.py @@ -0,0 +1,101 @@ +import numpy as np +from matplotlib.pyplot import plot, figure, subplot, grid, show, legend +from components.plant.pot import Pot +from components.pid.temp_controller import TempController + + +if __name__ == '__main__': + ctrl_params = { + "Hold": { + "kp": 0.4, + "ki": 0.0, + "kd": 0.0, + "kt": 0.0 + }, + "Heat": { + "kp": 0.08, + "ki": 0.008, + "kd": 0.0, + "kt": 1.5 + } + } + + plant_params = { + "theta" : 20, + "C" : 4190, + "M" : 20, + "L" : 0.2, + "Td" : 30, + "kn" : 0.2, + "gain" : 1.0 + } + + dt = 1.0 + k_noise = 0.001 + + temp_ist = 0 + temp_soll = 20 + heatrate_soll = 1.25 + ctrl = TempController(dt, ctrl_params) + plant = Pot(dt, plant_params) + _y = np.empty(0) + _fb = np.empty(0) + _t = np.empty(0) + _temp_soll = np.empty(0) + _heatrate_ist = np.empty(0) + _temp_ist = np.empty(0) + a = 0.5 + fb = 0 + rho = 0.02 + temps = [{'Temp': 20, 'Duration': 1000}, {'Temp': 40, 'Duration': 1000}, {'Temp': 50, 'Duration': 1000}, {'Temp': 60, 'Duration': 1000}, {'Temp': 70, 'Duration': 1000}, {'Temp': 80, 'Duration': 1000}, {'Temp': 78, 'Duration': 1000}] + + t = 0 + for temp in temps: + temp_soll = temp['Temp'] + hold_counter = temp['Duration'] + hold = False + ctrl.set_theta_soll(temp_soll) + ctrl.set_heatrate_soll(heatrate_soll) + + while True: + if hold: + if hold_counter == 0: + break + hold_counter -= 1 + plant.process() + temp_ist = plant.get_temperature() + k_noise * np.random.randn() + ctrl.set_theta_ist(temp_ist) + ctrl.process() + + y = 3500*ctrl.get_power() + power = max(0, y) + plant.set_power(power) + fb = plant.get_power() + if abs(temp_ist - temp_soll) < 0.1: + hold = True + + temp_ist -= rho + _temp_ist = np.append(_temp_ist, temp_ist) + _temp_soll = np.append(_temp_soll, temp_soll) + _y = np.append(_y, y) + _fb = np.append(_fb, fb) + _t = np.append(_t, t) + _heatrate_ist = np.append(_heatrate_ist, max(-1, min(3, ctrl.heatrate_ist))) + t += 1 + + figure(1) + subplot(3, 1, 1) + plot(_t, _temp_ist, _t, _temp_soll, 'r-', linewidth=1) + legend(["ist", "soll"]) + grid(True) + subplot(3, 1, 2) + plot(_t, _y, '-b', _t, _fb, '-r', linewidth=1) + legend(["y", "pot"]) + grid(True) + subplot(3, 1, 3) + plot(_t, _heatrate_ist, '-b', linewidth=1) + legend(["heatrate"]) + grid(True) + show() + + print("End of program") diff --git a/scripts/demos/pid/demo_temp_controller_smith.py b/scripts/demos/pid/demo_temp_controller_smith.py new file mode 100644 index 0000000..31e7027 --- /dev/null +++ b/scripts/demos/pid/demo_temp_controller_smith.py @@ -0,0 +1,111 @@ +import numpy as np +from matplotlib.pyplot import plot, figure, subplot, grid, show, legend +from components.plant.pot import Pot +from components.pid.temp_controller_smith import TempController + + +if __name__ == '__main__': + ctrl_params = { + "Hold": { + "kp": 0.4, + "ki": 0.0, + "kd": 0.0, + "kt": 0.0 + }, + "Heat": { + "kp": 0.08, + "ki": 0.008, + "kd": 0.0, + "kt": 1.5 + } + } + + plant_params = { + "theta" : 20, + "C" : 4190, + "M" : 20, + "L" : 0.2, + "Td" : 30, + "kn" : 0.2, + "gain" : 1.0 + } + + dt = 1.0 + k_noise = 0.001 + + temp_ist = 0 + temp_soll = 20 + heatrate_soll = 1.25 + ctrl = TempController(dt, ctrl_params, plant_params) + plant = Pot(dt, plant_params) + _y = np.empty(0) + _fb = np.empty(0) + _t = np.empty(0) + _temp_soll = np.empty(0) + _heatrate_ist = np.empty(0) + _temp_ist = np.empty(0) + _temp_ist_model = np.empty(0) + _temp_ist_model_delay = np.empty(0) + a = 0.5 + fb = 0 + rho = 0.02 + temps = [{'Temp': 20, 'Duration': 1000}, {'Temp': 40, 'Duration': 1000}, {'Temp': 50, 'Duration': 1000}, {'Temp': 60, 'Duration': 1000}, {'Temp': 70, 'Duration': 1000}, {'Temp': 80, 'Duration': 1000}, {'Temp': 78, 'Duration': 1000}] + + t = 0 + for temp in temps: + temp_soll = temp['Temp'] + hold_counter = temp['Duration'] + hold = False + ctrl.set_theta_soll(temp_soll) + ctrl.set_heatrate_soll(heatrate_soll) + + while True: + if hold: + if hold_counter == 0: + break + hold_counter -= 1 + plant.process() + temp_ist = plant.get_temperature() + k_noise * np.random.randn() + ctrl.set_theta_ist(temp_ist) + ctrl.process() + + y = 3500*ctrl.get_power() + power = max(0, y) + plant.set_power(power) + ctrl.set_model_power(power) + fb = plant.get_power() + if abs(temp_ist - temp_soll) < 0.1: + hold = True + + temp_ist -= rho + _temp_ist = np.append(_temp_ist, temp_ist) + _temp_soll = np.append(_temp_soll, temp_soll) + _y = np.append(_y, y) + _fb = np.append(_fb, fb) + _t = np.append(_t, t) + _heatrate_ist = np.append(_heatrate_ist, max(-1, min(3, ctrl.heatrate_ist))) + _temp_ist_model = np.append(_temp_ist_model, ctrl.theta_ist_model) + _temp_ist_model_delay = np.append(_temp_ist_model_delay, ctrl.theta_ist_model_delay) + t += 1 + + figure(1) + subplot(3, 1, 1) + plot(_t, _temp_ist, _t, _temp_soll, 'r-', linewidth=1) + legend(["ist", "soll"]) + grid(True) + subplot(3, 1, 2) + plot(_t, _y, '-b', _t, _fb, '-r', linewidth=1) + legend(["y", "pot"]) + grid(True) + subplot(3, 1, 3) + plot(_t, _heatrate_ist, '-b', linewidth=1) + legend(["heatrate"]) + grid(True) + + figure(2) + plot(_t, _temp_ist, '-b', _t, _temp_ist_model, '-g', _t, _temp_ist_model_delay, '-r', linewidth=1) + legend(["plant", "model (fast)", "model (delayed)"]) + grid(True) + show() + + print("End of program") diff --git a/scripts/demos/plant/__init__.py b/scripts/demos/plant/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/scripts/demos/plant/demo_pot.py b/scripts/demos/plant/demo_pot.py new file mode 100644 index 0000000..21ee5a3 --- /dev/null +++ b/scripts/demos/plant/demo_pot.py @@ -0,0 +1,59 @@ +import numpy as np +from matplotlib.pyplot import plot, figure, subplot, grid, show, legend +from components.plant.pot import Pot + + +if __name__ == '__main__': + dt = 1.0 + power_on = 3500 + steps_heat = 4000 + steps_cool = 40000 + + pot_params = { + "theta" : 20, + "C" : 4190, + "M" : 20, + "L" : 0.4, + "Td" : 30, + "kn" : 0.2, + "gain" : 1.0 + } + + pot = Pot(dt, pot_params) + + _t = np.empty(0) + _temp = np.empty(0) + _p_pot = np.empty(0) + _power = np.empty(0) + + t = 0 + pot.set_power(power_on) + for _ in range(steps_heat): + pot.process() + _t = np.append(_t, t) + _temp = np.append(_temp, pot.get_temperature()) + _p_pot = np.append(_p_pot, pot.get_p_pot()) + _power = np.append(_power, pot.get_power()) + t += 1 + + pot.set_power(0) + for _ in range(steps_cool): + pot.process() + _t = np.append(_t, t) + _temp = np.append(_temp, pot.get_temperature()) + _p_pot = np.append(_p_pot, pot.get_p_pot()) + _power = np.append(_power, pot.get_power()) + t += 1 + + figure(1) + subplot(2, 1, 1) + plot(_t, _temp, '-b', linewidth=1) + legend(["temp"]) + grid(True) + subplot(2, 1, 2) + plot(_t, _power, '-g', _t, _p_pot, '-r', linewidth=1) + legend(["power", "p_pot"]) + grid(True) + show() + + print("End of program")