Analytical, Experimental, and Numerical Investigation of Energy in Hydraulic Cylinder Dynamics of Agriculture Scale Excavators

: This paper discusses energy behaviors in hydraulic cylinder dynamics, which are important for model-based control of agriculture scale excavators. First, we review hydraulic cylinder dynamics and update our physical parameter identiﬁcation method to agriculture scale experimental excavators in order to construct a nominal numerical simulator. Second, we analyze the energy behaviors from the port-Hamiltonian point of view which provides many links to model-based control at laboratory scale at least. At agriculture scale, even though the nominal numerical simulator is much simpler than an experimental excavator, the analytical, experimental, and numerical energy behaviors are very close to each other. This implies that the port-Hamiltonian point of view will be applicable in agriculture scale against modeling errors.


Introduction
Unmanned hydraulic machinery is gaining popularity and a new generation of modelbased control of hydraulic machinery is coming (e.g., [1][2][3][4][5]). In this paper, we analyze energy behaviors of hydraulic cylinder dynamics, which are important for model-based control of agriculture scale excavators. At the agriculture scale, the most widely used controls are PID or PID-like control since the physical parameters are absent and we are free from any physical parameter identification. However, model-based control has the possibility to achieve better performance than PID control [6][7][8]. Moreover, good PID control sometimes has a special structure which is meaningfully clarified in the context of modelbased control. Especially, in order to prevent a serious accident or to explain the accident procedure, model-based control would play the key role. Even against modeling errors, the robust control theory such as the port-Hamiltonian theory guarantees robustness [9,10]. By the passivity [9,11], port-Hamiltonian systems and control can achieve a large region of attraction.
Regarding model-based control of hydraulic machinery, the paper analyzes energy behaviors in hydraulic cylinder dynamics from the port-Hamiltonian point of view, which provides many links to model-based control at laboratory scale at least [12,13]. Note that the analyzed energy in this paper is not the global energy (e.g., [14][15][16][17]) defined for the whole hydraulic machinery but just one of the local energies defined for each joint whose analysis is still open [18] and critical for better understanding of the whole hydraulic machinery eventually. To the best of our knowledge, however, on agricultural scales that should be larger than the laboratory scale, such local energy is mainly analyzed in terms of statics [1] and is not fully analyzed in terms of dynamics. On the agricultural scale, numerical energy behavior and experimental energy behavior in hydraulic cylinder dynamics have not been analyzed yet from the port-Hamiltonian point of view [12,13,18,19]. In the paper, on the agricultural scale, we discuss the analytical, experimental, and numerical energy behaviors of a Hamiltonian in hydraulic cylinder dynamics.
The paper has two main results. As the first main result, our previous physical parameter identification method [20,21] for laboratory scale hydraulic arms with closed-center valves is updated for agriculture scale excavators with open-center valves. The updated physical parameter identification method is valuable since no trial and error is needed and the physical parameters are uniquely identified by using experimental input and outputs only. As the second main result, we analyze the analytical, experimental, and numerical energy behaviors in hydraulic cylinder dynamics from the port-Hamiltonian point of view and show that these behaviors are very close to each other. This is an unexpected result since the hydraulic cylinder dynamics in the nominal numerical simulator are much simpler than the experimental dynamics or Computational Fluid Dynamics (e.g., [22][23][24]). This implies that the port-Hamiltonian point of view would be applicable on the agricultural scale against modeling errors as well as on the laboratory scale.
The paper is organized as follows. In Section 2, we review hydraulic cylinder dynamics and the pump pressure dynamics. In Section 3, we update our physical parameter identification method to the experimental excavators in order to construct the nominal numerical simulator. In Section 4, we discuss the cross validation in which we use a new experimental input and outputs that are never used in the physical parameter identification. In Section 5, we discuss, experimentally, numerically, and analytically, the energy behaviors from the port-Hamiltonian point of view on the agricultural scale instead of on the laboratory scale. In Section 6, the paper is concluded.

Nonlinear Nominal Model
Let us consider a hydraulic arm whose nonlinear nominal model of Figures 1 and 2 is a coupling of the hydraulic cylinder dynamics [18] and a pump pressure dynamics defined later. First, the hydraulic cylinder dynamics is described by where the joint angle q i (t) (rad), the joint angle velocityq i (t) (rad/s), the cap side pressure p +i (t) (Pa), the rod side pressure p −i (t) (Pa), and the orifice area (the input) u i (t) (m 2 ) are the functions of time t [s]. The drift term f 0 (x) and the input term g o (x)u are defined as follows: where the subscript i = 1 or 2 is the joint number and the subscript + or − denotes the cap side and the rod side, respectively. The piston displacement s i (t) (m), the pump pressure p p (t) (Pa), and the tank pressure p t (t) (Pa) are the functions of time t (s). The cylinder Jacobian matrix G ii (q i ) (m) is the function of the joint angle q i (t) (rad). The external force F ei is reaction force from the piston end walls. The physical parameters are the link mass m i (kg), the link moment of inertia around the center of the gravity I i (kg · m 2 ), the stroke L i (m), the link length i (m), the distance between the joint center and the piston endpoint ij (m), the distance between the joint center and the center of the gravity gi (m), the piston minimal length ci (m), the piston areas A ±i (m 2 ), the damping coefficient F di (N · s/m), the Coulomb friction coefficient F ci (N), the gravity acceleration g (m/s 2 ), the gravity coefficient W i (m/s 2 ), the bulk modulus b (Pa), the flow gain coefficient C f ±i (m/(s · √ Pa)), the internal leakage coefficient C ii (= C i1 or C i2 ) (m 3 /(s · Pa)), the external leakage coefficient C e±i (m 3 /(s · Pa)), and the pipe volume V min ±i (m 3 ). The nominal model introduces the restricted domain s i (t) ∈ (−L i /2, L i /2). The origin of the piston displacement s i (t) and the cylinder pressure p ±i (t) are the middle point of the stroke range, air pressure, respectively. C f +i h +i (p +i , p p , p t ) and C f −i h −i (p −i , p p , p t ) are approximated by Bernoulli's principle.  Second, the pump pressure dynamics are described bẏ where the pump bulk modulus b p (Pa) and the pipe volume V p (m 3 ) are the physical parameters. The pump flow rate Q p (m 3 /s) is the function of time. The input flow rate Q i (m 3 /s) and the tank flow rate Q t (m 3 /s) are approximated by Bernoulli's principle as hydraulic cylinder dynamics: where the tank flow coefficient C f t (m/(s · √ Pa)) and the maximum orifice area u max i (m 2 ) are the physical parameters. (1) ignores the piston mass, the cylinder mass (cf. [4]), and the nonlinear friction effect at least. Equation (2) assumes the steady flow and the negligible servo dynamics of an open-center valve [25,26]. See [18] for more details of the nonlinear nominal model.

Physical Parameter Identification
In this section, we update our previous physical parameter identification method [20,21] for agriculture scale excavators with open-center valves. The previous method cannot be applicable to any open-center valves but the updated method is applicable to closed-center valves as a special case. More precisely, the updated method can discuss the pump pressure dynamics as well as the hydraulic cylinder dynamics. In the rest of this paper, just for simplicity, without loss of generality, we show only the case of joint 2 as shown in Figure 3.

Experimental Condition (Identification)
From Equations (1) and (2), the equations of motion for the mechanical and the continuity equations for fluid systems are as follows: and <STEP1> It is possible to see that Equations (4) and (5) can be linear equations with respect to the unknown physical parameters.
Note that the right-hand sides are not zero and known.
is obtained at each time t = t 1 , · · · , t N . From the project theorem in Hilbert space [27], the set of the unknown physical parameter a 1 and a 2 in Equations (6) and (7) are identified uniquely as the following values: Remark 2. Apparently the obtained equation is similar to that in our previous paper [20,21]. However, as stated above, the pump pressure dynamics are not counted in the previous paper but are taken into account in this paper only. In this sense, the above method in the paper is updated from the previous method. Figure 4 shows an appearance of the experimental setup. The experimental setup consists of an agricultural scale excavator (HITACHI, ZX30U-5A), pipes (BRIDGESTONE, PFH06, 3/8 inch), and a servo open-center valve only. In the paper, the second joint of the excavator is focused precisely. The excavator is driven by the servo open-center valve of Figure 5 whose orifice area u 2 is measured by a displacement sensor (MIDORI PRECISIONS, LP-20F) and a DC servo motor (SANYO, T720-0.12EL8, 0.2 kW), a rack and a pinion (reduction ratio 1/25), an AD board (interface, PCI-3155, 16 bit), a control PC (Linux, 3.8 GHz, 4.0 GB) [28], and a DA board (interface, PCI-3325, 12 bit). In the control PC, a 2DOF local controller is designed and implemented to make the orifice area u 2 into the actual control input against the large overlapping in the valve. The open-center valve has a direct relief valve. The outputs of the excavator are a piston displacement s 2 (t) measured by a piston displacement sensor (MIDORI PRECISIONS, CPP-45-150LS), the cap, rod, pump, and tank pressures p +2 (t), p −2 (t), p p (t), and p t (t) measured by the pressure sensors (KEYENCE, AP-16S) and the oil temperature measured by a thermocouple (NIHONDENSOKU, TNI-3.2-10-EXA4M). The oil temperature is 27 • C. Figure 6 shows the system configuration in which an ideal sinusoidal signal r 2 and a control input v 2 to the DC servo motor. Table 1 shows the value of the known physical parameters before the physical parameter identification.
External leakage coefficient Unknown Since the physical parameter identification is in the off-line case, the time derivative of the piston displacementṡ 2 (t) and the cap and rod sides' pressuresṗ ±2 (t) are sufficiently estimated via the first difference approximation and the standard moving average (43order). In order to identify the physical parameters, the ideal sinusoidal signal r 2 (t) = 100 sin(2π0.5t) (m 2 ) is applied and we cannot provide slower cases any more because the piston collides with the cylinder ends. The maximum input u 1max and u 2max is not used in the experiment.      Figure 7 shows the experimental results from the ideal sinusoidal signal r 2 (t) = 100 sin(2π0.5t) (m 2 ). The input u 2 was not equal to the ideal sinusoidal signal r 2 around the zero mark. This is due to the large overlapping. However, except for this overlapping effect, the input u 2 was equal to the almost sinusoidal signal r 2 . The time responses of the pressures and the piston displacement p +2 (t), p −2 (t), p p (t), p t (t), and s 2 (t) were measured and the time responses of the piston velocityṡ 0 (t) was estimated. According to the input u 2 , every output responded sufficiently. The resolution and ranges of all sensors have no problems. In these senses, against actual noise and disturbances, the experimental setup is precisely justified to analyze the hydraulic cylinder dynamics on the agriculture scale instead of the laboratory scale. The relief valve did not work and the pressures did not saturate. Table 2 shows the result of the physical parameter identification. In spite of the fact that the proposed method does not impose any restrictions on the sign or value range of the parameters and does not require any initial estimates of the parameters, all the identified physical parameter values were positive. The condition numbers of X 1N were 26 and 220 where the number of the row was 17,500. The condition numbers were sufficiently small in this experiment.

Experimental and Numerical Condition (Cross Validation)
First, the nominal numerical simulator is constructed with the physical parameters in Table 2 on the laptop (2.2 GHz, 4.0 GB, MATLAB R2015a, Simulink. Second, by the well-known Adams-Bashforth method (ode113) the numerical outputs are generated by an experimental input u 2 (t) that is never used in the physical parameter identification in Section 3. Third, the experimental outputs and the corresponding numerical outputs are generated by the same experimental input and are compared with each other. In order to discuss the model validation, the ideal sinusoidal signal r 2 (t) = 100 sin(2π1.0t) (m 2 ) is applied and we cannot provide faster cases any more because not only the cylinders but also the whole hydraulic machinery is vibrated. Note that, in general, a frequency to discuss the model validation should be different from a frequency to identify the physical parameters and also should be higher to investigate the modeling error or uncertainty. Here the initial state values of the nominal numerical simulator are given from the measured values except the initial piston velocityṡ 2 (0), which is not measured but estimated by just first order differentiation. Figure 8 shows the results of the cross validation. The black dots denote the experimental input and outputs from the experimental excavator. The red lines denote the numerical outputs from the nominal numerical simulator. In this section, we compare the experimental outputs and numerical outputs via the fit ratio [29]:

Experimental and Numerical Results and Discussion (Cross Validation)
where y is the output of the experimental excavator,ŷ is the output of the nominal numerical simulator, and E[x] is the average value of x. Table 3 shows the fit ratio. Every ratio was positive, though a fit ratio can be negative very easily in general. As mentioned in Section 3, the pump pressure dynamics are not discussed in the previous paper [20,21]. In Figure 8, successfully, the output of the pump pressure from the experimental excavator and that from the nominal numerical simulator were close to each other. The amplitude errors can be observed but the phase and the number of the peaks were very similar. From a model-based control viewpoint, the phase-delay matching is often more important than the gain matching. Against modeling errors, not only the cap and rod sides pressures p +2 (t) and p −2 (t) but also the pump pressure p p (t) are generated sufficiently. One of the guidelines of the acceptable modeling error is a phase delay condition. From the viewpoint of the passivity [9,11] by which port-Hamiltonian systems can achieve the large region of attraction, the phase delay between the experimental and numerical outputs should be less than 90 degrees. In Figure 8, the maximum phase delay of the model was 33.5 degrees in the rod pressure. In this sense, the modeling error is acceptable.

Energy Behavior Analysis
In this section we discuss the analytical, experimental, and numerical energy behaviors of the excavator on the agricultural scale. First, we review the port-Hamiltonian system and apply the port-Hamiltonian system to the experimental excavator. In particular, we derive the time derivative of the Hamiltonian analytically. Second, the analytical, experimental, and numerical energy behaviors from the port-Hamiltonian point of view in agriculture scale are compared.

Analytical, Experimental, and Numerical Condition (Energy Behavior Analysis)
Let us review the port-Hamiltonian system. One of the most famous representations is the following input-state-output representation [9,10]: where z ∈ R n is the state,ū ∈ R m is the input andȳ ∈ R m is the corresponding output. J(z) = −J T (z) ∈ R n×n is a skew symmetric matrix and R(z) = R T (z) ∈ R n×n is a symmetric matrix. The notation ∇ z is the partial derivative with respect to z. In case of hydraulic cylinder dynamics, z,ū, J(z), R(z), g(z), and H(z) are given as: which is one of the Hamiltonians and not unique due to the existence of the Casimir function [12,13]. For the Hamiltonian in the paper, the first term is the kinetic energy. The second term is the potential energy. The third and fourth terms are the internal energy. This choice of the Hamiltonian is meaningful from a class of model-based control of view [9,10,13,30,31]. The time derivative of the HamiltonianḢ(z) is given as: The dissipation term and the supply rate term in the time derivative of the Hamiltonian −∇ T z HR∇ z H andȳ Tū are as follows: Figure 9 shows the energy analysis procedure. First, we compare the analytical, experimental, and numerical energy behavior of the time derivative of the Hamiltonian. The analytical and numerical energy behaviors of the time derivative of the Hamiltonian are obtained indirectly by Equations (10)- (12), while the experimental energy behavior is derived directly from the measured Hamiltonian and its derivative. Second, we compare the analytical and numerical energy behaviors of the dissipation term and supply rate term. The analytical and numerical energy behaviors of the dissipation term and supply rate term are obtained by Equations (11) and (12), respectively. Finally, we compare the analytical and numerical energy behaviors of the damping term, Coulomb friction term, internal leakage term, and external leakage term in the dissipation term. The analytical and numerical energy behaviors of these four terms are obtained by Equation (11), respectively.
In the paper, we analyze the two cases, that is, the slowest case and the fastest case. In the slowest case, the control input u 2 has the frequency f = 0.5 [Hz]. In the fastest case, the control input u 2 has the frequency f = 1.0 [Hz]. Figure 10 shows the analytical, experimental, and numerical time responses of the time derivative of the Hamiltonian in the slow case and fast case. The analytical, experimental, and numerical energy behaviors were close to each other. One of the guidelines of the acceptable error is a phase delay condition. From the viewpoint of the passivity [9,11] by which port-Hamiltonian systems can achieve the large region of attraction; the phase delay between the analytical, experimental, and numerical outputs should be less than 90 degrees. In Figure 10, the maximum phase delay of the time derivative of Hamiltonian was 57.6 degrees. In this sense, the error is acceptable. One of the most important features is the asymmetric behaviors of the time derivative Hamiltonian. Note that the input u 2 is symmetric in the sense of the standard sinusoidal signal A sin(2π f t). In our study, this asymmetrical behavior comes from the drift of the piston displacement s 2 (t) in Figures 7 and 8. This drift is observed even in the absence of nonlinear frictions [18,32] and now observed experimentally in the presence of the actual frictions. In terms of energy, when the input u 2 is positive, the pump provides the cylinder with little energy to operate the piston, which is stored as internal energy in the fluid, and the piston is operated passively, mainly by converting the stored internal energy in the fluid into kinetic energy. Another important feature is the difference between the slow case and the fast case. This difference comes from the control limitation in the servo open-center valve to overcome the overlapping. In the period u 2 ≈ 0, the overlapping effect in the fast case is stronger than that in the slow case. This is the nature of the servo open-center valve.    Figure 11 shows the numerical and analytical time responses of the dissipation term and the supply rate term of the time derivative of the Hamiltonian. The numerical and analytical time responses of the dissipation term and the supply rate term were close to each other again. However, Figure 11 is valuable because the close behaviors in Figure 10 are just a necessary condition and not a sufficient condition for the close signals in Figure 11. Every response of the time derivative of the Hamiltonian in Figure 10 can be calculated only by the responses in Figure 11 but the opposite is also true. Indeed, it is observed that the main peak in Figure 10 in the fast case comes from the main peak in the supply rate and not from the dissipation. In addition, since the value of the input u 2 is less than 1.0 × 10 −4 and also the values of part exp(p ±2 /b) are almost 1 in our study, the peak in the supply rate comes from the pressure peak eventually. Figure 12 shows the components of the time responses of the dissipation terms (the damping term, the Coulomb friction term, the internal leakage term, and the external leakage term). For every term, the numerical and analytical energy behaviors were close to each other again. Figure 12 is valuable because the close behaviors in the dissipation terms in Figure 11 are just a necessary condition and not a sufficient condition for the close behaviors in the four terms in Figure 12. Every response of dissipation in Figure 11 can be calculated only by the responses in Figure 12 but the opposite is also true. Indeed, except when the peaks in the supply rate are observed, the internal leakage as well as external leakage can be almost zero. These zeros correspond to the fact that these pressures are close to the tank pressure p t in our study. On the other hand, the Coulomb friction term is much smaller than the other three terms and the damping term strongly depends on the piston velocity causing the piston displacement s 2 (t) drift observed in Section 3 and 4. These observations in Figures 11 and 12 support the results in Figure 10 eventually.

Analytical, Experimental, and Numerical Results and Discussion (Energy Behavior Analysis)
In all, via the analytical, experimental, and numerical energy behaviors, the applicability of the port-Hamiltonian theory is confirmed on the agricultural scale against modeling errors.

Conclusions
In this paper, a method for identifying physical parameters of an open-centered multi-degree-of-freedom hydraulic arm is proposed, and the effectiveness of the proposed method is verified by model validation and energy analysis for a one-DOF real construction machine arm. The pump pressure is formulated and a method to uniquely identify the unknown parameters of the open-centered hydraulic arm is given. The agreement between the nonlinear nominal model constructed by the identification and the output signal of the experimental excavator was confirmed under the experimental condition of different frequency from that of the identified input. In other words, we succeed in updating the physical parameter identification method to the pump pressure dynamics. The nominal model [20,21] based on many physical assumptions would be valid only in a limited scale range of time, length, mass and so on. On the other hand, the scale range in industrial excavators (e.g., a mini-shovel mining shovel) is not narrow and our scale is different from mining machines scale or the laboratory machines scale. It would be possible in future work to discuss such other scales. By calculating the output of the Hamiltonian, the validity of updating the model to the port-Hamiltonian system of an open-centered hydraulic arm was evaluated. In all, via the analytical, experimental, and numerical energy behaviors, the applicability of the port-Hamiltonian theory is confirmed in agriculture scale against modeling errors. Future work will be to verify the model and with simultaneous motion of two axes, and to analyze the energy.