Omawianym procesem jest reaktor fermentacji drożdży (Saccharomyces cerevisiae) omówiony w przykładzie 1. Choć z punktu widzenia regulacji predykcyjnej reaktor fermentacji jest procesem jednowymiarowym (wejściem jest natężenie przepływu chłodziwa \(F_{\mathrm{ag}}\), zmienną wyjściową jest temperatura w reaktorze \(T_{\mathrm{r}}\)), do celów optymalizacji punktu pracy uwzględnia się również główne zakłócenie procesu, którym jest temperatura surowca \(T_{\mathrm{in}}\). Zmiany wartości temperatury surowca oznaczają konieczność cyklicznego obliczania puktu pracy dla algorytmu regulacji, czyli wartości zadanej temperatury [Ł11a]. Struktura układu optymalizacji punktu pracy i regulacji reaktora fermentacji jest schematycznie przedstawiona na rys. 30. Nieliniowa charakterystyka procesu \(T_{\mathrm{in}}(F_{\mathrm{ag}},T_{\mathrm{in}})\) została przedstawiona na rys. 31.
W porównaniu z przykładem 1, modele muszą mieć dwa wejścia (\(F_{\mathrm{ag}}\) i \(T_{\mathrm{in}}\)) oraz jedno wyjście (\(T_{\mathrm{r}}\)). Do celów optymalizacji punktu pracy potrzebny jest model statyczny, do regulacji predykcyjnej potrzebny jest model dynamiczny. Ponieważ zakresy zmienności sygnałów są różne, zostają one przeskalowane: \(u=0{,}01(F_{\mathrm{ag}}-F_{\mathrm{ag,nom}})\), \(h=0{,}2(T_{\mathrm{in}}-T_{\mathrm{in,nom}})\), \(y=0{,}1(T_{\mathrm{r}}-T_{\mathrm{r,nom}})\), gdzie wielkości \(F_{\mathrm{ag,nom}}=18\) l/h, \(T_{\mathrm{in,nom}}=25\) \(^{\circ}\mathrm{C}\), \(T_{\mathrm{r,nom}}=29.573212\) \(^{\circ}\mathrm{C}\) odpowiadają nominalnemu punktowi pracy.
Na podstawie dynamicznego fizykalnego modelu procesu otrzymano fizykalny model statyczny. Służy on do generacji danych, wykorzystywanych do identyfikacji modeli statycznych. Każdy ze zbiorów danych (zbiór danych uczących, weryfikujących i testowych) zawiera 2000 punktów, które są równomiernie rozłożone na przedstawionej na rys. 31 charakterystyce statycznej.
Model statyczny można opisać ogólnym równaniem \begin{equation} y^{\mathrm{ss}}=f^{\mathrm{ss}}(u^{\mathrm{ss}},h^{\mathrm{ss}}) \label{w_modstatsn_przyk} \end{equation} Wyznaczono wiele różnych modeli. Do zastosowania w algorytmach optymalizacji punktu pracy wybrano statyczny model neuronowy (sieć perceptronowa) zawierający jedynie 3 neurony w warstwie ukrytej. Błąd modelu dla zbioru danych weryfikujących wynosi \(6{,}6431 \cdot 10^{-2}\), model ma 13 parametrów (wag). Model wielomianowy dziesiątego rzędu ma porównywalny błąd (\(6{,}7234 \cdot 10^{-2}\)), ale ma aż 121 parametrów.
W wyniku symulacji fizykalnego modelu procesu otrzymano dane wykorzystywane do identyfikacji modeli dynamicznych procesu. Zostały one przedstawione na rys. 32, każdy ze zbiorów zawiera po 4000 próbek.
Wybrano model o dynamice drugiego rzędu \begin{align} y(k)=f(&u(k-1),u(k-2),h(k-1),h(k-2),\nonumber\\ &y(k-1),y(k-2))\label{w_moddynsn_przyk} \end{align} zawierający 3 neurony w warstwie ukrytej. Wyznaczono również model liniowy drugiego rzędu \begin{align} y(k)=&b_1 u(k-1)+b_2 u(k-2)+c_1 h(k-1)+c_2 h(k-2)\nonumber\\ &-a_1 y(k-1)-a_2y(k-2)\nonumber \end{align} Choć model liniowy ma te same argumenty co model nieliniowy, z uwagi na nieliniowy charakter procesu, jest on bardzo niedokładny. Dla zbioru danych weryfikujących błąd dynamicznego modelu liniowego wynosi \(6{,}9683 \cdot 10^{1}\), podczas gdy błąd dla modelu neuronowego ma wartość \(2{,}3435 \cdot 10^{-1}\). Na rys. 33 przedstawiono wyjście modelu neuronowego i modelu liniowego dla zbioru danych weryfikujących.
Porównano działanie dwóch struktur sterowania
W pierwszej strukturze rozwiązuje się dwa nieliniowe zadania optymalizacji w każdej iteracji. W drugiej strukturze liniowe i kwadratowe zadanie optymalizacji jest rozwiązywane w każdej iteracji, co 100 iteracji wywoływana jest procedura nieliniowej optymalizacji. W regulacji predykcyjnej wykorzystano neuronowy model dynamiczny \eqref{w_moddynsn_przyk}, w optymalizacji punktu pracy zastosowano neuronowy model statyczny \eqref{w_modstatsn_przyk}.
W celu maksymalizacji zysku osiąganego z produkcji, w optymalizacji punktu pracy minimalizuje się funkcję celu \begin{equation} J_{\mathrm{E}}=-F_{\mathrm{ag}}\nonumber \end{equation} W optymalizacji punktu pracy i regulacji predykcyjnej uwzględnia się te same ograniczenia zmiennej wejściowej procesu \begin{equation} F_{\mathrm{ag}}^{\min}=0 \ \mathrm{l/h}, \ F_{\mathrm{ag}}^{\max}=100 \ \mathrm{l/h}\nonumber \end{equation} Dodatkowo uwzględnia się ograniczenie zmiennej wyjściowej procesu \begin{equation} T_{\mathrm{r}} \le 25 \ ^{\circ}\mathrm{C}\nonumber \end{equation} Scenariusz zmian zakłóceń, wymuszający cykliczną optymalizację punktu pracy, ma postać \begin{equation} T_{\mathrm{in}}(k)=\left\{ \begin{array}{ll} T_{\mathrm{in,nom}} & \quad k<100\\ T_{\mathrm{in,nom}}+5(\sin(2(k-100))) \ \mathrm{l/h} & \quad 100\leq k\leq 600 \end{array}\nonumber \right. \end{equation}
Oba algorytmy regulacji predykcyjnej (NO i NPL) mają takie same parametry: \(N=10\), \(N_{\mathrm{u}}=2\), \(\lambda=2\). Na rys. 34 przedstawiono wyniki symulacji obu struktur sterowania. Na początku warstwa nieliniowej optymalizacji punktu pracy (LSSO) wyznacza optymalny punkt pracy (\(F_{\mathrm{ag}}^{\min}=65{,}3764 \ \mathrm{l/h}\), \(T_{\mathrm{r}}^{\min}=25 \ ^{\circ}\mathrm{C}\)), różny od początkowego punktu pracy. Po około 70 iteracjach proces znajduje się w optymalnym punkcie pracy. W drugiej strukturze nieliniowa optymalizacja punktu pracy wywoływana jest dla \(k=100,200,300,400,500,600\), we wszystkich pozostałych iteracjach wykorzystywana jest pomocnicza optymalizacja punktu pracy z linearyzacją modelu statycznego.
Trajektorie otrzymane w strukturze z pomocniczą optymalizacją punktu pracy współpracującą z suboptymalnym algorytmem NPL są bardzo podobne do tych, które otrzymano w złożonej obliczeniowo strukturze z nieliniową optymalizacją punktu pracy stosowaną w każdej iteracji i nieliniowym algorytmem NO. Dla całego horyzontu symulacji (600 iteracji) wyznaczono ekonomiczny wskaźnik jakości \begin{equation} J_{\mathrm{E}}=\sum_{k=1}^{600} J_{\mathrm{E}}(k)=\sum_{k=1}^{600}(-F_{\mathrm{ag}}(k))\nonumber \end{equation} W obu przypadkach otrzymano zbliżone wartości: dla pierwszej struktury otrzymano \(J_{\mathrm{E}}=-2{,}2778\cdot 10^6\), dla drugiej struktury otrzymano \(J_{\mathrm{E}}=-2{,}3137\cdot 10^6\).
Różne scenariusze optymalizacji punktu pracy reaktora fermentacji przedstawiono w pracy [Ł11a].
Poprzednia strona: 7. Algorytmy optymalizacji punktu pracy | Następna strona: 9. Podsumowanie |