Experimental System Identification and Black Box Modeling of Hydraulic Directional Control Valve

Directional control valves play a large role in most hydraulic systems. When modeling the hydraulic systems, it is important that both the steady state and dynamic characteristics of the valves are modeled correctly to reproduce the dynamic characteristics of the entire system. In this paper, a proportional valve (Brevini HPV 41) is investigated to identify its dynamic and steady state characteristics. The steady state characteristics are identified by experimental flow curves. The dynamics are determined through frequency response analysis and identified using several transfer functions. The paper also presents a simulation model of the valve describing both steady state and dynamic characteristics. The simulation results are verified through several experiments.


Introduction
Directional control valves (DCVs) are indispensable components in nearly all hydraulic systems. They are used to control the direction and rate of flow between power sources and actuators. Their characteristics, both steady-state and dynamic, directly influence the behavior and performance of the system they are used to control. DCVs are used for a wide range of applications and therefore exist in many different variations which may be divided into different categories. One category is servo valves and high-performance proportional valves; either electrically actuated (proportional valves) or electro-hydraulically actuated with one or more pilot stages (both servo and proportional valves). They all include some kind of closed loop spool position control and are used in applications where fast response, high precision and repeatability are required. Another category of DCVs is the so-called mobile proportional valves (MPVs) used for mobile applications such as, tractors, wheel loaders and truck mounted cranes, i.e., open loop motion control systems where an operator closes the control loop. They are normally pressure compensated, making them load independent, and may be mechanically, hydraulically or electro-hydraulically actuated. In the latter case they usually include onboard electronics with closed loop spool position control and dither function for elimination of static friction (stiction) between spool and valve housing. MPVs with electro-hydraulic actuation (EHA) are increasingly being used in closed loop motion control systems and often replace servo valves in systems with relatively low performance requirements. The main motivation for doing so is that MPVs are much cheaper than servo valves. In addition, the MPVs are more flexible and easier to integrate in a system due to their modular design. However, in terms of response and precision, MPVs cannot be compared to servo valves. Therefore, to qualify the use of MPVs, it is often necessary to carry out dynamic simulations to analyze and verify the performance of the system where the valves are being used. This requires accurate valve models which usually need to be experimentally verified, because the required model data often is doi:10.4173/mic.2015.4.3 insufficient or not existing.
In literature there are three ways of modeling a proportional valve; using a black box model, gray box model and a white box model. The black box model consists of empirically determined functions that are found by matching test data with some function. This could be used to make a polynomial of the flow characteristics or a transfer function to characterize the dynamics of the valve. Black box models do not contain any reference to the physics of the plant. White box models are made according to physical laws of the plant and should therefore be an accurate model of the system. The gray box model represents a mix of both the black and the white box models. Valve dynamics are often modeled as second order systems as it closely resembles the dynamics up to its eigenfrequency Tørdal and Klausen (2013); Bak and Hansen (2012). Eryilmaz and Wilson (2006) developed a nonlinear mathematical model for a generic proportional control valve which should be usable for a variety of different control valves if the geometry is known. Valdés et al. (2008) made a CFD model of a ABS NO valve in order to extract valve flow during different pressures and openings. Further they created an analytic model of the geometry of the valve and discovered that they give comparable results during turbulent flow. Niksefat and Sepehri (1999) made a hydraulic force controller by using nonlinear Quantitative Feedback Theory (QFT). For this purpose they had to make a complete nonlinear mathematical model of the actuator system consisting of a valve and a cylinder, and use the Golubev method to derive a rational transfer function of the system. Xu et al. (2014) modeled a large flow capacity three-stage valve to be used in large scale hydraulic systems by using white box analytic equations. Amirante et al. (2013) made a CFD model of a valve and optimized the design to reduce the flow forces acting on the spool during a prescribed flow rate. Posa et al. (2013) studied the flow through a 4/3-directional control valve using a CFD analysis during different openings. They concluded that the axial forces (pressure forces) on the spool was roughly proportional to the pressure drop and the valve opening. Dasgupta and Murrenhoff (2011) modeled the dynamics of a servo valve using a bondgraph simulation technique based on a white-box model of the inner geometry of the servo valve.
In this paper an MPV-type of DCV is analyzed and it is shown that the dynamics can be more accurately described with a third order system. An approach for testing the valve and processing measured data in order to visualize the valve dynamics and carry out system identification is described. The approach is quite elaborate but by comparing it with other approaches to identify the valve dynamics it is shown that it is the only one yielding accurate results. Finally a valve model is proposed which is validated and verified with experimental results.

Considered Valve
MPVs are widely used in both onshore (off-highway) and offshore applications mainly for systems with variable supply pressure (load sensing systems) but also for systems with constant supply pressure. They are available from manufacturers such as Danfoss, Brevini, HAWE, Parker, Bosch Rexroth, HUSCO and others. Common for nearly all MPVs is the modular design which makes it possible to configure customized valve groups for different applications with relatively little effort. The considered valve is a Brevini HPV41, schematically illustrated in Figure 1, consisting of the following modules:  The pressure reducing valve in the supply module reduces the supply pressure to a level needed for the pilot stage, which is integrated in the EHA module. The pressure relief valve protects the pilot stage from pressure peaks.
Via the main spool the load pressure, p A or p B , is transmitted through the LS circuit to the pressure compensator; a pressure reducing valve with a variable pressure setting which is offset by a constant value corresponding to the pre-compression of the spring. Due to this offset the pressure, p C , between the compensator and the main spool is always reduced to give an approximately constant pressure drop across the main spool metering edge. Thereby the controlled flow is independent of the load pressure and proportional to the position of the main spool. For this reason pressure compensated MPVs are also referred to as load independent.
The main components of the EHA module, see Figure 2, are two proportional pressure reducing valves (PRVs) controlled by an integrated electronic circuit which includes a microprocessor and an LVDT for measuring of the main spool position. Based on the deviation between the actual spool position and the control signal, a spool position reference, the two PRVs are controlled to create an axial pressure difference on the main spool, causing it to move. The details of the control scheme for the valve are, however, out of the scope of this paper. Further information about the HPV41 and the EHA module (the MHPED) is given by Brevini Fluid Power (2012). Figure 2: EHA main components.

Frequency Response Analysis
The frequency response tests of the hydraulic valve have been carried out in the Mechatronics laboratory at the University of Agder. The testing procedure was carried out in accordance with the British Standard for hydraulic valve testing BS (2009). The MPV was instrumented and mounted in a test bench as shown in Figure 3. The pump connected to the P-port of the valve delivers a constant pressure of 100 bar and a maximum flow of 100 l/min. The testing procedure has been programmed using a LabVIEW program which easily communicates with the DAQ board. The LabVIEW program generates a voltage signal which is controlling the modules mounted on the valve. The valve was tested with three different amplitudes, the amplitudes are 5%,10% and 20%. The test offsets are ±50%.
The valve has to be tested with an offset since the valve is highly non-linear around the middle spool position. This non-linearity is caused by the deadband characteristics in the spool opening geometry. The frequency response analyses Nise (2011) is carried out using traditional methods. The input signal is compared with the output response of the system. The phase lag and the magnitude difference between the input and the output signal have been measured for each single test frequency. The test frequencies ranged in [0.1 · · · 4.9]Hz, and the step between each test was chosen to be 0.2Hz. A simplified illustration of the frequency response test is shown by Figure 4. Figure 4 shows that the input and output are equal at low frequencies. When the frequency is increased, both a phase lag and a magnitude difference is introduced. The input and output amplitude of each test Tørdal and Klausen (2013) frequency can be found and the magnitude M is calculated using Equation 1. (1) A is the measured amplitude at the given test frequency. The initial amplitude A @0.1Hz is the sine amplitude of the first frequency. This will normalize the Bode plot to begin at 0 dB. The phase lag between the input and output signals for each test frequency is needed in order to have a bode plot representation of the linearised valve dynamics. The phase lag φ is calculated using Equation 2.
T p is the period of the sine waves and δt is the time lag between the input and measured response. By testing the valve on several frequencies, the magnitude and phase can be calculated for each respective test frequency. The test results are stored in two vectors which can be illustrated by Equation 3.
Where M is the magnitude vector in decibel [dB] and φ is the phase vector in degrees. The result of the frequency response test can be represented using a bode plot which illustrates the change in magnitude and phase as the frequency increases.

Data Processing
The measured data from the flow sensor contains noise. Such noise could be removed using a zero-phase forward-backward lowpass filter, however the amplitude and phase needs to be extracted afterwards. Instead of filtering, fitting a sine wave to the measurements using an optimization algorithm will give satisfying results. Fitting a sine wave is also less time demanding since the fitted sine signal will describe the amplitude, offset and phase lag immediately. The approximated measured flowQ which is a close fit to the actual measured flow is given by Equation 4.
B is the offset, A is the amplitude, f is the test frequency and φ is the phase lag. These four parameters are identified using an optimization algorithm. The optimization algorithm tunes the four unknown parameters in order to minimize the error between the measured flow Q and the approximated flowQ. The absolute error vector between the approximation and the measured data is given in Equation 5.
This error should be minimized in order to have a close fit between the measurements and the approximation. The objective function which is used to minimize the error is described by Equation 6.
The optimization is carried out for each single test frequency and the four unknown parameters are used to represent the approximated flow for further analysis. The optimization algorithm is provided by Matlab's fmincon function Byrd et al. (1999). An illustration of the result of using the optimization in Matlab is shown in Figure 5. Figure 5 shows that the approximated flow signal (blue) represents the measured flow (red) in a satisfactory way. The amplitude of the approximated flow signal from Equation 4 is used to find the magnitude of the response using Equation 1, and the phase of the response is equal to the phase of the approximated flow signal.

Steady State Flow Characteristics
In order to investigate the steady-state flow characteristic of the valve, the flow is measured during a steady- In Figure 6 the x-axis is the reference flow signal sent to the valve, and the y-axis shows the measured flow response. The red dots represents the measured data and the blue line represents the fitted polynomial function. This plot clearly shows a deadband in the range Q ref ∈ [−20l/min, 20l/min] and a maximum flow reaching ±70l/min at maximum opening. The measured deadband is expected and the maximum flow is correct compared to the datasheet of the valve,Brevini Fluid Power (2012).

Dynamic Characteristics
Several tests of the valve were carried out. The valve was tested with three input signal amplitudes and in both negative and positive offset configurations. The first three tests represent the positive offset test configuration. The results are given by the bode plot in Figure 7  It can be seen that the results are quite similar when looking at the phase plot. But the magnitude of the three different test amplitudes deviates from each other. This is probably caused by the friction nonlinearities introduced by the sliding spool inside the valve housing. The valve has also been tested in a negative offset in order to compare the results with the positive offset test configuration. The negative offset test configuration is presented by the bode plot in  Figure 8 also indicates that the magnitude is more dependent to the test input amplitude. When comparing the two test configurations, the test results are not similar. The positive offset test indicates a faster dynamic behavior in terms of phase lag and different magnitudes at 5% and 10% sine amplitude input.

System Identification
The system identification is performed using Matlab's System Identification Toolbox, Garnier et al. (2003); Ljung (2009); Young and Jakeman (1980). The toolbox has a useful graphical user interface (GUI) which is used to import the frequency domain data from the lab test results presented in Section 5.2. The toolbox is used to estimate time continuous transfer functions from bode frequency response data, which should characterize the dynamics between the reference flow and the measured flow response. An overview of the system identification toolbox is shown in Figures 9 and 10.  1. The data type is set to be frequency domain data.
2. The data format of the experimental data is given as three vectors containing the magnitude M , phase φ and the frequency in rad/s.

3.
Here the variable reference to the three vectors stored in Matlab's workspace is specified.
4. The sampling interval is set to 0 to estimate a continuous Transfer Function. Since the sampling interval of the experimental data is much faster than the systems bandwidth, it should therefore give satisfactory results.
5. Press this button to import the experimental frequency domain data.
6. Choose "Transfer Function Models" from the drop down menu, then specify number of poles and zeros in your estimated transfer function. In this case 1 zero and 3 poles are used since this gives the closest fit to the experimental data.
7. The information about the approximated transfer function is found here.
To illustrate the result of using the system identification toolbox a bode plot of the estimated transfer function is compared with the input experimental data. A test performed in a negative offset and 10% sine input amplitude has been chosen to illustrate the fit between an estimated transfer function and the experimental data. The comparison is given by Figure 11. Lab Test TF Estimate Figure 11: Example in 50% negative offset.
In Figure 11 the red dots show the experimental frequency response and the blue line shows the estimated transfer function. As the figure shows, the fit between the estimated model and the experimental data is satisfactory. The system identification process can be repeated using the GUI, but the process could also be implemented using Matlab functions in a Matlab script. Such implementation would be favorable if large amounts of experimental data is to be processed. The resulting transfer function found by using the system identification toolbox can be expressed using a general form given by Equation 7.
Q is the valve flow and Q (ref ) is the reference flow. The coefficients [b 1 , b 0 ] and [a 2 , a 1 , a 0 ] are found from the system identification process.
Three transfer functions are estimated from test data sets at 50% positive offset, the results are illustrated in Figure 12. The fitted transfer functions shown in Figure 12 represents the experimental test data given by Figure 7. The coefficients describing Equation 7 which are found using the system identification method are given in Table 1. The fit to the experimental data is given in percent to the right in the table. The negative offset tests is also represented using three fitted transfer functions. Figure 13 shows the bode diagrams of the three fitted transfer functions. The accompanying coefficients describing the three bode plots shown in Figure 13 are given by Table 2.
Because b 1 is negative in Tables 1 and 2, it gives that the systems described by these coefficients have a zero in the right-half-plane. A right-half-plane zero produces a +20dB/decade slope on the magnitude, while reducing the phase by 90 • . This behavior is most likely produced by the transfer function estimator to create a total of 360 • phase lag while limiting the total magnitude slope to -40dB/decade. However a right-halfplane zero produces an unfavorable behavior in the time domain which is easiest shown via a step-test. A step-test of the system given in Equation 8 is shown in Figure 14 G(s) = −187.3s + 8069 s 3 + 37.7s 2 + 1308s + 7789 Step response Figure 14: Step test response of a system with a righthalf-plane zero.
The step response in Figure 14 shows a small dip at the start, i.e. the response starts negative before it goes positive and settles at a final value. This behavior is due to the right-half-plane zero. On the valve system this would result in a small negative opening in the beginning of a step input, but this is not the case in the real world. This is clearly a disadvantage of using the transfer function estimator as it does not factor in the time domain response.

Transfer Function Estimation Comparison
To estimate a third order transfer function using Matlab's toolbox, the user has to perform a frequency response test and analysis of the valve. This is a costly and time consuming method, and one should know the benefits of performing a frequency response analysis compared to a simple step response test. A step response should include the dynamic responses from all frequency ranges and calculating a first or second order transfer function from this response is faster than performing a frequency response test. In this section a first and second order transfer function are calculated from a step test and compared with a frequency response analysis. A step test to 80% opening of the valve with a measured flow response is shown in Figure  15. Step test with annotated settling and rise time.
In Figure 15 the blue curve is the step signal sent to the valve converted into reference flow, and the red curve shows the measured flow response. The double arrow annotated T r shows the rise time which is between 10% and 90% of final value and T s shows the time where the flow signal is 98% of the final value. The error between the step input and flow measurements origins from the non-linear opening area of the spool. The measured rise and settle times are given by Equations 9 and 10: From this information two first order transfer functions are estimated. A first order transfer function is given by: The constant a is calculated with the following equations: Since the two constants calculated in Equations 12 and 13 are not equal, the resulting transfer functions may not be accurate. A better approximation of the system comes from estimating a second order transfer function. A second order system is given by: where ω n is the natural frequency of the system and ζ is the damping ratio. In order to determine these two constants, the two Equations 15 and 16 are used: T r · ω n = 1.76ζ 3 − 0.417ζ 2 + 1.039ζ + 1 where Equation 16 is a polynomial approximation for the normalized rise time for an underdamped second order system Nise (2011) (page 181). These two equations are solved by substituting Equation 15 into Equation 16 and solving for ζ: The roots are solved numerically and the available values for ζ are shown in Table 3. The four transfer functions and the actual frequency response is shown in the bode plot in Figure 16. In Figure 16, the red and blue curve show the two first order transfer functions given in Equations 19 and 20. Compared with the test data shown as red dots, these two first order systems are close on the magnitude and phase up to 1 Hz, but diverge from the test data thereafter. The green curve, which shows an estimated second order transfer function which is given in Equation 21, manages to get a larger phase on the system at higher frequencies, but it is not a good fit with the test data. On the magnitude it is a good fit up to roughly 1Hz, but on the phase plot it diverges from the test data quite quickly at roughly 0.3Hz. The best approximation of the test data comes from an estimated third order transfer function which is shown as a cyan-colored curve. It is practically equal to the test data set on all the tested frequencies.
The four shown transfer functions have different accuracies when compared to a measured frequency response, but the time needed to acquire the necessary information and to calculate a transfer function is often the inverse of the accuracy. In this case, the first order transfer function approximated from a step response is a fine approximation up to 1Hz and the time required is small compared to performing a frequency response analysis. However the performance of an estimated transfer function is much better when fitted to a frequency response rather than a step response. Therefore the user should decide which accuracy is required of the transfer function in his model before spending too much time on tests.

Modeling, Simulation and Verification
In this section the valve is modeled and simulated in Matlab, and the time response is compared with a measured flow response from the experimental setup. The model structure is shown in Figure 17.  Figure 17 introduces a mean transfer function G M (s) describing the dynamics of the valve. It is identified by calculating a mean phase and gain from the six experimental tests shown in Figures 7 and 8. The mean transfer function is fitted using the system identification toolbox described in Section 6. The flow curve is extracted from Figure 6 and implemented using a look up table.
The first verification test is carried out by testing the valve with a sine input flow reference Q (ref ) with an amplitude of 70l/min and a frequency of 0.1Hz. The result is shown in Figure 18. The measured flow signal (blue) and the simulated flow signal (red) match satisfactory. Also the valves dead band is included in the simulation results. Two additional tests have been carried out at both 0.3Hz and 0.5Hz in order to validate the model further. The results are shown in Figures 19 and 20. The results from testing the valve at both 0.3Hz and 0.5Hz also state the model gives a quite accurate description of both the steady-state and the dynamics.

Conclusion
In this paper a mobile type of directional control valve (DCV) -a Brevini HPV41 -has been investigated and a method for analyzing it's dynamic behavior has been presented. The model was identified by measuring the steady state and dynamic characteristics of the flow generated by movements of the spool position. The flow curves in Section 5.1 describe how much flow the spool lets through at a certain opening, e.g. the steady state characteristics. The dynamic characteristics were identified by performing a frequency response analysis of the flow response shown in Section 5.2. The dynamic characteristics were described as continuous transfer functions by using a system identification toolbox in Matlab as described in Section 6. These transfer functions had a fit percentage of roughly 90% which is a good fit. Finally a complete black box model of the valve with an MHPED module was simulated in Section 8. The comparison figures in this section show that the black box model almost matches the measured flow data and therefore the model is verified. We have with this result shown that the presented methods are effective of making a black box model of a hydraulic valve.