• Introduction
• Loading the AIsys package
• Example 4: Hydraulic System

Hydraulic system - Example of a Non-Linear Affine-Input System

by   A. Kugi, K. Schlacher and R. Novak, Department of Automatic Control, Johannes Kepler University, Linz, Austria,

kugi@mechatronik.uni-linz.ac.at

NOTE: This worksheet demonstrates the use of the AIsys package to solve the problem of a Hydraulic System.

Introduction

The AIsys package contains several algorithms for the analysis and synthesis of nonlinear Affine-Input control systems, or AI-systems for short. This class of systems describes a large number of physical systems in many engineering applications. However, only the increasing availability of low cost digital signal processors in combination with the increasing power of computer algebra programs like Maple enable the practical use of these nonlinear control strategies.  

The AIsys package contains four examples which prove that the proposed algorithms are also feasible for practical problems. This mechanical example is the ball and wheel experiment. Other examples in the Maple Application Center include benchmark tests for nonlinear control systems, namely the chemical stirrer vessel, the squirrel cage induction motor and a hydraulic system.

Loading the AIsys package

The AIsys package has been archived into a Maple Library. To install it download the file AIsys2.zip and extract the contents into a directory.

Note: Do not extract these files into your Maple "lib" directory since you could overwrite your main Maple library.

>

restart;

>	libname := "C:/mylib/aisys", libname;

>	with(AIsys);

To find out more read this worksheet or the help file:

>

#?AIsys

Example 4: Hydraulic System

                                          

The above figure presents the schematic diagram for a hydraulic ram. For more details of the plant and an industrial application, see also [Kugi2000] A. Kugi, Non-linear Control Based on Physical Models, Lecture Notes in Control and Information Sciences 260, Springer-Verlag, 2000. The continuity equations of this system read as

>	diff(P1(t),t) = Eoil (Q1 - vp*A1 - Cleak_int*(P1-P2)-Cleak_ext1*P1)/(V01+xp*A1); diff(P2(t),t) = Eoil (Q2 + vp*A2 + Cleak_int*(P1-P2)-Cleak_ext2*P2)/(V02-xp*A2); diff(xp(t),t) = vp; diff(vp(t),t) = (P1*A1 - P2*A2)/m;

Here P1, P2 denote the pressures in the chambers, xp is the piston displacement measured from the center, vp is the piston velocity. Additional parameters are the cylinder volumes V01, V02 (measured at xp = 0), the piston area A, the effective bulk modulus Eoil and the valve flows Q1 and Q2. The supply pressure is given by Ps, the tank pressure by Pt, the mass of the moving parts is m.  The external and internal leakage flows are assumed to be laminar with the leakage coefficients Cleak_ext1, Cleak_ext2 and Cleak_int. The valve flows for positive and negative spool displacements xs are given by the equations below, Kd is a valve parameter.

If xs > 0 the flows are:
Q1_pos := Kd*xs*sqrt(Ps-P1);
Q2_pos := -Kd*xs*sqrt(P2-Pt);
If xs < 0 the flows are
Q1_neg := Kd*xs*sqrt(P1-Pt);
Q2_neg := -Kd*xs*sqrt(Ps-P2);
If xs = 0 the flows are
Q1_zero := 0;
Q2_zero := 0;

Two different models for positive and negative valve displacements are necessary, the plant input is the valve position xs and the plant output is the piston position xp.

>

f:=[Eoil/(V01+xp*A1)*(-vp*A1- Cleak_int*(P1-P2)-Cleak_ext1*P1),Eoil/(V02-xp*A2)*(vp*A2+ Cleak_int*(P1-P2)-Cleak_ext2*P1),vp,(P1*A1-P2*A2)]; g_pos:=[Eoil/(V01+xp*A1)*Kd*sqrt(Ps-P1),-Eoil/(V02-xp*A2)*Kd*sqrt(P2-Pt),0,0]; g_neg:=[Eoil/(V01+xp*A1)*Kd*sqrt(P1-Pt),-Eoil/(V02-xp*A2)*Kd*sqrt(Ps-P2),0,0]; xx:=[P1,P2,xp,vp];

Often in hydraulic systems the assumption Ps-P1 = P2 - Pt is made, in this case the flow Q = Q1 = Q2 can be used as input and the vector field g becomes

>	g_Q:=[Eoil/(V01+xp*A1),-Eoil/(V02-xp*A2),0,0];

>	`AIsys/Global/simp`  := proc(x)  simplify(x); end:

Setup the mathematical model for the AIsys package

>	f:=M_AIsys(f); g_pos:=M_AIsys(g_pos); g_neg:=M_AIsys(g_neg); g_Q:=M_AIsys(g_Q);

Computation of the reachability distribution, to check if the systems hangs up we want more outputs

>	infolevel[AIsys]:=2:

>	print(`For positive valve displacement`); RDistr(f,g_pos,xx):

`Initialization: delta:={}`

`delta := delta + g1`

`Computing Lie bracket: [f,g1]`

`delta := delta + [f,g1]`

`Computing Lie bracket: [delta2,delta1]`

`delta := delta + [delta2,delta1]`

`Computing Lie bracket: [f,delta2]`

`delta := delta + [f,delta2]`

`Dimension of reachability distribution:`, 4

More outputs can be obtained by setting infolevel[AIsys]:=3: or infolevel[AIsys]:=4:

What happens if the flow Q is used as input?

>	RDistr(f,g_Q,xx):

`Initialization: delta:={}`

`delta := delta + g1`

`Computing Lie bracket: [f,g1]`

`delta := delta + [f,g1]`

`Computing Lie bracket: [delta2,delta1]`

`Computed Lie bracket is linear dependent.`

`Computing Lie bracket: [f,delta2]`

`delta := delta + [f,delta2]`

`Computing Lie bracket: [delta3,delta1]`

`Computed Lie bracket is linear dependent.`

`Computing Lie bracket: [delta3,delta2]`

`delta := delta + [delta3,delta2]`

`Dimension of reachability distribution:`, 4

Is the system still reachable without any leakages? To answer this question the side-relations in the simplifiers routines are adjusted in the following way:

>	`AIsys/SIMP`[2]:={Cleak_ext1=0,Cleak_ext2=0,Cleak_int=0};

>	RDistr(f,g_Q,xx):

`Initialization: delta:={}`

`delta := delta + g1`

`Computing Lie bracket: [f,g1]`

`delta := delta + [f,g1]`

`Computing Lie bracket: [delta2,delta1]`

`Computed Lie bracket is linear dependent.`

`Computing Lie bracket: [f,delta2]`

`delta := delta + [f,delta2]`

`Computing Lie bracket: [delta3,delta1]`

`Computed Lie bracket is linear dependent.`

`Computing Lie bracket: [delta3,delta2]`

`Computed Lie bracket is linear dependent.`

`Computing Lie bracket: [f,delta3]`

`Computed Lie bracket is linear dependent.`

`Dimension of reachability distribution:`, 3

This result coincides with the analysis of the linearized system.

Next the observability is checked, readjusting the simplifier routines.

>	`AIsys/SIMP`[2]:=0;

>	print(`For negative valve displacement`); ODistr(f,g_neg,xp,xx):

`Initialization: delta := {}`

`Computing differential dh`

`delta := delta + dh`

`Computing Lie derivative: [delta1,f]`

`delta := delta + Lie derivative`

`Computing Lie derivative: [delta1,g1]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta2,f]`

`delta := delta + Lie derivative`

`Computing Lie derivative: [delta2,g1]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta3,f]`

`delta := delta + Lie derivative`

`Dimension of observability codistribution:`, 4

What happens if Q is used as input?

>	print(`With flow as input`); ODistr(f,g_Q,xp,xx):

`Initialization: delta := {}`

`Computing differential dh`

`delta := delta + dh`

`Computing Lie derivative: [delta1,f]`

`delta := delta + Lie derivative`

`Computing Lie derivative: [delta1,g1]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta2,f]`

`delta := delta + Lie derivative`

`Computing Lie derivative: [delta2,g1]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta3,f]`

`delta := delta + Lie derivative`

`Dimension of observability codistribution:`, 4

What happens without any leakages?

>	`AIsys/SIMP`[2]:={Cleak_ext1=0,Cleak_ext2=0,Cleak_int=0}:

>	print(`With flow as input and no leakages`); ODistr(f,g_Q,xp,xx):

`Initialization: delta := {}`

`Computing differential dh`

`delta := delta + dh`

`Computing Lie derivative: [delta1,f]`

`delta := delta + Lie derivative`

`Computing Lie derivative: [delta1,g1]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta2,f]`

`delta := delta + Lie derivative`

`Computing Lie derivative: [delta2,g1]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta3,f]`

`Computed Lie derivative is linear dependent.`

`Computing Lie derivative: [delta3,g1]`

`Computed Lie derivative is linear dependent.`

`Dimension of observability codistribution:`, 3

Next some other properties are checked

>	`AIsys/SIMP`[2]:=0:

Computing relative degree with MIMO and SISO algorithms

>	MRelDeg(f,[g_pos],[xp],xx): print(`SISO Algorithm`); RelDeg(f,g_pos,xp,xx,true):

`Computing Lie derivative: L[g1]h1`

`Computed Lie derivative is 0`

`Computing Lie derivative: lie := L[f]^1h1`

`Computing Lie derivative: lie := L[g1]lie `

`Computed Lie derivative is 0`

`Computing Lie derivative: lie:=L[f]^2h1`

`Computing Lie derivative: lie := L[g1]lie `

`Lie derivative <> 0 -> r1 :=3`

`Checking the rank of the decoupling matrix`

`Decoupling matrix seems to be regular.`

`The relative degree vector is `, [3]

`Computing Lie derivative : L[g]L[f]h`

`Computing Lie derivative of output h along f: L[f]^2h`

`Checking column dimension of jacobian A = d(t_vec)/dx`

`Computed relative degree is 3`, ` jacobian A is regular`

Computing the maximum relative degree

>	MaxRelDeg(f,[g_pos],xx);

`Initialization: delta := {}`

`delta := delta + g1`

`Computing involutive closure ...`

`Rank of involutive closure of delta is 1`

`Compute Lie bracket: [f,g1]`

`delta := delta + [f,g1]`

`Computing involutive closure ...`

`Rank of involutive closure of delta is 3`

`Computing Lie bracket [f,delta2]`

`delta := delta + [f,delta2]`

`Computing involutive closure ...`

`Rank of involutive closure of delta is 4`

So the maximum relative degree is 3

For control purposes the nonlinear feedback for the Input/Output Linearization can be computed by UTrans

If xs > 0
print(`For positive valve displacement use transformation:`);
nlp:=UTrans(f,g_pos,xp,xx,v,true);
If xs < 0
print(`For negative valve displacement use transformation:`);
nln:=UTrans(f,g_neg,xp,xx,v,true);

`Computing Lie derivative : L[g]L[f]h`

`Computing Lie derivative of output h along f: L[f]^2h`

`Checking column dimension of jacobian A = dt/dx`

`Computed relative degree is 3 and jacobian of A is regular`

`Used coefficients: a0, a1, a2, `

`Coefficient: a3 set to 1`

`Computing Lie derivative : L[g]L[f]h`

`Computing Lie derivative of output h along f: L[f]^2h`

`Checking column dimension of jacobian A = dt/dx`

`Computed relative degree is 3 and jacobian of A is regular`

`Used coefficients: a0, a1, a2, `

`Coefficient: a3 set to 1`

The problem in here is that the information of the whole state is necessary to implement this feedback law. Usually, only the piston position xp and the pressures P1 and P2 are measured. The velocity signal vp can only be obtained by approximate differentiation of xp. But this procedure often fails, since in many industrial environments the transducer noise is considerable high. Therefore, we pose the question if it is possible to find an exact linearizing output such that the static feedback law is independent of the velocity vp? This question can be answered by the algorithm ElimSt , see, e.g., [Schlacher2001] K. Schlacher, A. Kugi, R. Novak, Input to Output Linearization with Constrained Measurements, In Preprints of the 5th IFAC Symposium on "Nonlinear Control Systems", Saint Petersburg, Russia, July 4-6, 2001, pp.471-476, 2001.

To get less information we set

>	infolevel[AIsys]:=0:

The measurements are

>	meas:=[P1,P2,xp]:

The computation of the restrictions gives

>	pdes:=ElimSt(f,g_pos,xx,meas,h);

>	pde1:=pdes[1]: pde2:=pdes[2]:

The solution of the first restriction is clear

>	restr1:=pdsolve(pde1,h(op(xx)));

Substituting this result in the second one gives the desired output

>	restr2:=simplify(subs(pde1=0,h(op(xx))=h(op(meas)),pde2));

>	output:=pdsolve(restr2,h(op(meas)));

For instance the output

>	z:=P1*A1-P2*A2+Eoil*ln((V01+A1*xp)^A1/(V02-A2*xp)^A2);

gives the desired result.

>	simplify(subs(h(op(xx))=z,pde1)); simplify(subs(h(op(xx))=z,pde2));

The nonlinear feedback for positive valve displacement is now independent of the velocity and becomes

>	UTrans(f,g_pos,z,xx,v);

Sometimes only the load pressure P1 - P2 can be measured, what happens in this case?

>	meas:=[P1-P2,xp]:

>	pdes:=ElimSt(f,g_pos,xx,meas,h);

"Problem has no solution."

Thus, in the general case the problem has no solution. But what happens without leakages and with the flow Q as input?

>	`AIsys/SIMP`[2]:={Cleak_ext1=0,Cleak_ext2=0,Cleak_int=0}:

>	pdes:=ElimSt(f,g_Q,xx,meas,h);

>	pde1:=pdes[1]: pde2:=pdes[2]: pde3:=pdes[3]:

The influence of the first and the second restriction is clear

>	pdsolve(pde1,h(op(xx)));

>	pdsolve(pde2,h(op(xx)));

Thus, we get the same output as before

>	pdsolve(subs(h(op(xx))=h(P1,P2,xp),pde3),h(P1,P2,xp));

>

Disclaimer:  While every effort has been made to validate the solutions in this worksheet, Waterloo Maple Inc. and the contributors are not responsible for any errors contained and are not liable for any damages resulting from the use of this material.