Full text
Embedded Implementation of a Neural Network emulating Nonlinear MPC in a process control application* Sebastian Leonow1and Raphael Dyrska1and Martin M¨ onnigmann1 Abstract— We present the design, training, and implementation of a nonlinear autoregressive neural network for the control of a multi-input, multi-output hydraulic plant. The network mimics the optimal control signals of a nonlinear model predictive controller and is implemented on a low-level microcontroller. While trained with simulation data only, experiments on the real plant show that not only the setpoint tracking, but to some degree also the constraint satisfaction and unmeasured disturbance rejection are adapted by the neural network. In contrast to the optimization-based predictive controller, the neural network easily runs on an ESP32 microcontroller and Micropython with guaranteed evaluation time and still achieves similar control performance as the predictive controller. I. INTRODUCTION Nonlinear model predictive control (NMPC) is a powerful method for the control of constrained nonlinear systems with multiple inputs and outputs. It is based on periodically solving a nonlinear program (NLP) for a finite prediction horizon under consideration of the system dynamics (see, e.g., [1], [2]). While powerful from a theoretical point of view, the complexity of the underlying optimization problem is often a limiting factor for its application. Several approaches exist in order to reduce the computational effort, either from an algorithmic point of view (see, e.g., [3]), or by exploiting structural information of the solution (see, e.g., [4], [5]). Especially challenging for an application on embedded hardware with strongly limited resources of computational power and memory, the embedded implementation of NMPC also received a lot of attention during the last years (see, e.g., [6], [7], [8]). Besides research focusing on NMPC directly, various approaches exist that combine machine learning methods with linear model predictive control (MPC) and NMPC. An overview of contributions combining MPC and machine learning methods is given in [9]. In [10], constraint satisfaction of a learning-based controller for a linear system is provided by introducing safety sets and a backup controller in case of constraint violation. To speed up the solution of the underlying optimization problem, the authors in [11] use classification methods to predict the set of active constraints expected for the MPC solution for the current state to warmstart the underlying active-set algorithm. *This paper is funded by the European Union’s Horizon Europe under grant no. 101079342 (Fostering Opportunities Towards Slovak Excellence in Advanced Control for Smart Industries) 1All Authors are with Department of Mechanical Engineering, Automatic Control and System Theory, Ruhr-Universit¨ at Bochum, 44801 Bochum, Germany sebatian.leonow, raphael.dyrska, [email protected] Instead of assisting in the solution process itself, machine learning, and especially neural networks are also used to replace the MPC and NMPC controller completely. In [12], [13], [14], neural networks were used to replace explicit MPC controllers, i.e., multi-parametric solutions of a linear MPC problem. In [15] and [16], the mixed integer optimizations within an MPC controller for domestic heating systems were mimicked by a neural network. The authors in [17] learned the control laws resulting from MPC for controlling the temperature of a six-zone building in a simulation case study. Another application of replacing optimal control with a neural network is given in [18] for the control of resonant power converters in a hardware-in-the-loop setup. In this paper, we present the design of a nonlinear autoregressive network with exogenous inputs that mimics the optimal control of an NMPC and discuss its implementation on an ESP32 microcontroller using Micropython as online interpreter. We compare the performance of the neural network based control to the original NMPC by applying both variants to a laboratory hydraulic plant. Besides the tracking of reference values, we also aim at learning some degree of constraint satisfaction by applying constraints to the NMPC control and using the resulting control action as target data for the neural network training. Although this procedure does not yield a guaranteed constraint satisfaction as discussed in [10], it showed satisfying results in a way that the neural network based control performs similarly to the NMPC in reference tracking, disturbance rejection, and constraint satisfaction, in the sample application (see Sec. IV). II. NEURAL NETWORK BASED CONTROL The aim of our approach is to mimic the control actions of an NMPC by a neural network, by providing control references r(t)and plant outputs y(t)as features to the network and receiving the resulting control action u(t)as target. Providing r(t)and y(t)instead of the control error e(t) = r(t)−y(t)is crucial for the network to inherit a observer (EKF) NMPC neural network Fig. 1. The neural network based controller covers observer and NMPC.
layer 1 layer 2 inputs outputs output feedback Fig. 2. NARX architecture, depicted for one neuron in each layer, with weights wand d, and biases b. The superscript indices iand jcorrespond to the enumerated neurons in each layer. The first subscript index denotes the corresponding layer, the second subscript index, if applicable, denotes the assignment of a weight (wor d) to the plant outputs yor the (fed back) network outputs u. The unit delay is denoted by z−1. The neurons from layer 1 have a nonlinear activation function h(σi), where σidenotes the sum of all weighted inputs and the bias of the current neuron i. certain degree of constraint satisfaction with respect to y(t) (see Sec. III-B). The network and NMPC controller structure, and plant interfaces are outlined in Fig. 1. We chose a nonlinear autoregressive neural network with exogeneous inputs (NARX) as fundamental neural network architecture. This choice is based primarily on previous experiments that showed promising results, in particular compared to recurrent networks, where the NARX architecture performed significantly better. We specifically used a fully connected NARX with two layers and one delay stage as depicted in Fig. 2. Layer 1 consists of 10 neurons, while layer 2 consists of 2 neurons. In the fully connected NARX all inputs r(k)y(k)and (fed back) outputs u(k−1) from the previous time step are connected to all layer 1 neurons through the corresponding weights wi 1,u and wi 1,u, respectively. One delay stage is included for all inputs and the delayed inputs are also connected to all neurons from layer 1 through the weights di 1,y and di 1,u, respectively. A total number of 2·(|r|+|y|+|u|)inputs to each neuron from layer 1 results, where |·| denotes the size of the corresponding vectors. The second layer neurons correspond with the number of outputs |u|and are connected to all layer 1 neurons. Table I summarizes the network architecture. For the weights wi 1,y ∈R(|r|+|y|)×1,wi 1,u ∈R|u|×1, di 1,y ∈R(|r|+|y|)×1, and di 1,u ∈R|u|×1follows. The two second-layer neurons are connected to all first layer neurons, resulting in wj 2∈R10×1. All biases bare scalar. The firstlayer neurons have a sigmoid activation function to cover the nonlinear NMPC control with h(σi) = 1 1 + exp(−σi).(1) The second layer neurons have a linear activation function. The output feedback is closed when the NARX is applied in the closed control loop, such that the current output of the network u(k)becomes a network input (u(k−1)) in the next time step. Since the control actions u(k)and u(k−1) are already known during training (from the NMPC), the NARX feedback loop is opened and u(k)and u(k−1) are used as both, target values and features, respectively, so that the NARX becomes an easier to train feedforward network (see [19]). TABLE I SPECIFIC NARX ARCHITECTURE. size delays connection activation inputs |r|+|y|+|u|- - - layer 1 10 1 full sigmoid (1) layer 2 |u|- full linear outputs |u|- - - A. Nonlinear model predictive controller We use a generic NMPC controller provided by the Matlab Model Predictive Control Toolbox 1and solve the nonlinear optimization problem VN(x) = min y(·),u(·),ε Jy+J∆u+Jε(2a) subject to x(k+ 1) = f(x(k), u(k)), k = 0, ..., N −1(2b) y(k) = g(x(k)), k = 0, ..., N (2c) y−ε≤y(k)≤y+ε, k = 0, ..., N, (2d) ε≥0(2e) u≤u(k)≤u, k = 0, ..., N −1.(2f) for a current state x(0) in every time step, and apply the optimal input value u∗(0) to the plant. The cost functions are defined as Jy= N X k=1 kr(k)−y(k)k2 Q, Jε=kεk2 S, J∆u= N−1 X k=0 ku(k)−u(k−1)k2 R. Cost function Jypenalizes the deviation of output y(k)from the given reference value r(k)for time step k. The input rate u(k)−u(k−1), i.e., the change in the applied input value between two time steps, is penalized in J∆uto avoid fastchanging actuator use. The slack variable is introduced to ensure feasibility by softening the output constraints on y(k). The use of the slack variable εis penalized by the the cost function Jε. The underlying optimization problems were solved using the SQP algorithm as part of the Matlab function fmincon. An Extended Kalman Filter (EKF) was used within the Nonlinear Model Predictive Control Toolbox to reconstruct the state vector xfrom plant outputs yand controller outputs u. III. APPLICATION: EMBEDDED PROCESS CONTROL We chose a laboratory-scale hydraulic plant with combined pressure and flow rate control as a demonstrator for the outlined control concept. The plant is depicted in Fig. 3 and consists of a variable-speed centrifugal pump and a controllable discharge valve. Pump and valve are the actuators required to set a desired pressure and flow rate in an enclosed chamber. The plant inputs are u= (n, v)T 1https://de.mathworks.com/help/mpc/index.html
M ~ ~ pump motor inverter valve outlet inlet PLC PC with MATLAB OPC controlled chamber ESP32 analog interface Fig. 3. Laboratory scale hydraulic plant. for pump speed and valve opening, respectively. The plant outputs are y= (p, q)T, i.e. the measured pressure and flow rate. The plant configuration resembles a frequent real-world control task that can be found, e.g., in reverse osmosis plants for seawater desalination [20]. The actuators and sensors are connected to a PLC and a standard PC with Matlab / Simulink, which are used to implement the NMPC. The PC is explicitly not required to run the NARX-based control, but performs data logging for the evaluation. A. Plant model A plant model is required for the NMPC and consists of the component models for pump, valve, and sensors. The model is organized in a Hammerstein structure with pump and valve as static, nonlinear models and the sensors as linear dynamic models. 1) Pump model: We assume quasi-static conditions for the pump and chose a standard centrifugal pump model to compute the steady state discharge pressure ϕfrom the steady state flow rate ψand the rotational speed n(see [21]), and adapted the model to the current pump: ϕ=cϕ,1·ψ2+cϕ,2·ψ·˜n+cϕ,3·˜n1.7,(3) with parameters cϕ,1=−8.726·10−6,cϕ,2=−1.691·10−6, cϕ,3= 201.906 ·10−6and a rotational speed scaling ˜n= 0.6·n+ 40. 2) Valve model: As for the pump we assume quasi-static conditions also for the valve and use the generic orifice equation (see [22] p. 18 ff) ψ=pϕ·cv(v),(4) with the valve coefficient cv(v)implemented as piecewise linear function. 3) Sensor models: The pressure and flow rate sensors are represented as linear dynamic, discrete-time state space models with state vector x(k) = xp(k)xq,1(k)xq,2(k)T and a sampling time TS= 0.5s. TSwas chosen relatively large to meet the cycle time of the online NMPC evaluation (see Fig. 9 in Sec. IV). The plant model x(k+ 1) = f(x(k), u(k)) results with f(x(k), u(k)) = ap01×2 02×1Aq·x(k) + I2×2 01×2·ϑ(k), (5) Fig. 4. The full operating range of the plant is inside of the solid black border, while we restrict the admissible range for pand qwithin the highlighted gray area by the controller constraints (6). The blue trajectory depicts data of three open-loop step responses, where blue dots correspond to sampled time points and lines are added for convenience. with 0and Idenoting zero and identity matrices, respectively, and ap= 0.6065, Aq=0.6065 0 0.3033 0.6065. The vector ϑ(k) = ϕ(k)ψ(k)Tinvokes the nonlinear models (3) and (4), and therefore connects (5) to the plant inputs u(k). The measured plant outputs are y(k) = p(k)q(k)Twith p(k) = xp(k)and q(k) = xq,1(k). The plant model (5) is used in the NMPC for closed-loop control of the real plant and for training data generation as described in Sec. II-A. Fig. 4 depicts the full operating range of the plant, bounded with a solid line, over the input range n, v ∈ [0,100]%. We deliberately limit the admissable operating range for the plant outputs to yand ¯y, highlighted by the rectangular area in Figure 4, with 0.25 30 ≤y(k)≤0.35 50 .(6) The lower and upper bounds on ycorrespond to the NMPC constraints (2d). The blue trajectory depicts the open loop operation of the plant with data from three step responses, from u= 47,45Tat point (1) to u=67,45Tat point (2) to u=47,54Tat point (3), to illustrate the open-loop characteristics. B. NMPC tuning and training data generation The NARX training data was generated in an extensive simulation study performed with the NMPC introduced in Section II-A running against the plant model (5) over a 5 hour interval, resulting in 36000 data points. The references rfor the training data generation were taken randomly from a set r∈[rmin, rmax]with rmin = 0.22,20Tand rmax =0.42,70T, which is also depicted in Fig. 6. Since we aimed at learning constraint satisfaction as well, reference values outside of the output constraints
were internally limited to the corresponding bound during the simulation study, to not include the effect of constraint softening in the training data. Therefore a modified reference value ˜r(k)was forwarded to the NMPC such that ˜r(k) = min(r(k), y), ˜r(k) = max(r(k), y), holds. The training data consists of data sets rand y, thus contains the original, unmodified reference values to train references outside the bounds and learn constraint satisfaction by the network. Note that the distinction between rand ˜ris not possible when using the control error e=r−y as feature input instead of independent data sets rand y, since econtains relative information only. The different combinations of references, as well as the distinction between control of the plant model and the real plant, are summarized in Fig. 5. The input values, i.e., the pump speed and the valve opening, are both scaled to a range from 0to 100%, which thus define the lower and upper bounds on u(k)as in (2f). Control and prediction horizons were both chosen to N= 5. The weighting matrices on outputs, inputs, and slack variables have the obvious dimensions. They were tuned such that the NMPC leads to the desired control performance and read Q=1040 0 1, R =0.25 0 0 0.25, S = 105. The EKF was implemented using a Gaussian white noise covariance matrix for the process noise as well as for the measurement noise QEKF = 10 ·I3×3, REKF = 0.01 ·I2×2. C. NARX training and embedded implementation We implemented the NARX structure in Keras and applied the non-normalized training data generated by the NMPC as outlined in Sec. III-B. Training was performed with the simulation EKF NMPC plant model EKF NMPC NARX real plant real plant evaluation NARX training ∈ , ∈ min , max ∈ min , max Fig. 5. The training data for the NARX was generated from a simulation run of the NMPC (+EKF), where the references were internally kept within the limits ˜r∈[y, ¯y], while a larger reference set r∈[rmin, rmax]was used as feature input to the NARX, as well as for the evaluation. Fig. 6. Left: Training data with steady-state plant outputs yand references r as a result of the NMPC control. Right: Performance throughout the training. NADAM solver, a batch size of 32 and the mean squared error (mse) as loss metric. Figure 6 depicts the training data set in the left diagram, where only the steady state values from the 5 hour data set are plotted. The NMPC controls the plant within the constrained operating area, while the references rare deliberately chosen from a larger, random set, for an efficient training of constraint satisfaction. The right diagram in Fig. 6 depicts the training progress with the mse as loss measure over the training epochs. We stopped the training after 200 epochs since the performance stalled at an mse = 0.6045. We exported the trained network parameters from Keras into a Python script and ran Micropython on the ESP32 with a cycle time of TSto evaluate the network as controller. IV. RESULTS We implemented the NMPC in Matlab / Simulink and used an OPC connection to the PLC to exchange measurements and control actions with the plant. The sampling times for the NMPC and NARX evaluation equal the sampling time of the plant model TS. The NARX-based controller on the ESP32 microcontroller was connected via an analog interface to Matlab / Simulink and from there to the PLC (cf. Fig. 3), which only served as an interface to the real plant. This allowed us to use Matlab / Simulink for data logging for both controllers. A. Time series results Figure 7 depicts the time series results for both controllers at the real plant. The upper two diagrams depict the plant outputs pand q, respectively, plotted together with the respective reference values rand the bounds from (6). The lower two diagrams depict the plant inputs nand v. As during the training, we again deliberately chose references rthat violate the bounds from (6), to demonstrate the constraint satisfaction ability of the neural network. It is obvious from the upper two diagrams in Fig. 7 that both, the NMPC and the NARX, respect the constraints. However, the NMPC allows some offset violation due to the slack variable (see Sec. II-A). Since the NMPC used only admissible references ˜rduring training data generation (to
Fig. 7. Time series results of a 25 minute closed-loop control with NMPC and NARX at the real plant. The references are rpand rqfor pressure and flow rate, respectively. not trigger constraint softening, cf. Fig. 5), the NARX control actually respects the limits more strictly than the NMPC control. The general control quality is sufficient for both controllers. The NARX controller shows a slightly higher steady-state offset. Table II summarizes the quantitative control quality. We use the mean squared error mse := K−1X K (r(k)−y(k)) ◦(r(k)−y(k)) , where ◦denotes the elementwise product, over K= 1600 samples between t= 600s and t= 1400s (cf. Fig. 7), i.e. excluding the (deliberate) references outside the constraints. Table II also lists the maximal absolute constraint violation max(|(y−¯y)|,∀y > ¯y, |(y−y)|,∀y < y). The numbers in parenthesis are relative to the maximal plant operating values max(p)=0.5bar and max(q) = 80 l/min. The quantitative measures underline the observations from the time series results in Fig. 7, namely that the reference tracking is comparable between both control variants, but the constraint satisfaction is better with the NARX-based control. Figure 8 depicts a time series of a disturbance rejection test, where we confronted both controllers with unmeasured disturbances on the two plant inputs. The lighter gray bars TABLE II QUANTITATIVE CONTROL QUALITY MEASURES. variant variable mse max. violation NARX p4.96 ·10−4(0.1%) 0.0053 (1.1%) q16.55 (16.5%) 0.684 (0.85%) NMPC p2.32 ·10−4(0.05%) 0.03 (6.4%) q16.12 (20.16%) 3.3 (4.11%) mark timespans where nis disturbed by +10% for the first bar, followed by −10% for the second bar. The following darker gray bar marks a timespan where the valve vwas disturbed by +10%. Both controllers react in a desired way by reducing the disturbance, however, both fail to reject it completely. For the NMPC, the plant-model mismatch increases significantly due to the unmeasured disturbance, leading to false predictions of the trajectory and an under-compensation of the disturbance. The constraint violation during the second disturbance (for t∈[90,120]) leads to a more rigorous control action of the NMPC. For the NARX-based control, the disturbance rejection is comparable to the NMPC results. This is remarkable because it demonstrates that the neural network was able to adapt the NMPC controller behavior instead of just memorizing control actions, since the disturbances were not part of the training data. Fig. 8. Time series of the disturbance rejection test for both control variants. During the timespans marked by a light gray color bar the pump speed n was disturbed by +10% and −10% (in that order), and during the timespan marked by a darker gray bar the valve opening vwas disturbed by +10%.
Fig. 9. Profiling results recorded during the test run from Fig. 7. A boxplot is included for both measurements, where the boxes cover 75% of the data, the whiskers cover 99.3%, and the magenta crosses mark the remaining data points as outliers. B. Profiling The two objectives of the NARX-based control are to reduce the computational effort and to provide a strict realtime guarantee. We quantified these with a profiling during the controller evaluation from Fig. 7. The profiling results are depicted in Fig. 9. The NPMC was evaluated on a standard PC with a 3.20 GHz CPU, while the NARX was evaluated on the ESP32 set to 240 MHz. We stress the evaluation times are not scaled by, or weigthed with, the CPU frequencies. Despite its drastically lower CPU frequency, the NARX evaluation on the ESP32 outperforms the NMPC evaluation on the PC, reaching nearly a factor of 10. Both variants stay within the sampling time of TS= 0.5s, however, the NMPC shows higher fluctuations in the evaluation time which suggests to keep a greater margin to the cycle time to not fail in updating the control action. The NARX-based control runs with very constant evaluation time, making it favorable for real-time conditions. V. CONCLUSIONS We outlined the mimicking of a nonlinear model predictive controller by a small-scale autoregressive neural network with the aim of efficient real-time implementation on a low level hardware. The results demonstrate a comparable control quality between both, the computationally expensive NMPC and the computationally much cheaper NARX-based controller, when applied to a hydraulic plant. Remarkably, the small-scale neural network was able to learn the NMPC behavior instead of memorizing the control actions, as demonstrated by the rejection of unmeasured disturbances. Moreover, the NARX-based controller inherits a certain degree of constraint satisfaction, again comparable to the NMPC. REFERENCES [1] D. Q. Mayne, J. B. Rawlings, C. V. Rao, and P. O. Scokaert, “Constrained model predictive control: Stability and optimality,” Automatica, vol. 36, no. 6, pp. 789–814, 2000. [2] L. Gr¨ une and J. Pannek, Nonlinear model predictive control. Springer, 2017. [3] M. Diehl, R. Findeisen, F. Allg¨ ower, H. Bock, and J. Schl¨ oder, “Nominal stability of real-time iteration scheme for nonlinear model predictive control,” IEE Proceedings-Control Theory and Applications, vol. 152, no. 3, pp. 296–308, 2005. [4] T. A. Johansen, “Approximate explicit receding horizon control of constrained nonlinear systems,” Automatica, vol. 40, pp. 293–300, 2004. [5] R. Dyrska and M. M¨ onnigmann, “Accelerating nonlinear model predictive control by constraint removal,” in Proc. of the 7th IFAC Conference on Nonlinear Model Predictive Control, 2021, pp. 278– 283. [6] R. Quirynen, K. Berntorp, and S. Di Cairano, “Embedded optimization algorithms for steering in autonomous vehicles based on nonlinear model predictive control,” in 2018 Annual American Control Conference (ACC). IEEE, 2018, pp. 3251–3256. [7] S. Adhau, S. Patil, D. Ingole, and D. Sonawane, “Implementation and analysis of nonlinear model predictive controller on embedded systems for real-time applications,” in 2019 18th European Control Conference (ECC). IEEE, 2019, pp. 3359–3364. [8] A. Raha, A. Chakrabarty, V. Raghunathan, and G. T. Buzzard, “Embedding approximate nonlinear model predictive control at ultrahigh speed and extremely low power,” IEEE Transactions on Control Systems Technology, vol. 28, no. 3, pp. 1092–1099, 2019. [9] L. Hewing, K. P. Wabersich, M. Menner, and M. N. Zeilinger, “Learning-based model predictive control: Toward safe learning in control,” Annual Review of Control, Robotics, and Autonomous Systems, vol. 3, pp. 269–296, 2020. [10] K. P. Wabersich and M. N. Zeilinger, “Linear model predictive safety certification for learning-based control,” in 2018 IEEE Conference on Decision and Control (CDC). IEEE, 2018, pp. 7130–7135. [11] M. Klauˇ co, M. Kal´ uz, and M. Kvasnica, “Machine learning-based warm starting of active set methods in embedded model predictive control,” Engineering Applications of Artificial Intelligence, vol. 77, pp. 1–8, 2019. [12] A. Domahidi, M. N. Zeilinger, M. Morari, and C. N. Jones, “Learning a feasible and stabilizing explicit model predictive control law by robust optimization,” in 2011 50th IEEE conference on decision and control and European control conference. IEEE, 2011, pp. 513–519. [13] S. Chen, K. Saulnier, N. Atanasov, D. D. Lee, V. Kumar, G. J. Pappas, and M. Morari, “Approximating explicit model predictive control using constrained neural networks,” in 2018 Annual American control conference (ACC). IEEE, 2018, pp. 1520–1527. [14] K. Kiˇ s, M. Klauˇ co, and M. Kvasnica, “Explicit MPC in the form of sparse neural networks,” in 2021 23rd International Conference on Process Control (PC). IEEE, 2021, pp. 163–168. [15] Y. L¨ ohr, M. M¨ onnigmann, M. Klauˇ co, and M. Kal´ uz, “Mimicking predictive control with neural networks in domestic heating systems,” in 2019 22nd International Conference on Process Control (PC19). IEEE, 2019, pp. 19–24. [16] Y. L¨ ohr, M. Klauˇ co, M. Fikar, and M. M¨ onnigmann, “Machine learning assisted solutions of mixed integer mpc on embedded platforms,” IFAC-PapersOnLine, vol. 53, no. 2, pp. 5195–5200, 2020. [17] J. Drgoˇ na, D. Picard, M. Kvasnica, and L. Helsen, “Approximate model predictive building control via machine learning,” Applied Energy, vol. 218, pp. 199–216, 2018. [18] S. Lucia, D. Navarro, B. Karg, H. Sarnago, and O. Lucia, “Deep learning-based model predictive control for resonant power converters,” IEEE Transactions on Industrial Informatics, vol. 17, no. 1, pp. 409–420, 2021. [19] K. S. Narendra and K. Parthasarathy, “Learning automata approach to hierarchical multiobjective analysis,” IEEE Transactions on Systems, Man, and Cybernetics, vol. 21, no. 1, pp. 263–272, 1991. [20] M. Mohammed and M. E. A. Mohamed, “Control of reverse osmosis process at a brackish water desalination station,” in International Conference on Electrical Systems & Automation. Springer, 2022, pp. 143–154. [21] S. Leonow, F. Wollenhaupt, and M. M¨ onnigmann, “Combined flow and pressure control for industrial pumps with simple adaptive MPC,” in Proceedings of the 21st International Conference on Process Control (PC17), High Tatras, 2017, pp. 315–320. [22] J. F. G¨ ulich, Centrifugal pumps. Springer, 2010, vol. 2.