16th IFAC Symposium on Automation in Mining, Mineral and Metal Processing August 25-28, 2013. San Diego, California, USA
Parameter estimation for a flotation process tracking simulator Janne Pietilä ∗ Jani Kaartinen ∗∗ Ari-Matti Reinsalo ∗∗∗ ∗
Aalto University Department of Automation and Systems Technology, Helsinki, Finland (e-mail:
[email protected]). ∗∗ (e-mail:
[email protected]) ∗∗∗ (e-mail:
[email protected])
Abstract: A parameter estimation method for an online flotation process simulator is described. The applications of an online process simulator include soft sensor implementations, process trajectory predictions and advisory feedback to the operator, which have potential to improve the process efficiency and minimize the detrimental effect of disturbances on the process or the environment. The online simulator uses a detailed and dynamic model of an actual industrial flotation process, and therefore accurately corresponds to the process phenomena present at the plant. Parameter estimation is required for the flotation kinetics in order for the simulator to adapt to changes in the process conditions. A relatively generic parameter estimation algorithm is developed and tested with a dual simulator setup. The particular requirements and limitations of adapting an online simulator are discussed, and modifications to well-known estimation algorithms are presented as a possible method of meeting the requirements and overcoming the limitations. The results show that online parameter estimation and simulation model tracking is possible with the chosen method, and point out areas of further development for application when the simulator is used alongside the real process. Keywords: estimation algorithms, adaptation, process simulators 1. INTRODUCTION Mineral flotation processes are often vast and complex, containing numerous subprocesses and units. Even with modern instrumentation and automation, the process state may not be completely known or the behaviour may not be completely predictable. As the trend is towards ever larger processes operated by fewer personnel, more tools are needed to handle the process operation and to obtain more information about the process state. Simulators have been an important aid in process control for quite a while now. Recently there has been interest to utilize dynamic simulators in place of steady-state process simulation. Also, these simulators are no longer restricted to process planning and design phases, but are seeing applications also in online operation. An overview of simulator integration to a process, along with a summary of the state of the art is given in Šindelář and Novák (2011). An online tracking simulator implements a dynamic model of the process and runs alongside the actual process, mirroring its behaviour. Essential to a tracking simulator is the ability to adapt to changes in the process conditions, so that the simulator outputs continue to match those of the actual process over an extended period of time. This would be easy to achieve were all the process conditions readily and accurately measurable. However, usually changes in the process parameters are only detected indirectly through measurements of the process outputs, or by laboratory analysis after a significant delay. The challenge is to find a method for adapting a complex simulation model to 978-3-902823-42-7/2013 © IFAC
122
changes in the process with a preferably simple algorithm, and information about process conditions gathered only through indirect online measurements. What is needed, then, is a method to estimate the unmeasurable parameters of the process and use this estimate in the dynamic online process simulator to ensure the continued precision in the simulator predictions. Of primary concern is that the outputs of the simulator match those of the process, making the simulator behave outwardly identically to the process. Second, it is also desirable that the internal parameters of the simulator remain plausible, assuming these parameters have some physical interpretation relating to the process. If this goal is achieved, more relevant information about the process behaviour can be obtained in real time as the process is operated. It should be noted, however, that when the simulator is operated alongside a real process, reference information about the internal process parameters is rarely available. A more easily measurable goal in this case would be to have a stable estimate of the process parameters, with gradual and smooth responses to changes in process conditions. Once an estimation method is implemented, the applications of an online tracking simulator include soft-sensor implementations and operator support. A comprehensive review of dynamic process simulators and their applications, including a detailed discussion of dynamic operator training simulators, can be found in Ahmad et al. (2010). Safety improvement and operational enhancement are highlighted as particular advantages of simulator integration. Soft sensor techniques and applications are re10.3182/20130825-4-US-2038.00048
IFAC MMM 2013 August 25-28, 2013. San Diego, USA
viewed in e.g. Kadlec et al. (2009), and a special emphasis on soft sensor adaptation is given in Kadlec et al. (2011). A comprehensive discussion of soft sensing in mineral processing can be found in Sbarbaro and Villar (2010, Chap. 4). 2. METHODS 2.1 Tracking simulators Recent advances in online process tracking simulators include the studies reported by Nakaya and Li (2013); Nakaya et al. (2011); Ishimaru et al. (2010); Nakaya et al. (2008). They utilize a very detailed dynamical model of the process to implement a simulation model that is constantly updated to match the real process. The simulator parameters are adjusted proportionally to the prediction error. However, in later studies this was found to be insufficient and the parameter identification is augmented with a dynamic data reconciliation method. The application examples are industrial processes, including a polymer electrolyte membrane fuel cell (Nakaya et al., 2008), a depropanizer process (Ishimaru et al., 2010) and a large petrochemical process (Nakaya and Li, 2013). The adaptation problem of tracking simulators has also been recently approched with traditional control theoretic methods by Friman and Airikka (2012). The problem with adapting a dynamic simulation model to the process is that the current process measurements alone do not determine the process state, but rather the whole process history must be considered. Friman and Airikka state that the traditional method of updating the parameters based on the prediction error alone is slow and difficult to tune for the simulator. They propose a method which uses a PI controller to adjust the simulator parameters, coupled with autotuning for the PI controller parameters. The proposed method is demonstrated on a simulation example and shown to function also with data from an industrial power plant. Besides parameter tracking, the estimation of the process state is also important in tracking simulators. After all, if the process and the simulator states differ, the output values will also be different despite identical parameter values (or, alternatively, the simulator parameters will diverge as the simulated outputs are driven to match the output values of the process). A review of some popular online state estimation methods is given by Schei (2008). Parameter estimation is also discussed, and in the context Schei remarks that: Usually it can not be assumed that the process excitations fulfil certain persistency requirements in order to ensure convergence of parameter estimates. Hence, the parameter vector θ should normally be a set of parameters which is identifiable from stationary data, and which do not require any particular excitations in order to obtain convergence. Typical choices are heat and mass transfer coefficients, kinetic parameters, etc. By carefully selecting a number of parameters nθ equal to the number of output measurements ny , zero steady-state deviations
in all predicted output measurements can normally be achieved. Some general points about process modeling and simulation, and parameter estimation are examined next. Assume that the dynamic behaviour of the actual process follows x(t) ˙ = f (x(t), u(t), w(t), θ(t)) (1) y(t) = g (x(t), u(t), w(t), θ(t)) where x(t) represents the vector-valued system state, u(t) is the similarly vector-valued input to the system, w(t) represents the unknown disturbances entering the system, and the system parameters are denoted by θ(t). This state-space equation govering the system behaviour is generally not known. Rather, by various means of system identification, it can be approximated by a dynamic model, e.g. in the form given by ′ x˙ (t) = f ′ (x′ (t), u(t), θest (t)) (2) yest (t) = g′ (x′ (t), u(t), θest (t)) The objective of system identification and parameter estimation in tracking simulators is to develop a model that approximates the behaviour of system outputs y(t) sufficiently accurately with the estimated outputs yest (t). The application of the model determines the sufficient accuracy. In many cases the dimension of the model state x′ (t) is less than the system state x(t). Similarly, only a subset of the process outputs y(t) may be of interest in the model application, and thus the dimension of the model outputs yest (t) may be lower as well. Of particular note are also the disturbances w(t): as they are unknown, they naturally can not be included in the model (2). (Were some of the disturbances known and measurable, they could be included in the model inputs u′ (t) even if the disturbances themselves could not be controlled.) Given the differences between the actual process and the simulation model, it follows that the parameters of the model θest (t) do not necessarily coincide with the parameters of the actual process θ(t). Although this is not strictly necessary for the model outputs yest (t) to follow the process outputs y(t) (or the application-specific subset of those), it is desirable as the simulator could then be utilized also in the soft sensor role, providing information about the time-varying system parameters θ(t) which would otherwise be time-consuming, difficult or even impossible to acquire. Whether the model parameters can be made to coincide with the actual parameters of the process depends strongly on the chosen model structure and foundation; a sensible assumption is that a model that is based on physical phenomena or phenomenological relationships is more likely to have the parameters coincide than, say, a regression or time-series model obtained from process data. Other requirements for the model in a tracking simulator include its applicability over the entire range of operating conditions in the process. There are various ways to achieve this. A global model that handles the variations in operating conditions by only changing the values of the model parameters demands a lot of knowledge about the actual process in the model development stage. Simple loworder models are less demanding, but only cover a small area of the operating conditions unless the process itself is sufficiently simple and linear. However, multiple linear
123
IFAC MMM 2013 August 25-28, 2013. San Diego, USA
models can be used in parallel and switched between, so that each model covers a different operating point.
θĞƐƚ͘;ƚнϭͿ
2.2 Dynamic HSC simulator
,^^ŝŵʹŵŽĚĞů
This research utilizes the Outotec HSC simulator as the simulation model. The HSC simulator is used to model the copper flotation process (Figure 1) of Pyhäsalmi Mine in central Finland. In the HSC simulator, unit processes – mainly flotation cells in this case – are modelled with dynamic mass balance calculations and first order flotation kinetics, with three floatability groups (fast, slow, and nonfloating) per mineral. The overall flotation kinetics are in addition affected by the bubble surface area flux (Gorain et al., 1997), and an entrainment factor (Savassi, 2005) is included in the froth recovery. A brief description of the simulation model can be found in Lamberg (2010). The input and output streams of the process units are connected to each other in an identical arrangement to the process flowsheet to achieve a faithful reproduction of the actual flotation process. The modelling approach used in HSC makes it easy to get the process engineers at the plant involved with the modelling effort, as the simulation model is represented in a similar flowsheet as the actual process is in the control room displays of the automation system. Also the modelled phenomena are ones that the process engineer is familiar with. More details about the HSC simulator can be found in an application example reported by Roine et al. (2011) and references therein. The use of flotation kinetic models means that the model parameters have a physical interpretation or are closely related to e.g. the ore properties in the actual process. It also allows the calibration of the model parameters by laboratory experiments with samples taken from the actual process. As the simulation model corresponds to the process representation in the automation system of the plant, the model is also readily integrated with the measurements provided by the automation system – an important feature for an online tracking simulator. However, the model structure in the HSC simulator is very different from the mathematical models on which the theory of system identification and parameter estimation is based. This makes the development of a suitable parameter estimation algorithm for the tracking simulator a challenging task. Also, due to the complexity of the model, the simulation speed is not very high which limits the choice of adaptation techniques. 2.3 Parameter estimation and simulator adaptation The general structure and information flow of the parameter estimation is presented in Figure 2. The Master process is either the actual, physical process or a HSC simulator emulating the actual process. The Slave process is always the simulator which follows the behaviour of the Master by changing the adaptable parameter values. The Parameter Estimator reads any available measurements from both the Master and the Slave processes, and implements the adaptation algorithm. The signal w(t) represents any unknown disturbances entering the actual process. A thorough discussion of implementing this online tracking simulator can be found in Kaartinen et al. (2013). 124
^>s
LJĞƐƚ͘;ƚͿ
θĞƐƚ͘;ƚͿ DĂƚůĂď
WZDdZ ^d/DdKZ
LJ;ƚͿ Ƶ;ƚͿ ǁ;ƚͿ
,^^ŝŵŵŽĚĞů ŽƌĐƚƵĂůƉƌŽĐĞƐƐ D^dZ
θ;ƚͿ
Fig. 2. Parameter adaptation Table 1. Output variable and kinetic parameter pairs for adaptation. Output variable
Kinetic parameter
Cu recovery (%) Zn recovery (%) S recovery (%) Pb recover (%) Cu grade in concentrate (%)
Chalcopyrite (Ccp) k-fast Sphalerite (Sp) k-fast Pyrite (Py) k-fast Galena (Gn) k-fast Gangue (Gang) k-fast
The adaptation parameters are a set of kinetic coefficients of flotation, listed in Table 1. Also listed are the corresponding model output variables for each parameter. For simplicity, only the fast floating of each mineral was considered. Correspondingly, only this one kinetic coefficient per mineral is used in the entire process for all size fractions, although the simulator allows coefficients to be defined separately for individual size fractions (and even for individual flotation cells). The principle is to use the observed difference in the output values of the process and the simulator to manipulate the corresponding adaptation parameter, and thus make the simulator track the process. Reference values for the kinetic coefficients have been determined experimentally with laboratory analysis; these values are used internally in the HSC simulator, while the adaptation is carried out by scaling the reference values. The prediction error between each output of the actual process (process simulator) y(t) and the simulation model (model simulator) yest (t) is calculated according to e(t) = y(t) − yest (t).
(3)
IFAC MMM 2013 August 25-28, 2013. San Diego, USA
Scavenging
Roughing feed
tailings Middlings flotation
Cleaning
concentrate Fig. 1. Schematic of the copper flotation circuit at Pyhäsalmi Mine Table 2. Weighting factors in the error cost function. Key: + for a positive value, - for a negative one. e(t)
ˆ e (t) D
ke
kd
+ + -
+ + -
1 1 1 1
1 0.5 1 0.5
ˆ e (t) of the rate of change (time Additionally, an estimate D derivative) in the output error is calculated by fitting a first order polynomial to Ne = 10 most recent observations of the error. The estimated derivative is then used to calculate the error cost function as a weighted sum: ˆ e (t). ǫ(t) = ke (t)e(t) + kd (t)D (4) The modified error function takes into account the rate of change in the model outputs. A similar approach, but with PID control as the adaptation algorithm, has been proposed and demonstrated by Friman and Airikka (2012). The weighting factors ke (t) and kd (t) have different values depending on whether the estimation error is decreasing or increasing, according to Table 2. The parameter update step is based on well-known recursive estimation techniques (see e.g. Ljung (1998)). The current estimate is recursively updated proportionally to the residual ǫ(t): θest (t + 1) = θest + Γ(t)ǫ(t) (5) The parameter update gain Γ(t) is usually chosen based on the model structure and equations. For the ith output it is determined by wi (t)Pi ϕi (t) (6) Γi (t) = wi (t)ϕi (t)P ϕi (t) + λ(t) ϕi (t) denotes the predictor for the ith output variable. Usually the predictor would be derived by linearizing the 125
model equation (2). However, since the HSC simulator uses local unit models which are not readily expressable in a state-space representation (2), the predictors were instead selected based on a simplified model of flotation kinetics. The recovery of a mineral to tailings of a flotation cell is determined by ctail Ftail Qtail = (7) Rtail = Qf eed cf eed Ff eed where Q(∗) denotes the mass flow of the mineral, c(∗) the concentration, and F(∗) the overall mass flow in stream (∗) (see e.g. Wills and Napier-Munn (2006)). As the flotation kinetics are directly affecting the mineral grades in the tailings stream, these grades were selected as predictors for mineral recoveries: ϕi (t) = ctail,i (t) (8) Pi is used to denote the inverse of the ith predictor (co)variance. It is updated each step according to Pi (t + 1) = 1/λ(t) (Pi (t) − Γi (t)ϕi (t)Pi (t)) (9) The update gain (6) includes the forgetting factor λ(t), which enables more weight being placed on more recent observations of the residual and thus allows the estimation algorithm to better adapt to changes in the system parameters. Here the forgetting factor is time-variant, and decreases if the estimation error increases: λ(t) = 1 − 1/Ω (1 − ϕi (t)Γi (t)) e2i (t), (10) where Ω = σi2 /(1 − λ0 ) and λ0 = 0.995 is the nominal value for the forgetting factor. Noise variance for the ith residual is represented by σi2 , and is estimated directly from the measurement. If the measurements are suspected to contain outliers, a more robust estimate can be obtained from the median absolute deviation (MAD) (Ljung, 1998): median (|ei − median(ei )|) , (11) σi = 0.7
Model parameter Ccp k−fast
Scavanger tailings Cu recovery (%)
IFAC MMM 2013 August 25-28, 2013. San Diego, USA
fast and there is practically no oscillation present. These results could indicate that the algorithm is sensitive to the operating point of the process, which could be solved with some sort of diagnostic online tuning of the algorithm, for example.
Model Process
4 3.9 3.8
The example shows how the chosen adaptation algorithm performs when a step change in one of the Master process parameters takes place. The primary objective of tracking the process output is reasonably well achieved. More development is needed for the Slave process parameters to accurately follow the Master parameters, as now there is a slight bias in the estimated parameter value. This bias seems to be the result of a discrepancy in the values of the state variables between the two simulation models, brought about by the change in a flotation parameter value. As the Slave simulator adapts to the new conditions, this bias will decrease, but currently the convergence is quite slow. The oscillation in both the simulator outputs and the estimated parameter values is also an issue meriting further attention.
3.7 3.6 3.5
0
5
1.2
10 15 20 Simulation time t [h]
25
Model Process
1.1 1 0.9 0.8 0
5
10 15 20 Simulation time t [h]
25
4. CONCLUSIONS AND DISCUSSION
Fig. 3. Results of the adaptation example with changes in the relative floatability of chalcopyrite. where ei represents the time series of the ith residual. Finally, to ensure that the effect of outliers on the residual is minimized, an IGGI estimator (Bao et al., as cited in Chao et al. (2008)) is used with (6): |ǫi (t)| ≤ k1 σi wi = 1 (12) wi = k1 / |ǫi (t)| k1 σi < |ǫi (t)| ≤ k2 σi wi = 0 |ǫi (t)| ≥ k2 σi 3. RESULTS The functioning of the parameter estimation algorithm is here illustrated with an example. The chalcopyrite fast flotation factor, Ccp k-fast, is altered in the Master simulation, and the behaviour of the respective parameter in the Slave simulator is then observed, along with the corresponding output variables. The Master and Slave simulators start from identical states, setpoint values and feed properties. At time t = 11 h53 min, the Ccp k-fast parameter value in the Master simulation is increased by 10 %, from 1 to 1.1 relative value (Figure 3). This causes more copper to float to concentrate, decreasing copper recovery to tailings. As the copper recoveries between the Master and the Slave simulators begin to diverge, the adaptation algorithm reacts by increasing the Ccp k-fast parameter in the Slave simulation. Initially the reaction is quite fast, leading to overshoot and oscillation. The parameter value in the Slave simulator stabilizes at t = 12 h 32 min . Similar step changes take place at t = 14 h 32 min when the Master simulator Ccp k-fast parameter returns to its original value, and again at t = 24 h 34 min when there is a 10 % decrease in the chalcopyrite floatability parameter. The reaction of the adaptation algorithm to the second step change is almost ideal, as the response is very 126
A parameter adaptation method and an estimation algorithm for a flotation process tracking simulator has been presented. The dynamic simulation model and the HSC simulator accurately reproduce the structure of an actual industrial mineral concentration process, and the behaviour of the output variables is reasonably well replicated. The tracking simulator is expected to run online with the process, and form the basis for advanced soft sensor and operator advisory systems. The presented estimation and adaptation method is able to track changes in the process so that the outputs of the Slave simulator quite accurately follow those of the Master. More development is still needed to improve the accuracy of the parameter estimate and to eliminate the oscillation observed immediately after step changes in the parameters. The presented estimation and adaptation method is fairly generic in the sense that it does not utilize information about the structure of the process plant in question. However, the algorithm is tuned for the process by various tuning constants, for which the values have so far been selected manually in an ad-hoc fashion. In the future, sophisticated tuning algorithms and less complex estimation methods will be seeked to decrease the reliance on manual tuning parameters. So far the parameter-output pairs (Table 1) have been considered independent of each other. However, simulations indicate that there are interdependencies, with one flotation parameter affecting several output variables. Future development of the estimation algorithm will take these crosscorrelations into account. Alternatives for the chosen adaptation method include a trajectory matching method, where the Slave simulator is iteratively run for a fixed time interval with different parameter values, the output trajectory being compared after each iteration with the corresponding trajectory of the actual process. The parameter values could then be selected with a range of optimization algorithms. The sim-
IFAC MMM 2013 August 25-28, 2013. San Diego, USA
ulation speed is the primary limitation of this approach; it would need to be hundreds of times faster than realtime in order to allow enough iterations to take place. In contrast, with the proposed algorithm a realtime simulation speed is sufficient, and there is in fact no advantage in a higher speed. A downside of the presented algorithm is the number of tuning parameters, for which it is not simple to determine appropriate values. An adaptation approach based on PID control will be tested in the future, the advantage being that numerous tuning methods for PID control are known and can be tested. The adaptation performance was here tested with step changes in the process parameters. However, it is not very likely that such abrupt changes occur in the real process. It would be useful for the future development of the adaptation method to implement tests with e.g. gradual drift in the kinetic flotation parameters. So far, development of the tracking simulator has been done in a laboratory environment with a process simulator performing the role of the Master process. In the future, the tracking simulator will be implemented alongside the real process plant. This brings a few changes to the overall setup; for example, sample times in the real process are longer, and not equal between different measurements. ACKNOWLEDGEMENTS The work reported in this paper is done in the SOREX project, which is a part of the Green Mining program funded by TEKES –– the Finnish Funding Agency for Technology and Innovation. The authors would like to thank Outotec Oyj and Pyhäsalmi Mine Oy for the fruitful collaboration. REFERENCES Ahmad, A., Low, E., and Shukor, S.A. (2010). Safety improvement and operational enhancement via dynamic process simulator: A review. Chemical Product and Process Modeling, 5(1). DOI:10.2202/1934-2659.1502. Chao, Z., Hua-sheng, H., Wei-min, B., and Luo-ping, Z. (2008). Robust recursive estimation of auto-regressive updating model parameters for real-time flood forecasting. Journal of Hydrology, 349(3–4), 376–382. Friman, M. and Airikka, P. (2012). Tracking simulation based on PI controllers and autotuning. In Proceedings of the IFAC Conference on Advances in PID Control (PID’12). Gorain, B.K., Franzidis, J.P., and Manlapig, E.V. (1997). Studies on impeller type, impeller speed and air flow rate in an industrial scale flotation cell. Part 4: Effect of bubble surface area flux on flotation performance. Minerals Engineering, 10(4), 367–379. Ishimaru, S., Nakaya, M., and Ohtani, T. (2010). An application of tracking simulator to depropanizer process. In Proceedings of the SICE Annual Conference, 1486–1489. Kaartinen, J., Pietilä, J., Remes, A., and Torttila, S. (2013). Using a virtual flotation process to track a real flotation circuit. In Proceedings of the 15th IFAC Symposium on Control, Optimization and Automation in Mining, Mineral and Metal Processing. In print. 127
Kadlec, P., Gabrys, B., and Strandt, S. (2009). Datadriven soft sensors in the process industry. Computers and Chemical Engineering, 33(4), 795–814. Kadlec, P., Grbić, R., and Gabrys, B. (2011). Review of adaptation mechanisms for data-driven soft sensors. Computers and Chemical Engineering, 35(1), 1–24. Lamberg, P. (2010). Structure of a property based simulator for minerals and metallurgical industry. In Proceedings of SIMS 2010, the 51st Conference on Simulation and Modelling. Ljung, L. (1998). System Identification: Theory for the User. Pearson Education. Nakaya, M., Ikegaya, Y., Nakabayashi, A., Ootani, T., Chen, D., Li, X., Wang, D., and Ogata, M.G. (2011). Online process simulator with hybrid model of physical model and just-in-time model. In Proceedings of the 18th IFAC World Congress, 1640–1645. Nakaya, M. and Li, X. (2013). On-line tracking simulator with a hybrid of physical and just-in-time models. Journal of Process Control, 23(2), 171–178. Nakaya, M., Seki, T., Kawaguchi, K., Onoe, Y., and Ootani, T. (2008). Model parameter estimation by tracking simulator for the innovation of plant operation. In Proceedings of the 17th IFAC World Congress, 2168– 2173. Roine, T., Kaartinen, J., and Lamberg, P. (2011). Training simulator for flotation process operators. In Proceedings of the 18th IFAC World Congress, 12138–12143. Savassi, O.N. (2005). A compartment model for the mass transfer inside a conventional flotation cell. International Journal of Mineral Processing, 77(2), 65–79. Sbarbaro, D. and Villar, R. (eds.) (2010). Advanced control and supervision of mineral processing plants. Springer. Schei, T.S. (2008). On-line estimation for process control and optimization applications. Journal of Process Control, 18(9), 821–828. Šindelář, R. and Novák, P. (2011). Framework for simulation integration. In Proceedings of the 18th IFAC World Congress, 3569–3574. Wills, B.A. and Napier-Munn, T. (2006). Wills’ Mineral Processing Technology. Elsevier, 7th edition.