Multiple Property Cross Direction Control of Paper Machines

Cross direction (CD) control in sheet-forming process forms a challenging problem with high dimensions. Accounting the interactions between different properties and actuators, the dimensionality increases further and also computational issues arise. We present a multiple property controller feasible to be used especially with imaging measurements that provide high sampling frequency and therefore enable short control interval. The simulation results state the benefits of multiple property CD control over single property control and single property control using full feedforward compensation. The controller presented may also be tuned in automated manner and the results demonstrate the effect of tuning on input saturation.


Introduction
Continuous film and sheet-forming processes, such as paper manufacturing, set a challenging problem for process control.The properties to be controlled vary temporally along the production line (i.e. machine direction, MD) and spatially across the machine (i.e.cross direction, CD), where the overall control objective is to produce an as flat surface as possible.Typically, the control problem is separated into MD control and CD control, where the former can be considered as a single-input-single-output (SISO) control loop for setpoint following and disturbance rejection.However, the CD process consists of hundreds of measurement locations and a very large number of actuators distributed across the machine, hence forming a high dimension, non-square multi-input-multi-output (MIMO) system to be used for disturbance rejection in CD.Another problem comes from the actuator constraints, which need to be considered in control design.Figure 1 shows the schematic picture of a paper machine with different actuator positions and the mea-surement position.The sources of disturbances come from the processes before the headbox as well as different processing stages indicated in Figure 1.There are significant interactions between different properties and actuators in the paper machine (Backstrom and He, 2004).Hence actuations performed to one property may have unfavorable effects to another property.For instance, the slice lip actuators primarily used for controlling the basis weight variations will also affect on moisture and caliper, which are controlled with downstream actuators.Further, there might also be several actuators for single property, e.g.steam box and rewet shower profilers for moisture.
Interactions are recognized and taken into account in MD control, where usually multivariable model predictive controller (MPC) is a straightforward solution (Lang et al., 2001;Kuusisto et al., 2002;Nyuan, 2007).Novel methods to account for interactions also in CD control have been reported in 1990s (Hall, 1991;Heaven et al., 1994), but only recently the interactions have been properly handled with multiple property controllers (Backstrom and He, 2004;Backström et al., 2001;Haznedat and Arkun, 2002;Duncan and Heath, 2008;Fan, 2003;Shakespeare et al., 2003;Fu et al., 2006;Lahouaoula and Gheoghe, 2010;Duncan and Heath, 2010).The early papers use feedforward type compensation (Hall, 1991) or frequency separation (Heaven et al., 1994) to prevent competing control actions.Backström et al. (2001) implemented the CD-MPC to a linerboard machine with two properties and four actuator sets.Only one actuator set was treated as manipulated variable(s) in MPC while three other actuator sets were maintained under conventional control and treated as measurable disturbances to MPC.This approach used downsampled profiles with dimensions equal to highest actuator resolution.In a later paper (Backstrom and He, 2004) the measurements were mapped to a common resolution greater than three times the highest actuator resolution.The controller was extended to account for MD profile control, too.Both these implementations use an efficient QP solver to meet the online requirements, where the scan time sets the upper limit for computational time.Haznedat and Arkun (2002) introduced a CD-MPC where computational issues were solved by using a reduced order subspace in control calculations.The reduction of the dimensionality of the problem was done with Karhunen-Loeve transformation.Another paper (Fan, 2003) using model reduction techniques applied the CD-MPC in wavelet domain.As a result, the computational time was nearly halved from 32.55 seconds to 17.18 seconds when using a full model and reduced model for the system with three properties and four actuator sets.The performance degradation due to reduced model was quite modest, only few percents.Shakespeare et al. (2003) presented a high resolution robust optimal controller to maintain full control capability.Robustness against mapping errors was achieved with a set of penalty functions applied to steady-state optimal controller for multiple profiles and multiple actuator sets.The results presented were, however, from a paper mill with two actuator sets for single profile control.Duncan & Heath extended a robust constrained IMC CD controller to multi-array processes.The interactions were first (Duncan and Heath, 2008) accounted only in spatial domain, but later (Duncan and Heath, 2010) also the different dynamics of actuator sets were considered.Again, the example was for two actuator sets and single profile control.In (Fu et al., 2006), both high resolution profiles and dynamics were accounted for.The computational issues were attacked with sparse matrix structure and gradient type optimization.The example given was from a paper mill with three properties and three actuator sets.A very recent paper (Lahouaoula and Gheoghe, 2010) reports improvements of σ-values of 10%, 60%, 32%, 33% and 22% for dry weight, moisture, caliper, top coat weight and bottom coat weight, respectively, when changing the traditional CD control to multivariable multiple CD-MPC controllers using a combined feedback-feedforward control strategy.These numbers clearly show the potential and importance of multiple property control.
Another important aspect of CD control is its poor dynamic performance, which is limited by the long dead-time between the actuating and sensing locations and even more profound by the sensing and estimation with the scanning sensor.The dynamic bandwidth of the control system could be increased with full web-width measurement arrays.The effect of different measurement scenarios and scanner path arrangements to estimation error are given in (Tyler and Morari, 1995;Chang et al., 2001).Efforts on finding a feasible and cost-effective solution for full web width measurement systems have been made (Francis, 1991;Williams et al., 1996;Ferguson, 1997;Raunio et al., 2010;Söderberg et al., 2010).One early method (Francis, 1991) combines the light transmittance information from a stationary array of sensors and from a scanning optical sensor at the wet end of the paper machine with the correct information from the dry end scanning sensor.Another system (Williams et al., 1996;Ferguson, 1997) uses CCD cameras to capture the two-dimensional variations within a sampling interval of few seconds.Other attempts of utilizing imaging measurements in paper and board machine estimation and control are presented in Raunio et al. (2010) and Söderberg et al. (2010).This emerging technology is expected to collect the CD profile data with substantially higher frequency than with the traditional scanning sensor.Hence more information on the dynamic variations of cross direction profiles would be available, the CD control interval could be decreased significantly and more of the two dimensional variation could be attenuated.However, the computational demand increases and MPC may fail to find a solution within the sampling time less than one scan.
The outline of this paper is following.First the principles of a cross direction process are given and ex-tended to handle multiple properties, followed by the derivation of the controller, both to spatial and temporal part.Also a novel tuning method is introduced.Simulation results are presented in a separate section and discussion based on the results is given.Finally, conclusions are drawn in the last section.

Modelling and control
Here, an option to multiple property CD-MPC is presented.The controller follows the idea of Chen (1997Chen ( , 1999)), where the spatial and temporal effects were controlled in sequence, but now the idea is lifted to handle multiple property control.Optimal control performance is maintained using a high-resolution profile data.An additional feature of the controller is the actuator constraint handling, where hard constraints are expected to be satisfied by re-tuning the spatial quadratic controller in an automated manner.The simulated examples given in Section 3 describe a process with three actuator sets and two properties.

Process model
The starting point for the multiple property model is the mathematical presentation for a single property (e.g. the basis weight in the paper machine).We use the general assumptions of the input-output relationships in the sheet forming process: the temporal part of the model can be separated from the spatial part of the model, the responses of each actuator are identical in the time domain, and the spatial response between N outputs and M manipulated variables is described with a diagonal interaction matrix G of dimensions NxM.For a single property, the process is described with a discrete-time model eq.(1).
Here, y(k) is the vector of the property values in CDdirection and u(k) the corresponding vector of the manipulated variable.In the equation above, the temporal part of the model g q −1 is usually described with the first order and delay transfer function eq.( 2).
The interaction matrix G consists of row vectors, each of them describing the effect of a single actuator to the width of the sheet.In practice, the responses of the actuators near the edges are different from those near the center, but here similar response shapes throughout the sheet are assumed and the responses near the edges are truncated.Most of the entries in interaction matrix are zero and there are 2p-1 non-zero elements in each row.Hence the matrix can be written as eq.( 3).
The vector d(k) describes the process disturbances, whose effect on the control actions are expected to be compensated for.For the multiple property process with I properties and J actuator sets, the models from different properties are simply stacked together and interaction matrices from different inputs to various outputs are expressed with one block-diagonal matrix.Hence, for I properties the model becomes eq.( 4).In eq.( 4) G ij describes the relationship between property i and actuator set j. Notice that the matrix G ij might have only zero entries.

Control approach
Due to the non-square MIMO characteristics, i.e. the number of inputs being (much) less than the number of outputs, it is impossible to produce a flat surface and for one property the objective is to minimize the variance of the profile across the web eq.( 5). (5) The optimal solution to the quadratic problem can be achieved with the following conventional control law eq.(6).
In eq.( 6), the optimal controller K is an inversion of the process model g q −1 G.However, an inversion does not exist, since the delay in the dynamic part of the model cannot be inverted and the system is not square.The simplest solution to the problem is to assume steady-state process (i.e.g q −1 = 1 and d(k+1) = d(k) ) and use a pseudo-inverse of the matrix In sheet forming processes, it is usually of interest to weight the profile deviations across the sheet and restrict the movement of actuators.Hence, the controller objective is extended to penalize the control actions and profile deviations in quadratic manner as eq.( 7).
In eq.( 7), weighting matrices Q and R are diagonal or band-diagonal matrices with proper dimensions.With the steady-state assumption above, the optimal solution can be calculated from eq.( 8).
The derivation of the controller is similar for multiple property process once the inputs u i and outputs y i are collected to vectors U and Y and the interaction matrices G ij are expressed as one block-diagonal matrix.The dimension expansion is also necessary for the weighting matrices Q and R. The controller now gives the optimal steady-state control moves with interactions between different properties and actuator sets incorporated.However, the steady-state assumption requires a long enough sampling interval that the time delay and actuator dynamics will die out and will not cause instability.This interval will most likely become too long to achieve the best possible control performance.
The assumption on the separable temporal and spatial parts of the model allows us to continue the control synthesis by designing a temporal controller to be used in sequence with the spatial, or steady-state, controller.An efficient option is internal model controller (IMC), which can be set to each actuator zone mainly to compensate for the dead time and modeling errors.With this arrangement, the sampling interval can be set to any value yet ensuring stable control performance.This controller structure separating the spatial and temporal parts was selected over any non-square dynamic MIMO controller structure (MPC or IMC) due to its computational performance.In general, IMC can handle any system by factorizing the process model eq.(9).
In eq.( 9), g + (q −1 ) contains the non-minimum phase characteristics.The best possible controller is then attained by using the realizable inverse g −1 − (q −1 ) as a controller.The controller will be designed more robust by introducing an exponential filter F (q −1 ), which also smoothens the overall control actions.For more comprehensive description of IMC, the reader is referred to Garcia and Morari (1985).The design applies also for systems with different dynamics and dead times, but here we consider a case where the filtering time constant f j is the same within every actuator set j.
Since the actuator dynamics were also assumed similar within every actuator set j, all controllers have the same form of eq.( 10) and the control action can be calculated from eq.( 11).
Now, the IMC for each actuator set can be coded efficiently as time series objects.The controller derived is almost similar to one proposed in Chen (1997Chen ( , 1999)).
The spatial controller finds the instantaneous changes to the optimal controller profile ensuring optimal unconstrained steady-state performance with respect to multiple properties.The temporal controller dynamically filters those control actions to ensure stability.Current control structure does not, however, account the dynamical interactions.Additionally, as recognized e.g. in Backström et al. (2001), the spatial controller sets penalties only to incremental control actions and hence the risk of actuator saturation is still evident.The temporal controller also restricts only the incremental control actions by filtering.Therefore, an anti wind-up solution is needed in the temporal controller to take into account the hard constraint, i.e. actuator saturation limits.

Input saturation and automated tuning
As demonstrated by (Garcia and Morari, 1985), the IMC remains stable while the information of input saturation is also fed to the internal model.However, in our system, clipping the inputs may produce a suboptimal performance and CD control actions will be aliased to MD variation.There exist some tailored anti wind-up solutions also to cross direction problems with linear controllers (Rojas et al., 2002), as well as tuning procedures ensuring the robustness under input constraints (Duncan andHeath, 2008, 2010;Heath and Wills, 2004), but we propose a different procedure.
Our novel solution is to automatically tune the spatial controller by decreasing the aggressiveness of control actions i.e. we are avoiding excessive control signals and severe input clippings.The automated tuning may incorporate a pre-defined weighting matrix for control actions.In addition, the weighting matrix R is re-evaluated at each sampling time, where additional weight is given depending on the current distance to the saturation limits.Figure 2 shows the relation between the actuator saturation limits and the input weighting.
In Figure 2, the relation comes from a second-order equation eq.( 12).tween the distance to the actuator saturation limit and the output penalty.
In eq.( 12), r m is the input weighting, r 0 is the predefined weighting factor, or the lower limit for the weighting, r max is the upper limit for the additional penalty, u max is the maximum (minimum) allowable control input i.e. saturation limit and u d = u max − u is the current distance to the saturation limit.In addition, one tuning parameter comes from the timeconstant of the filter in IMC.In Heath and Wills (2004), a rule for tuning the filter time constant was given as a relation between the allowable rate of control movement and the range of control inputs eq.( 13).
The final control structure is presented in Figure 3. Notice that it is assumed in Figure 3 that the control target (reference) is zero for all outputs.

Simulation results
The controller developed was coded with Matlab R and implemented to the paper machine simulator presented elsewhere (Ylisaari et al., 2009;Ohenoja et al., 2010).It is assumed in our simulations that the profile data is updated every second based on imaging measurements.Information about simulating the imaging measurements and estimation can also be found from Ohenoja et al. (2010).The control interval is also set to one second.For the simulated system with dimensions 800x164, it takes 0.002 seconds to perform the calculations for one control move with a dual core desktop computer with CPU speed of 2.2 GHz and a memory of 4 GB.
The performance of the control system can be evaluated in numerous ways.One feasible performance indicator for multiple property control was given in Haznedat and Arkun (2002), see eq.( 14).In eq.( 14), y p is the controlled profile and d p is the disturbance, or open-loop profile for property p.The first term describes the performance for single property and the second term ensures the balance between different properties.If no control actions are applied, the standard deviations std(y p ) and std(d p ) will be equal, the second term will be zero and hence the maximum value of κ is the number of properties evaluated.Additionally, we will evaluate the control performance also in terms of actuating energy.
We simulate two properties and three actuator sets.Both properties have 400 measurement positions and there are 67, 40 and 57 actuators in different actuator sets.The spatial and dynamic responses of each type of actuators are presented in Figure 4. Table 1 explains the effect between different actuators and properties.These responses may be considered to describe the effect of dilution headbox, steam box and rewet shower profilers to basis weight and moisture.Other simulation settings are presented in Table 2.With the presented system, the input-output model can now be written as eq.( 15) and the weighting matrices Q and R as eq.( 16) and eq.( 17), respectively.
Table 1: Simulated interactions between the different actuator sets and properties.

Steady-state performance
The steady-state disturbance profiles for Property 1 and Property 2 are presented in Figure 5.We illustrate the performance of the multiple property control (MP) presented in Section 2 by comparing the results with a set of solutions from single property control (SP).The latter one uses a feedforward compensation between Actuator 2 and Actuator 3 to avoid unnecessary competing control actions for Property 2. Single property control is also tested using a full feedforward compensation (SPFF) i.e. compensation from Actuator 1 to Property 2 as well.The weights and input saturation limits used in single property control are similar to those for multiple property control.Finally, we also present optimal constrained solution with respect to κ, found using an optimization function 'fmincon' in Matlab R .The results of the steady-state performance are listed in Table 3.

Performance under input saturation
Next figures show simulation results with saturating actuators.The pure CD variations are similar as in previous example, but during the simulation period, a dynamical profile variation is added which will cause  input saturation.The dynamical variation also represents a disturbance that can be assumed to be detected with imaging measurements, but may not reliably be detected with scanner based methods due to sparseness of the measurements and heavy filtering usually applied.The maximum disturbance profiles are presented in Figure 6. Figure 7 presents the dynamical and random variation part of the simulated disturbances for Property 1.The dynamical variation for Property 2 comes from the relation of 0.7•Property 1.The simulation results with and without the automated tuning are presented in Figures 8 and 9.

Discussion
The results of the steady-state performance, listed in Table 3, show how multiple property control outperforms the single property control (SP) as expected.
Multiple property control allows more variance in Property 1 to achieve better overall control result (κ) with minimum actuating energy.Single property control with full feedforward compensation (SPFF) gives slightly better result with respect to κ, but the balance between the two properties and actuating energies are worse than with multiple property control.The optimal constrained solution with respect to κ shows substantial performance increment, but with the cost of increased actuating energy and the number of actua-tors lying in their saturation limits.This would mean very limited actuating possibilities for further profile changes.
The effect of tuning is demonstrated by simulating dynamical variations that will cause input saturation.Figure 8 shows the temporal value of the performance indicator κ throughout the simulation period.It can be observed that we encounter slight performance degradation when using automated tuning.From Figure 9 we observe the actuating energies and saturation violations, respectively.The solution with automated tuning clearly uses the actuators in a different way compared to the case without automated tuning.We observe that there are no saturation violations for Actuator 3 when using automated tuning, whereas without the tuning, some inputs of Actuator set 3 lie on their saturation limits constantly.There are total 5025 saturation violations for the case without the automated tuning, 3947 of them lying in dynamical region (MD positions 201...800).For the case with automated tuning, there are 2860 saturation violations, all of them lying in dynamical region.The results show that this kind of novel tuning method could be applicable.However, the tuning should also include output weighting to avoid unnecessary sluggishness of control.Also the second-order function applied for additional penalty (Figure 2) is quite aggressive as it penalizes all control actions.For example, an "S-shaped function" would give no or only moderate penalty addition when actu-   ator is operating far from the saturation limit.

Conclusions
The results presented above show the benefits of multiple property CD control; we can expect more uniformity when observing all properties of interest and with less actuating energy.The controller presented is a computationally efficient option for the problem described and is expected to be a feasible control method to be used with higher sampling frequency provided by imaging measurements.The automated tuning procedure presented offers a way to avoid severe input saturation, but further development and analysis of the method is required.

Figure 2 :
Figure 2: Second-order curve showing the relation between the distance to the actuator saturation limit and the output penalty.

Figure 4 :Figure 5 :
Figure 4: Actuator response shapes (on top) and dynamics (on bottom) for three types of actuators.Dotted line is for Actuator 1, dashed line for Actuator 2 and solid line for Actuator 3.

Figure 6 :
Figure 6: The effect of dynamical disturbances to profile information for Property 1 (solid line on top) and Property 2 (solid line on bottom).The dotted lines are the corresponding steady-state disturbance profiles plotted also in Figure 5.

Figure 7 :Figure 8 :
Figure 7: The dynamical and random part of the simulated variations.

Figure 9 :
Figure 9: Actuating energies and cumulative sum of saturating actuators.On top, results without automated tuning of spatial controller.On bottom, results with automated tuning of spatial controller.The dotted line is for Actuator 1, dashed line for Actuator 2 and solid line for Actuator 3.

Table 3 :
Steady-state performance of the controllers.MP is multiple property control, SP is single property control, SPFF is single property control with full feedforward compensation.