Sage Journals: Discover world-class research

Abstract

Uncertainty exists in many industry fields and needs to be dealt properly to avoid unexpected failure. This article proposes a new approach to deal with the uncertain problems encountered by the mathematical modeling of an active hydraulically interconnected suspension system. As the need for both riding comfort and the controllability is soaring nowadays, the traditional passive and semi-active suspension system could barely keep up with the pace, and the proposed active hydraulic system could be one of the solutions. In order to deal with the uncertain factors in the hydraulic system, an interval analysis method for the dynamic responses of nonlinear systems with uncertain-but-bounded parameters using Chebyshev polynomial series is introduced. The comparisons conducted in this article demonstrate the accuracy and computational efficiency of the proposed uncertain problem solver and reveal the influences of uncertain parameters in fluid and mechanical components on the dynamic responses of active hydraulically interconnected suspension.

Keywords

Active suspension hydraulic system vehicle dynamic response interval uncertain analysis Chebyshev series

Introduction

Vehicle has become inalienable in modern society as it extends people’s home range with incredible convenience. As a result, the reliability and the performance of riding attract enormous attention from manufactures, and to some extent decide the popularity of one vehicle. The key in improving the controllability of a vehicle and the riding quality largely depends on suspension which connects a vehicle body to its wheels and isolates the vibration from rough ground,¹ and it also accounts for most fatal crashes, which is around 33%,² related to uncontrollable motions (roll, pitch, bounce, and articulation).³ The situation of rollover, which is the most dangerous, even gets worse with the increasing popularity of sport utility vehicles as they have higher centers of gravity (CGs).⁴ To solve the aforementioned problems, an adequate adjusted suspension system is imperative.⁵

Modern suspension system could be classified into three categories which are passive suspension, semi-active suspension, and active suspension. The prevailing type nowadays goes to passive suspension⁶ giving credit to its cost-effectiveness and reliabilities. However, the drawbacks of this kind of system are also evident, one of which is the compromise between riding comfort and handling stability,⁷ because the increase in suspension roll stiffness and damping will inevitably result in abrupt bumping over rough geography. The other drawback could be illustrated by anti-roll bar, which would undesirably stiffen the suspension warp mode and weaken the road holding ability. As a complement to the passive suspension, semi-active suspension⁸ was brought out. It utilizes the adjustability of the area of orifice in the damper to control the damping force in order to achieve the control over the suspension. But it fails to introduce external energy into the suspension system which differentiates it from the active suspension.^9,10 The most remarkable feature of an active suspension system is that the suspension force is provided,¹¹ at least a portion, by active power sources. It enables the active activities based on the signals acquired by various sensors attached to the key points of the vehicle to measure the motion of the body and suspension, and a control unit will process these data and then adjust the parameters according to the driving condition.¹² There are extensive surveys about the conventional active suspension, and many of them achieved great improvements compared with passive and semi-active ones. However, the conventional active suspensions have their own limitations which constrain the application within a narrow scope such as sport vehicles and luxury vehicles. Typical conventional active suspension has independent structures including four independent controlled actuators, attributing to soaring cost, reduced reliability, increased power consumption, and inherent complexity. The desired control force would be directly applied to the vehicle chassis to achieve superior performance only if the system could satisfy the significant power requirement. Instead of adopting independent layout of four actuators, the idea of introducing compact interconnected fluid circuits into suspension design helps solve the aforesaid problems properly. A good example could be the dynamic ride control (DRC) sport suspension system. The diagonally interconnected mechanical structure adopts a pump as the external pressure source to provide force into the diagonally linked shock absorbers to achieve a stable motion during cornering, avoiding rollover.

To take the progress one step further, a more sophisticated and robust active hydraulically interconnected suspension (HIS) model is proposed;¹³ it cooperates four hydraulically interconnected actuators, which are interconnected by two circuits in roll plane, into each wheel station and a control unit comprised a motor, a bump, a servo valve, a tank, an accumulator, and so on, achieving the goal of structure simplicity and reduced overall cost. It could actively tilt the vehicle against the uncontrollable rollover by providing desired restoring forces,¹⁴ while not consuming much energy, because the variable here is the vertical force, not the suspension deflection, which means not much work should be done.¹⁵

Due to the presence of fluid flow and the companying uncertain parameters such as pipe friction and working fluid damping coefficient,¹⁶ large uncertainty is introduced into the mathematical modeling. The nonlinear nature of the hydraulic system weakens the accuracy of the mathematical model of the actuator that is crucial in analyzing. In this article, to solve the uncertain problems encountered during mathematical modeling, an interval analysis method for dynamic response of nonlinear systems with uncertain-but-bounded parameters is proposed using Chebyshev expansion series.

Uncertainties are inherent in nearly all the real-world problems¹⁷ such as loads, material properties, boundary conditions, fraction, geometry, and the problem of uncertainty especially looms large because of the increase in complexity and precision of modern systems.¹⁸ Unfortunately, even small uncertainties of parameters may attribute to large variation in system dynamic change due to the characteristics of enlargement of propagation.¹⁹ Many methods have been brought out mainly fitting into two categories: probabilistic methods²⁰ and non-probabilistic methods.²¹ The probabilistic methods have acquired decent achievements in practical engineering, but it is limited by the needs to express uncertain parameters as stochastic variables with precise probability distribution, which means to know the complete information.

Interval method,²² as a non-probabilistic method, taking advantages of the ease of acquiring the boundary of the uncertainties, shows great potential in dealing with this problem and has been verified by a range of engineering design problems.²³ The minimal and maximal responses of the uncertain objective and constraints are the corresponding bounds of an interval function; in other words, the interval method calculates the upper bound and lower bound of the true solution. The bounds can be obtained by using the interval arithmetic or optimization method. However, encountered by most interval arithmetic, overestimation,²⁴ which is caused by the so-called “wrapping effect,”²⁵ compromises the accuracy, especially in the process of numerical iterations where it would be accumulated. Although the Taylor interval method^26–28 and Chebyshev interval method are proposed to control overestimation, it is still hard to eliminate the overestimation completely. The optimization method may obtain higher accurate bounds, but its efficiency would be quite low, especially for complicated engineering problems.

In this article, Chebyshev series expansions are introduced to build a surrogate model of the active HIS model. Chebyshev surrogate model achieves a higher numerical accuracy when dealing with uncertain problems. The scanning method, which is more accurate than interval arithmetic, is used to compute the bounds based on the Chebyshev surrogate mode, and it is much more efficient than the scanning method based on original complicated model.

Modeling of active hydraulic interconnected suspension

A representative model is the fundamental of precise and thorough investigation. Due to the characteristic of symmetry of a conventional vehicle chassis, a half-car model is adopted here to simplify the analysis process without the risk of losing generality.

Simplified half-car model

Figure 1 shows the proposed half-car model.

Figure 1.

Half-car model.

As a number of dynamic systems are governed by ordinary differential equations, by using Newton’s second law, the motion equations of this model could be derived as follows

\begin{array}{l} (m h^{2} + I) \overset{\cdot\cdot}{θ} | = - 2 l^{2} k_{s} θ + (Z_{u 1} - Z_{u 2}) l k_{s} - 2 l^{2} C_{s} \dot{θ} + ({\dot{Z}}_{u 1} - {\dot{Z}}_{u 2}) L C_{s} + F_{1} l - F_{2} \\ m_{s} {\overset{\cdot\cdot}{Z}}_{s} | = - k_{s} (Z_{s} + θ l - Z_{u 1}) - C_{s} ({\dot{Z}}_{s} + \dot{θ} l - {\dot{Z}}_{u 1}) - k_{s} (Z_{s} - θ l - Z_{u 2}) - C_{s} ({\dot{Z}}_{s} + \dot{θ} l - {\dot{Z}}_{u 2}) - F_{1} - F_{2} \\ m_{u} {\overset{\cdot\cdot}{Z}}_{u 1} | = - k_{t} Z_{u 1} + k_{s} (Z_{s} + θ l - Z_{u 1}) + C_{s} (Z_{s} + \dot{θ} l - {\dot{Z}}_{u 1}) + F_{1} + k_{t} Z_{g 1} \\ m_{u} {\overset{\cdot\cdot}{Z}}_{u 2} | = - k_{t} Z_{u 2} + k_{s} (Z_{s} + θ l - Z_{u 1}) + C_{s} ({\dot{Z}}_{s} + \dot{θ} l - {\dot{Z}}_{u 1}) + F_{1} + k_{t} Z_{g 2} \end{array}

(1)

where $Z_{g 1}$ and $Z_{g 2}$ are the road surface inputs at wheels, $Z_{s}$ is the vertical displacements of sprung mass, $Z_{u 1}$ and $Z_{u 2}$ are the vertical displacements of the unsprung mass, $F_{1}$ and $F_{2}$ represent the forces exerted on the car body by the cylinders, $I$ is the roll moment of inertia, $l$ is the distance from the CG of the vehicle, $θ$ is the roll angle, $m_{u}$ is the unsprung mass, $m_{s}$ is the sprung mass, $k_{s}$ is the suspension damper coefficient, and $k_{t}$ is the tire stiffness.

Modeling of hydraulic system

There are three main components which comprise the active suspension which are as follows: four double-direction hydraulic actuators, two interconnected hydraulic circuits, and a compact pressure control unit. Figure 2 is the schematic drawing of the active HIS. Hydraulic actuators possess the advantages of low cost and energy efficiency, which are crucial in commercial application. By adopting advanced control strategies, a mechanical hydraulic system could achieve even better performances. However, hydraulic actuators have the inherent drawback which is nonlinear and complicated dynamics. This, to a considerable extent, influences the dynamic response analysis of the whole system.

Figure 2.

Schematic diagram of active HIS.

The four actuators, which are controlled by the pressure control unit, are mounted according to the specific structure of a vehicle between the wheels and the chassis. The active suspension could promptly act to tilt the vehicle to prevent it from rollover. Take the condition in Figure 2 as an example, the vehicle body shows the trend to roll to the left side, the pressure of circuit A would rise, under the command of control unit, to stiffen left side suspension, while the pressure of circuit B would be reduced to keep the body in balance. The anti-moment could be expressed as follows

\begin{matrix} F_{l} = 2 a (P_{A} - P_{B}) \\ F_{r} = 2 a (P_{A} - P_{B}) \\ M_{c} = \frac{F_{l} l}{2} + \frac{F_{r} l}{2} \end{matrix}

(2)

where $P$ is the pressure in the cylinders, $a$ is the piston area, $F$ is the force provided by actuators, $l$ is the distance between the left and the right, and $M$ is the desired anti-moment.

The compact control unit energizes the actuators according to the angle sensors; in other words, it reacts to the demand of the actual road conditions. This feature could also be called demand dependence, which could help the system avoid oversensitive reactions. There are two reasons for avoiding oversensitive reactions; the first one is that “nulling” tilt would consume considerable energy while help little in handling the motion of the vehicle. The other reason is the feeling of acceleration also plays an important role in driving experience,²⁹ which would be compromised by oversensitivity. As a result, the active HIS would help vehicles maintain controllable motion while have comparable efficiency with much lower energy consumption.

Involving the hydraulic system, which possesses high-level uncertainties, makes the numerical analysis more sophisticated. In order to make the expression more concise, the uncertain parameter vector, which is k-dimensional, could be expressed as follows

\underline{x} \leq x \leq \bar{x}

(3)

where $\underline{x} = [{\underline{x}}_{1}, {\underline{x}}_{2}, \dots, {\underline{x}}_{k}]^{T}$ and $\bar{x} = [{\bar{x}}_{1}, {\bar{x}}_{2}, \dots, {\bar{x}}_{k}]^{T}$ denote the lower bounds and upper bounds of the uncertain parameters, respectively. So the interval notation could be defined as

[x] = [\underline{x}, \bar{x}] = {x_{i}, {\bar{x}}_{i} \leq x_{i} \leq {\bar{x}}_{i}, i = 1, 2, \dots, k}

(4)

Figure 3 shows the flow and input sign convention in the HIS.

Figure 3.

The flow and input sign convention in the HIS.

It can be seen that $P_{i}, i = 1, 2, \dots, 7$ are the pressures in the system, $q_{i}, i = 1, 2, \dots, 7$ are the aforesaid flows, and $A_{pu}$ and $A_{pd}$ are the areas of the upside and downside of the piston, respectively. As the force provided by the hydraulic system is directly affected by the areas of the piston, and due to the accuracy of manufacturing, $A_{pu}$ and $A_{pd}$ are considered to be the uncertain parameters, which are denoted as $[A_{pu}]$ and $[A_{pd}]$ . $P_{2}$ , $P_{4}$ , $q_{2}$ , and $q_{4}$ , which are provided by the control unit and controlled by the servo valves, are decided by the restoring force requirement and are given by the following equations

{\overset{\cdot}{P}}_{2} = \frac{E}{L [A_{L}]} (q_{2} - q_{1} - q_{7}) = \frac{E}{L [A_{L}]} (G P_{c 1} - G P_{2} - q_{1} - q_{7})

(5)

{\overset{\cdot}{P}}_{4} = \frac{E}{L [A_{L}]} (q_{4} - q_{3} - q_{5}) = \frac{E}{L [A_{L}]} (G P_{c 2} - G P_{4} - q_{3} - q_{5})

(6)

where $E$ is the bulk modulus which can be expressed as $E = ((- dP) / (dV)) V$ , $V$ is the chamber volume, $L$ is the pipe length from valve to actuators, $[A_{L}]$ is the pipe cross-sectional area, and $G$ is the proportional coefficient of servo valve controller and $P_{ci}$ ; $i = 1, 2$ is the reference pressure value of servo valve.

Take cylinder chamber $P_{1}$ as an example to calculate the pressure

{\dot{P}}_{1} = \frac{E}{s [A_{p d}]} ([A_{p d}] {\dot{x}}_{1} + q_{1}) = \frac{E}{s [A_{p d}]} ([A_{p d}] (- l \dot{θ} + {\dot{Z}}_{s} - {\dot{Z}}_{u 1}) + q_{1})

(7)

where ${\overset{\cdot}{x}}_{1} = - l \overset{\cdot}{θ} + {\overset{\cdot}{Z}}_{s} - {\overset{\cdot}{Z}}_{u 1}$ is the relationship between the piston movement and the vehicle body movement, $s$ is the height of the analyzed oil volume, and $[A_{pd}] {\overset{\cdot}{x}}_{1}$ is the volumetric flow due to piston motion. This equation can be applied to other chambers by substituting the subtitle of the parameters.

In the hydraulic system, the other influential uncertainties come from the pipe friction and the working fluid damping coefficient, and a little change happening in these parameters would lead to large fluctuations in the system dynamic responses. What makes the situation worse is that it is hard to acquire full information about these parameters. Express the uncertain parameters of pipe friction and working fluid damping coefficient with [Re] and [C], respectively, then the flow into the piston chamber could be given as follows

F - [C] q_{1} - [Re] {\overset{\cdot}{q}}_{1} = ma

(8)

(P_{2} - P_{1}) [A_{L}] - [C] q_{1} - [Re] {\overset{\cdot}{q}}_{1} = L [A_{L}] ρ {\overset{\cdot}{q}}_{1}

(9)

{\overset{\cdot}{q}}_{1} (L [A_{L}]) ρ + [Re]) = - P_{1} [A_{L}] + P_{2} [A_{L}] - [C] q_{1}

(10)

{\overset{\cdot}{q}}_{1} = - \frac{P_{1}}{L [A_{L}] ρ + [Re]} + \frac{P_{2}}{L [A_{L}] ρ + [Re]} - \frac{C}{L [A_{L}] ρ + [Re]} q_{1}

(11)

The flows into other chambers could be calculated similarly

{\overset{\cdot}{q}}_{3} = - \frac{P_{3}}{L [A_{L}] ρ + [Re]} + \frac{P_{4}}{L [A_{L}] ρ + [Re]} - \frac{C}{L [A_{L}] ρ + [Re]} q_{3}

(12)

{\overset{\cdot}{q}}_{5} = - \frac{P_{5}}{L [A_{L}] ρ + [Re]} + \frac{P_{4}}{L [A_{L}] ρ + [Re]} - \frac{C}{L [A_{L}] ρ + [Re]} q_{5}

(13)

{\overset{\cdot}{q}}_{7} = - \frac{P_{7}}{L [A_{L}] ρ + [Re]} + \frac{P_{2}}{L [A_{L}] ρ + [Re]} - \frac{C}{L [A_{L}] ρ + [Re]} q_{7}

(14)

Integrated modeling

After the mathematical modeling of the half-car model and the active HIS system, it is feasible to couple these two systems together. The roll angle acquired by the sensors could be used as the input of the active HIS, and the hydraulic forces, which are shown below, provided by the active HIS are the restoring forces for the car model

F_{1} = P_{1} [A_{pd}] - P_{3} [A_{pu}]

(15)

F_{2} = P_{5} [A_{pd}] - P_{7} [A_{pu}]

(16)

In order to obtain the assembled system, on the fundamental of the coupling interactions between the half-car model and the active HIS system, the equations concerning these two models should be combined together and be expressed in matrix form. The state vectors could be like

x_{HIS} = [\begin{matrix} P_{1} & P_{3} & P_{5} & P_{7} & q_{1} & q_{3} & q_{5} & q_{7} & P_{2} & P_{4} \end{matrix}]^{T}

(17)

x_{car} = [\begin{matrix} θ & Z_{s} & Z_{u 1} & Z_{u 2} & \overset{\cdot}{θ} & {\overset{\cdot}{Z}}_{s} & {\overset{\cdot}{Z}}_{u 1} & {\overset{\cdot}{Z}}_{u 2} \end{matrix}]^{T}

(18)

The state of the HIS system could be expressed like

x = [\begin{matrix} x_{HIS} & x_{car} \end{matrix}]^{T}

(19)

u = [\begin{matrix} u_{HIS} & u_{car} \end{matrix}]^{T}

(20)

A = [\begin{matrix} A_{HIS} & A_{HIS - car} \\ A_{car - HIS} & A_{car} \end{matrix}]

(21)

B = [\begin{matrix} A_{HIS} & 0 \\ 0 & A_{car} \end{matrix}]

(22)

where $A_{HIS - car}$ and $A_{car - HIS}$ are the coupling between the half-car and half-HIS systems, respectively. $A_{HIS - car}$ converts the pressure from the hydraulic system into the active force to vehicle system and $A_{car - HIS}$ transfers the vehicle body displacement to the pressure change in the hydraulic system

A_{HIS - car} = [\begin{matrix} \begin{matrix} gd \\ - hd \end{matrix} & \begin{matrix} - gu \\ hu \end{matrix} & \begin{matrix} \begin{matrix} - gd & gu \end{matrix} \\ \begin{matrix} - hd & hu \end{matrix} \end{matrix} \\ id & - iu & \begin{matrix} 0 & 0 \end{matrix} \\ 0 & 0 & \begin{matrix} id & - iu \end{matrix} \end{matrix}]

(23)

A_{car - HIS} = [\begin{matrix} \begin{matrix} - ee \\ ee \end{matrix} & \begin{matrix} ff \\ - ff \end{matrix} & \begin{matrix} \begin{matrix} - ff & 0 \end{matrix} \\ \begin{matrix} ff & 0 \end{matrix} \end{matrix} \\ ee & ff & \begin{matrix} 0 & - ff \end{matrix} \\ - ee & - ff & \begin{matrix} 0 & - iu \end{matrix} \end{matrix}]

(24)

where $gd = ([A_{pd}]^{*} l) / ({ms}^{*} h^{2} + I)$ , $gu = ([A_{pu}]^{*} l) / ({ms}^{*} h^{2} + I)$ , $hd = [A_{pd}] / ms$ , $hu = [A_{pu}] / ms$ , $id = [A_{pd}] / mu$ , $iu = [A_{pu}] / ms$ , $ee = (E^{*} l) / s$ , and $ff = E / s$ . $l$ is the distance from CG to a suspension spring/damper, $h$ is the height of the CG, $ms$ is the sprung mass, and $mu$ is the unsprung mass.

The complete model of the proposed system in the form of

\begin{matrix} \overset{\cdot}{x} = Ax + Bu \\ [\begin{matrix} {\overset{\cdot}{x}}_{HIS} \\ {\overset{\cdot}{x}}_{car} \end{matrix}] = [\begin{matrix} A_{HIS} & A_{HIS - car} \\ A_{car - HIS} & A_{car} \end{matrix}] [\begin{matrix} x_{HIS} \\ x_{car} \end{matrix}] \\ + [\begin{matrix} A_{HIS} & 0 \\ 0 & A_{car} \end{matrix}] [\begin{matrix} u_{HIS} \\ u_{car} \end{matrix}] \end{matrix}

(25)

The proposed active HIS-assisted half-car model is developed with Simulink.

Chebyshev surrogate model

As demonstrated in the previous sections, the high uncertainties introduced by the hydraulic system should be adequately dealt with to develop a robust and reliable analysis model. In engineering, the uncertainties induced by bounded parameters could be handled by convex model or interval model. Traditionally, low-level Taylor series–based interval method can only solve the problems whose uncertain levels are small; otherwise, the overestimation phenomenon would be frustrating. In this section, Chebyshev series expansion³⁰ is introduced to build a surrogate model of the hydraulic interconnected suspension system, and then it is combined with scanning method to compute the bounds.

Considering the responses of the model as a continuous function without knowing its analytical expression, let us assume that there exists a polynomial $p (x)$ which fits the criterion that it converges to $f (x)$ on $[a, b]$

‖ f (x) - p (x) ‖_{\infty} < ε, x \in [a, b]

(26)

where $f (x)$ is continuous over $[a, b]$ , and this expression holds for any $ε > 0$ .

Let $P_{n} (x)$ denote the set of polynomials whose degree is not bigger than $n$ , $n \in N^{0}$ , a unique polynomial $p_{n}^{*}$ exists in $P_{n} (x)$

‖ f (x) - p (x) ‖_{\infty} \geq ‖ f (x) - p_{n}^{*} ‖_{\infty} = E_{n} (f)

(27)

where $x \in [a, b]$ .

$E_{n} (f)$ is the infimum of maximum error, which holds for all $p (x) \in P_{n} (x)$ other than $p_{n}^{*}$ . It is defined as

E_{n} (f) = inf_{p_{n} (x) \in P_{n} (x)} ‖ f (x) - p_{n}^{*} ‖_{\infty}

(28)

where $x \in [a, b]$ . In this equation, $p_{n}^{*}$ is the best uniform approximation of degree $n$ to $f (x)$ on $[a, b]$ . But it is not feasible to calculate $p_{n}^{*}$ as the computation cost rises considerably when $n > 2$ . Fortunately, the truncate Chebyshev series could be used to approximate the original function, which could achieve a desirable result, because it is very close to the best uniform approximation polynomials.

The Chebyshev polynomial could be expressed by $T_{n}$ for $x \in [a, b]$ of degree $n$

T_{n} (x) = \cos n θ, where θ = a c r \cos (\frac{2 x - (b + a)}{b - a}) \in [0, π]

(29)

For the multi-dimensional problem, the concept on tensor product should be adopted. By generating the tensor products of the one-dimensional equations, the Chebyshev polynomial for a k-dimensional problem is expressed as

T_{n_{1}, n_{2}, \dots, n_{k}} (x_{1}, x_{2}, \dots, x_{k}) = \cos (n_{1} θ_{1}) \cos (n_{2} θ_{2}) \dots \cos (n_{k} θ_{k})

(30)

where $θ_{i} = acr \cos (x_{i})$ .

Because of the orthogonality of this series, the following equation holds

\int_{- 1}^{1} \frac{T_{n} (x) T_{m} (x)}{\sqrt{1 - x^{2}}} d x = \int_{0}^{π} (\cos n θ \cos m θ) = {\begin{matrix} π, & m = n = 0 \\ \frac{π}{2}, & m = n \neq 0 \\ 0, & m \neq n \end{matrix}

(31)

where $1 / (\sqrt{1 - x^{2}})$ is the weighting function.

The truncated Chebyshev series expansion of $f (x)$ which is continuous on $[a, b]$ could be expressed as

f (x) \approx p_{n} (x) = \frac{1}{2} f_{0} + \sum_{i = 1}^{n} f_{i} T_{i} (x)

(32)

where $f_{i}$ are the constant coefficients.

The truncated error³¹ is

e_{n} (f) = | f (x) - p_{n} (x) | \leq \frac{2^{- n}}{(n + 1)!} ‖ f^{(n + 1)} ‖_{\infty}

(33)

As to the calculation of the constant coefficients $f_{i}$ , it can be derived using the following equation

f_{i} = \frac{2}{π} \int_{- 1}^{1} \frac{f (x) T_{i} (x)}{\sqrt{1 - x^{2}}} dx = \frac{2}{π} \int_{0}^{π} (f (\cos n θ) \cos i θ) d θ

(34)

where $i = 0, 1, 2, \dots, n$ . Similarly, the concept of tensor product is used to generate the k-dimensional coefficients

f_{i_{1}} f_{i_{2}}, \dots, f_{i_{n}} = {(\frac{2}{π})}^{k} \int_{0}^{π} \dots \int_{0}^{π} (f (\cos θ_{1}, \dots, \cos θ_{k}) \cos i_{1} θ_{1}, \dots, \cos i_{k} θ_{k}) d θ_{1} \dots d θ_{k}

(35)

where $i_{1}, \dots, i_{k} = 0, 1, 2, \dots, n$ . In order to solve equations (34) and (35), Gaussian–Chebyshev³² integral is adopted; the integral formula of the above equation can be expressed as

\int_{- 1}^{1} \frac{1}{\sqrt{1 - x^{2}}} f (x) dx \approx \frac{π}{m} \sum_{j = 1}^{m} f (x_{j}) = \frac{π}{m} \sum_{j = 1}^{m} f (\cos θ_{j})

(36)

So the coefficients of Chebyshev polynomials are as follows

\begin{array}{l} f_{i} = \frac{2}{π} \int_{- 1}^{1} \frac{f (x) T_{i} (x)}{\sqrt{1 - x^{2}}} d x \approx \frac{2}{π} \frac{π}{m} \sum_{j = 1}^{m} f (x_{j}) T_{i} (x_{j}) \\ = \frac{2}{m} \sum_{j = 1}^{m} f (\cos θ_{j}) \cos i θ_{j} \end{array}

(37)

Equation (37) is a linear combination of the values of the function which makes the calculation Chebyshev coefficient easy to obtain.

For a k-dimension problem

\begin{matrix} f_{i_{1}, \dots, i_{k}} = {(\frac{2}{π})}^{k} \int_{- 1}^{1} \dots \int_{- 1}^{1} \frac{f (x) T_{i_{1}, \dots, i_{k}} (x)}{\sqrt{1 - {x_{1}}^{2}} \dots \sqrt{1 - {x_{k}}^{2}}} d x_{1} \dots d x_{k} \\ \approx {(\frac{2}{m})}^{k} \sum_{j_{1} = 1}^{m} \dots \sum_{j_{k} = 1}^{m} f (\cos θ_{j_{1}} \dots \cos θ_{j_{k}}) \cos i_{1} θ_{j_{1}} \dots \cos i_{1} θ_{j_{k}} \end{matrix}

(38)

In order to minimize the integral error, $m$ is usually set not less than $n + 1$ . After obtaining equation (32), we have built a surrogate model of original complex model, so the scanning method can be used to compute the bounds of responses.

The algorithm can be summarized as:

Input: n, m,

x = [a, b] \subset I R^{k}

1. Produce the interpolation points

θ_{j} = \frac{(2 j - 1)}{2 m}

;

x_{j_{i}} = \frac{a_{i} + b_{i}}{2} + \frac{b_{i} - a_{i}}{2} \cos θ_{j}

;

j = 1, 2, \dots, m, i = 1, 2, \dots, k

2. Run the Simulink model to compute

f (t, \cos θ_{j_{1}}, \dots, \cos θ_{j_{k}})

3. Compute the coefficient

f_{i_{1}, \dots, i_{k}}

of Chebyshev polynomials by using equation (38).

4. Construct the Chebyshev surrogate model

f (t, x) = \sum_{i_{1} = 0}^{n} \dots \sum_{k = 0}^{n} {(\frac{1}{2})}^{p} f_{i_{1}, \dots, i_{k}} \cos i_{1} θ_{j_{1}}, \dots, \cos i_{k} θ_{j_{k}}; θ \in [0, π]^{k}

5. Using the scanning method to compute the bounds

[\underline{f}, \bar{f}]

based on the surrogate model.

Output

[f (t, x)]

Interval uncertain analysis of active HIS

In this section, the proposed Chebyshev interval method is used in the active HIS system to deal with the uncertain parameters aforesaid to develop a robust and reliable model.

Parameter tuning

In order to guarantee the proposed model could adequately represent the real situation, parameter tuning should be conducted at first. Table 1 shows the parameters adopted in this article.

Table 1.

Integrated system parameters.

Parameters	Value	Parameters	Value
I: roll moment of inertia ( $kg \times m^{2}$ )	388	s: initial fluid height in the piston (m)	0.075
ms: sprung mass (kg)	950	G: valve P controller gain	$1.5 \times 10^{- 9}$
mu: unsprung mass (kg)	41	$[A_{pu}]$ : upper piston cross-sectional area (m²)	$[5.049 \times 10^{- 4}, 5.151 \times 10^{- 4}]$
l: distance from CG to a suspension spring/damper (m)	0.5	$[A_{pd}]$ : lower piston cross-sectional area (m²)	$[4.059 \times 10^{- 4}, 4.141 \times 10^{- 4}]$
k_t: tire stiffness (N/m)	$2.7 \times 105$	L: pipe length from valve to actuator (m²)	1.0
k_s: suspension damper coefficient (N s/m)	$44 \times 103$	$ρ$ : hydraulic oil density ( $kg / m^{3}$ )	870
E: bulk modulus (MPa)	1400	[C]: working fluid damping coefficient	$[4.5 \times 10^{6}, 5.5 \times 10^{6}]$
$[A_{L}]$ : pipe cross-sectional area	$[3.9196 \times 10^{- 5}, 3.9988 \times 10^{- 5}]$	$[Re]$ : pipe wall friction coefficient	$[0.9 \times 10^{4}, 1.1 \times 10^{4}]$

CG: center of gravity.

Most of the parameters are obtained through previous work,¹⁵ such as the suspension spring stiffness, the suspension damper coefficient, and the valve controller gain. In this article, the main concentration is focused on the uncertain parameters whose little change would introduce large fluctuation in the system dynamic responses. These parameters are separated into two sets according to their nature, which is hydraulic or not, for the sake of adequate demonstration. The first set includes upper piston cross-sectional area (m²), which can be expressed as $[A_{pu}]$ , the lower piston cross-sectional area as $[A_{pd}]$ (m²), and the pipe cross-sectional area $[A_{L}] (m^{2})$ . These parameters could directly influence the restoring force provided by the active HIS because the force is the product of hydraulic pressure and the area. Most importantly, these parameters could be different from the rating due to the accuracy of manufacturing and the working conditions. In this specific experiment, the range is set to 1%. As to the second set, they are the pipe wall friction coefficient [Re] and the working fluid damping coefficient [C]. Hydraulic characteristics are heavily dependent on these two parameters as they could influence the pressures’ rise time and the delivered pressure. What makes the situation worse is that these two parameters have large uncertain level and are difficult, if possible, to have an accurate value.

Uncertain analysis of fluid parameters

The anti-roll mode would be activated when the vehicle experiences rolling under the operation of turning. In this simulation, the vehicle is supposed to have a right turn; as a result, it would have a tendency to roll to the left side. In order to counteract this tendency, the HIS system would rise the pressure of the upper chamber on the left side and reduce the pressure on the right side to provide the entailed supportive force on the left. To simulate this action, a 25-bar signal and a 50-bar signal are sent to the lower chamber circle and upper chamber circle on the left side through the servo valve, respectively. The order of the Chebyshev inclusion function is set to $n = 2$ , so $m = 3$ . As a comparison to evaluate the effectiveness and the accuracy of the proposed Chebyshev interval method, scanning method, with 11 symmetrical sampling points is adopted. Figure 5 shows the dynamic responses of deflections on both sides with the uncertain parameters of [C] and [Re].

It can be seen in Figure 4(a) that the deflection experiences an oscillation at the beginning of the process and then a steep ascent until the predetermined value. It is evident that there exists a considerable area between the upper bounds and the lower bounds, which is caused by the uncertain parameters of [C] and [Re] even under the condition that the uncertain range was set to a relative low level. If the uncertainty rises, the differences between these two bounds would be even greater. As to the upper bound, it takes about 800 ms to reach the stable state which is expressed by the horizontal trend; the lower bound takes around 1200 ms to reach a relative stable state. It should be point out that there are actually two upper bounds and two lower bounds in Figure 4(a) which can only be seen in the amplified plot Figure 4(b) because the differences between them are very small. In Figure 4(b), the upper bound and the lower bound which enclose the other two bounds are generated by the proposed Chebyshev inclusion method, and the enclosed two bounds are the results of the comparison set named scanning method, which represent the relatively accurate real situation. The little difference between these two methods demonstrates the good approximation ability that the Chebyshev surrogate model possesses. Another merit that lies in the proposed method is the computational efficiency. The iterative numbers of the Chebyshev method and the scanning method are $m^{2} = 9$ and $11^{2} = 121$ . Achieving the great reduction in iteration times without losing calculation accuracy makes the proposed method more applicable and possesses the generalization capability. In Figure 4(c), the deflections are much the same as Figure 4(a) shows but in the opposite direction, and the oscillation in Figure 4(c) is more obvious than that of Figure 4(a). The amplified plot shows the details of the fluctuation, and it can be seen that the Chebyshev inclusion method also performs well to enclose the scanning method trajectory.

Figure 4.

Deflection responses of both sides with uncertain parameters of [C] and [Re]: (a) deflection on left side, (b) amplified plot on left side, (c) deflection on right side, and (d) amplified plot on right side.

Figure 5 shows the dynamic pressure responses in the four chambers of the two actuators. As to the left actuator, the stable pressure of the upper chamber is 25 bar and the stable pressure of the lower chamber is 50 bar. The pressure difference between these two chambers makes the actuator to extend or, in other words, to stiffen the suspension. To the right side, it is exactly the opposite situation. It can also be seen that, to the upper chambers, the transient pressure would reach around 50 bar and then decrease about 10 bar. The pressure of the lower chambers, at first, would only reach around 40 bar. The oscillation during the establishment of the stable pressure is obvious, which is the characteristics of hydraulic. The influences introduced by the uncertain parameters are much the same as the deflection situation demonstrated; the upper bounds and the lower bounds draw the possible region where the real pressure response trajectory might be. The effectiveness of the proposed Chebyshev internal method is illustrated in Figure 6.

Figure 5.

Pressure responses of the four chambers of two actuators with uncertain parameters of [C] and [Re]: (a) left upper chamber, (b) left lower chamber, (c) right upper chamber, and (d) right lower chamber.

Figure 6.

Amplified pressure responses of the four chambers of two actuators with uncertain parameters of [C] and [Re]: (a) amplified left upper chamber, (b) amplified left lower chamber, (c) amplified right upper chamber, and (d) amplified right lower chamber.

In Figure 6, four amplified plots about the detail comparison between the scanning method and the Chebyshev method are shown. During the fluctuation, the proposed Chebyshev method shows great capability of tracking the real value provided by the scanning method, even at the sharp corners of the trajectory.

Figure 7 is the car roll responses with uncertain parameters of [C] and [Re]. It is worth pointing out that the oscillations happened in the responses of deflections and pressures are absent in this situation, which means the passengers in the vehicle would feel the undesirable fluctuations. The reason is the vehicle’s original suspension and joints absorb the unstable energy and filter the oscillations.

Figure 7.

Car roll responses with uncertain parameters of [C] and [Re]: (a) car roll and (b) amplified plot of car roll trajectory.

Uncertain analysis of mechanical components parameters

As addressed before, the analysis conducted above are related to uncertain parameters of [C] and [Re]; another set of uncertain parameters in the proposed active HIS system includes $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ .Therefore, the Chebyshev interval method runs the Simulink model $m^{3} = 27$ times and the scanning method runs $11^{3} = 1331$ times. This set of parameters is concerned with the manufacturing accuracy and the working conditions such as pressure change and thermal unstable, which is very meaningful and could be found everywhere in the industrial fields. In order to make this article more concise and more heuristic, only the comparisons of the differences between the influences of different sets of uncertain parameters are presented. Figure 8 shows the influence comparison between the two sets of uncertain parameters in deflections.

Figure 8.

Influence comparison between the two sets of uncertain parameters in deflections: (a) the deflection on right side and (b) the deflection on left side.

In Figure 8, it can be seen that at the beginning, the influence of [C] and [Re] is larger than that of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ , which is illustrated by the discrepancy between the bounds. But as the process goes on, the influence concerning $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ becomes bigger and surpass the influence of [C] and [Re] at around 500 ms. Here, an interesting phenomenon happens; it can be seen that the influence of the introduced uncertainty is different to these two parameter sets. To [C] and [Re], the uncertainty only affects the deflection process, but not the stable value; in other words, the two bounds converge to a signal value in the end. As to $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ , the uncertainty affects not only the responding process but also the final stable value. It can be seen that the upper bounds and the lower bounds ended in a great discrepancy. The reason for this phenomenon is that the pressure here is constant as it is set, the pipe wall friction coefficient and the fluid damping coefficient only delay the pressure rise time, and the deflection will eventually reach the same value as it is determined by the supporting force. As to the $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ , although the pressure is the same, the restoring force is jointly decided by the pressure and the active area, which results in the final differences.

Figure 9 shows the pressure response in the two actuators under the two sets of uncertain parameters. Different from the final stable discrepancies in the deflection responses with $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ , responses with two sets of uncertain parameters all converged to the same pressures which are predetermined by the input signal. Although it converges to the same value, the response uncertain area or interval area of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ is obvious larger than that of [C] and [Re] which reveal the fact that the influence of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ is greater than that of [C] and [Re].

Figure 9.

Influence comparison between the two sets of uncertain parameters in pressures: (a) left upper chamber, (b) left lower chamber, (c) right upper chamber, and (b) right lower chamber.

Figure 10 illustrates the car roll responses under the two sets of uncertain parameters. The same with suspension deflection; although the pipe wall friction coefficient and the fluid damping coefficient are the interval parameters, the stable car roll of the two bounds is the same. On the contrary, the stable position of the car roll under the uncertain parameters of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ is different. One particular point which is worth pointing out is that the uncertain range of [C] and [Re] in this demonstration is set to 10%, while the uncertain range of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ is set to 1%. Although the uncertain range of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ is only 10% of the range of [C] and [Re], the influence of $[A_{pu}]$ , $[A_{pd}]$ , and $[A_{L}]$ is still larger. This result emphasizes the importance of the control of the active area.

Figure 10.

Influence comparison between the two sets of uncertain parameters in car roll.

Conclusion

This article presents a new approach regarding the mathematical modeling of an active HIS with uncertain parameters by introducing the Chebyshev inclusion method. The active HIS system overcomes the drawbacks encountered by traditional passive suspension system semi-active system of lacking flexibility by adopting active hydraulic control unit and overweighs the conventional active suspension system by employing interconnection structures. What accompanies the introducing of hydraulic system in the mathematical model is the problem of uncertainty. Uncertainty would considerably degrade the veracity of a mathematical model and even make the model invalid, if not properly handled. As the pipe friction coefficient and the fluid damping coefficient are of high-level uncertainty and are hard to acquire relative accurate values, taking the advantage of the concept of interval method, a non-probabilistic method called Chebyshev interval model is proposed. Although these parameters are uncertain, the boundaries of them are relative easy to acquire.

Two sets of uncertain parameters are chosen in this article, representing the uncertainty of hydraulic systems and the uncertainty of manufacturing accuracy, respectively. To evaluate the effectiveness and accuracy of the proposed method, a thorough comparison has been conducted. Three points have been revealed and verified by the demonstration.

The presented Chebyshev interval method could solve the uncertain problem with similar accuracy but higher efficiency compared with the standard scanning method.

The possible working area which is enclosed by the upper bound and lower bound sheds light upon the necessity to conduct uncertain analysis simulation prediction with a predetermined uncertain parameter may contradict the real situation.

Footnotes

Academic Editor: Weicha Sun

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This research was partially supported by the Australian Research Council (Discovery Projects) (DP150102751) and the National Natural Science Foundation of China (11472112).

References

Cao

Liu

. State of the art in vehicle active suspension adaptive control systems based on intelligent methodologies. IEEE T Intell Transp 2008; 9: 392–405.

National Highway Traffic Safety Administration (NHTSA). Initiative to address the mitigation of vehicle rollover. Docket No. NHTSA-2003-14622-1, June2003. Washington, DC: NHTSA.

Cao

Rakheja

C-Y.

Roll- and pitch-plane coupled hydro-pneumatic suspension: part 1—feasibility analysis and suspension properties. Vehicle Syst Dyn 2010; 48: 361–386.

National Highway Traffic Safety Administration (NHTSA). Traffic safety facts 1997. Washington, DC: Department of Transportation, NHTSA, 1998.

Zhu

Wang

Zhang

. H_∞ control of a novel low-cost roll-plane active hydraulically interconnected suspension: an experimental investigation of roll control underground excitation. SAE Int J Passeng Cars Mech Syst 2013; 6: 882–893.

Khajavi

Abdollahi

Comparison between optimized passive vehicle suspension system and semi active fuzzy logic controlled suspension system regarding ride and handling. In: Proceedings of the world academy of science, engineering and technology, vol. 1, 2007, pp.57–61. International Journal of Mechanical, Aerospace, Industrial, Mechatronic and Manufacturing Engineering.

Mantaras

Luque

Ride comfort performance of different active suspension systems. Int J Vehicle Des 2005; 40: 106–125.

Ivers

Miller

LR.

Semi-active suspension technology: an evolutionary view, vol. 40. New York: ASME, 1991, pp.327–346.

Iijima

Akatsu

Takahashi

. Development of a hydraulic active suspension. SAE technical paper 931971, 1993.

10.

Yokoya

Kizu

Kawaguchi

. Integrated control system between active control suspension and four wheel steering for the 1989 Celica. SAE technical paper 901748, 1990.

11.

Hrovat

Survey of advanced suspension developments and related optimal control applications. Automatica 1997; 33: 1781–1817.

12.

Karkoub

Zribi

Active/semi-active suspension control using magnetorheological actuators. Int J Syst Sci 2006; 37: 35–44.

13.

Smith

Zhang

Jeyakumaran

. Ride simulations of a half-car with a hydraulically interconnected passive suspension. In: Proceedings of the FISITA-2006 world automotive congress, Yokohama, Japan, 22–27 October 2006, pp.1–10.

14.

Sun

Gao

Kaynak

Finite frequency control for vehicle active suspension systems. IEEE T Contr Syst T 2011; 19: 416–422.

15.

Zhang

Wang

Modeling of a new active suspension for roll control, 2009, https://opus.lib.uts.edu.au/handle/10453/11092

16.

Wang

Zhang

. Experimental investigation of a hydraulically interconnected suspension in vehicle dynamics and stability control. SAE technical paper 2012-01-0240, 2012.

17.

Schuëller

Jensen

HA.

Computational methods in optimization considering uncertainties—an overview. Comput Method Appl M 2008; 198: 2–13.

18.

Hanss

The transformation method for the simulation and analysis of systems with uncertain parameters. Fuzzy Set Syst 2002; 130: 277–289.

19.

Cheng

Sandu

Uncertainty quantification and apportionment in air quality models using the polynomial chaos method. Environ Modell Softw 2009; 24: 917–925.

20.

Chen

Jin

Sudjianto

Analytical variance-based global sensitivity analysis in simulation-based design under uncertainty. J Mech Des: T ASME 2005; 127: 875–886.

21.

Moens

Vandepitte

A survey of non-probabilistic uncertainty treatment in finite element analysis. Comput Method Appl M 2005; 194: 1527–1555.

22.

Ishibuchi

Tanaka

Multiobjective programming in optimization of the interval objective function. Eur J Oper Res 1990; 48: 219–225.

23.

Jiang

Xie

. Interval uncertain optimization of vehicle suspension for ride comfort. CMES: Comp Model Eng 2014; 98: 443–467.

24.

Jackson

Nedialkov

NS.

Some recent advances in validated methods for IVPs for ODEs. Appl Numer Math 2002; 42: 269–284.

25.

Lin

Stadtherr

MA.

Validated solutions of initial value problems for parametric ODEs. Appl Numer Math 2007; 57: 1145–1162.

26.

Qiu

Wang

Parameter perturbation method for dynamic responses of structures with uncertain-but-bounded parameters based on interval analysis. Int J Solids Struct 2005; 42: 4958–4970.

27.

Qiu

Two non-probabilistic set-theoretical models to predict the transient vibrations of cross-ply plates with uncertainty. Appl Math Model 2008; 32: 2872–2887.

28.

Qiu

Non-probabilistic convex models and interval analysis method for dynamic response of a beam with bounded uncertainty. Appl Math Model 2010; 34: 725–734.

29.

Wang

Shen

Integrated vehicle ride and roll control via active suspensions. Vehicle Syst Dyn 2008; 46: 495–508.

30.

Luo

Zhang

. A new interval uncertain optimization method for structures using Chebyshev surrogate models. Comput Struct 2015; 146: 185–196.

31.

Jiang

Zhao

Shu

Numerical approximation. Shanghai, China: Fudan University Press, 2008.

32.

Abramowitz

Stegun

IA.

Handbook of mathematical functions: with formulas, graphs, and mathematical tables. Chelmsford, MA: Courier Corporation, 1964.

Interval uncertain analysis of active hydraulically interconnected suspension system

Abstract

Keywords

Introduction

Modeling of active hydraulic interconnected suspension

Simplified half-car model

Modeling of hydraulic system

Integrated modeling

Chebyshev surrogate model

Interval uncertain analysis of active HIS

Parameter tuning

Uncertain analysis of fluid parameters

Uncertain analysis of mechanical components parameters

Conclusion

Footnotes

Declaration of conflicting interests

Funding

References