Application Center - Maplesoft

App Preview:

Squirrel cage induction motor

You can switch back to the summary page by clicking here.

Learn about Maple
Download Application


 

AIsys_3.mws

Full and Transformed Model of a Squirrel Cage Induction Motor - 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 Full and Transformed Model of a Squirrel Cage Induction Motor.

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;

libname :=

> with(AIsys);

[ElimSt, GLie, Interior_prod, InvClos, Invar, Invol...
[ElimSt, GLie, Interior_prod, InvClos, Invar, Invol...

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

> ?AIsys

Example 3: Full and Transformed Model of a Squirrel Cage Induction Motor

[Maple OLE 2.0 Object]

The above figure presents the schematic diagram for the windings of an induction motor. The mathematical model of the induction motor takes the following form:

> diff(psi(t),t) = -R*inv(L(theta))*psi+u;
diff(omega(t),t) = 1/Theta*(-1/2*psi^T*diff(''inv(L(theta))'',theta)*psi-M_L);
diff(theta(t),t) = omega;

diff(psi(t),t) = -R*inv(L(theta))*psi+u

diff(omega(t),t) = 1/Theta*(-1/2*psi^T*diff('inv(L(...

diff(theta(t),t) = omega

Here, psi denotes the flux linkage vector and theta is the angle of the rotor measured in the stator frame of coordinates, omega is the angular velocity of the rotor, M_L is the load torque, u is the vector of stator voltages and Theta is the moment of inertia, R is the constant resistance matrix with the stator winding resistance R_1 and the rotor winding resistance R_2 and L(theta) denotes the inductance matrix with the main inductance L_M and the leakage inductances L_S1 and L_S2 of stator and rotor, respectively.

> with(linalg):

Warning, the protected names norm and trace have been redefined and unprotected

> R11 := evalm(R_1*array(identity,1..3,1..3)):
R22 := evalm(R_2*array(identity,1..3,1..3)):
R := stackmatrix(augment(R11,matrix(3,3,[seq(0,i=1..9)])),
augment(matrix(3,3,[seq(0,i=1..9)]),R22));
L11 := matrix(3,3,[L_M+L_S1,-L_M/2,-L_M/2,-L_M/2,L_M+L_S1,-L_M/2,
-L_M/2,-L_M/2,L_M+L_S1]):
L22 := matrix(3,3,[L_M+L_S2,-L_M/2,-L_M/2,-L_M/2,L_M+L_S2,-L_M/2,
-L_M/2,-L_M/2,L_M+L_S2]):
L12 := matrix(3,3,[cos(theta),cos(theta+2*Pi/3),cos(theta-2*Pi/3),
cos(theta-2*Pi/3),cos(theta),cos(theta+2*Pi/3),
cos(theta+2*Pi/3),cos(theta-2*Pi/3),cos(theta)]):
L21 := transpose(L12):
L := stackmatrix(augment(L11,L12),augment(L21,L22));

R := matrix([[R_1, 0, 0, 0, 0, 0], [0, R_1, 0, 0, 0...

L := matrix([[L_M+L_S1, -1/2*L_M, -1/2*L_M, cos(the...

> psi := [psi_S1,psi_S2,psi_S3,psi_R1,psi_R2,psi_R3];

psi := [psi_S1, psi_S2, psi_S3, psi_R1, psi_R2, psi...

> f1 := map(simplify,map(expand,evalm(-R&*inverse(L)&*psi),trig),trig):

> f1 := convert(f1,list):

> M_el := -1/2*simplify(evalm(psi&*map(simplify,map(expand,map(diff,
inverse(L),theta),trig),trig)&*psi)):

> f2 := 1/Theta*(M_el-M_L):

> f := [op(f1),f2,omega];

f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...
f := [-R_1*(3*psi_S1*L_M^2+2*psi_S1*L_M*L_S2+6*psi_...

> g := [[1,0,-1,0,0,0,0,0],[0,1,-1,0,0,0,0,0]];

g := [[1, 0, -1, 0, 0, 0, 0, 0], [0, 1, -1, 0, 0, 0...

Define simplifiers

> `AIsys/SIMP`[1] := {sin(theta)^2+cos(theta)^2 = 1};
`AIsys/IntProd/simp` := proc(x)
simplify(x);
end;

`AIsys/SIMP`[1] := {sin(theta)^2+cos(theta)^2 = 1}

`AIsys/IntProd/simp` := proc (x) simplify(x) end pr...

Setup mathematical model for AIsys package

> xx := [op(psi),omega,theta];

xx := [psi_S1, psi_S2, psi_S3, psi_R1, psi_R2, psi_...

> f := M_AIsys(f);
g := M_AIsys(g);

f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...
f := TABLE([compts = vector([-R_1*(3*psi_S1*L_M^2+2...

g := [TABLE([compts = vector([1, 0, -1, 0, 0, 0, 0,...

Check, if system is reachable

> RDistr(f,g,xx):

`Dimension of reachability distribution:`, 6

As one can see, there is a subsystem of order 2 which is not reachable. This is why, we perform the well known Blondel-Park transformation in order to obtain a reduced model of the induction motor, where the unreachable part is eliminated. For the subsequent investigations, we take the mathematical model of the induction motor benchmark of R. Ortega (see www.supelec.fr/invi/lss/fr/personnels/ortega/benchmi/benchmi.html for further details). Here, x1 is the angular velocity, x2 and x3 are the components of the transformed rotor flux and x4, x5 are the components of the transformed stator current, tau_L denotes the load torque and u[1] , u[2] are the components of the transformed stator voltage.

> restart;
diff(x(t),t) = f(x)+g[1](x)*u[1]+g[2](x)*u[2];

diff(x(t),t) = f(x)+g[1](x)*u[1]+g[2](x)*u[2]

> f := [n_p*L_sr/D_m/L_r*(x2*x5-x3*x4)-R_m/D_m*x1-tau_L/D_m,
-1/T_r*x2-n_p*x1*x3+L_sr/T_r*x4,
n_p*x1*x2-1/T_r*x3+L_sr/T_r*x5,
L_sr/T_r/sigma/L_s/L_r*x2+n_p*L_sr/sigma/L_s/L_r*x1*x3-gamma*x4,
-n_p*L_sr/sigma/L_s/L_r*x1*x2+L_sr/T_r/sigma/L_s/L_r*x3-gamma*x5];
g := [[0,0,0,1/sigma/L_s,0],[0,0,0,0,1/sigma/L_s]];
xx := [x1,x2,x3,x4,x5];

f := [n_p*L_sr/D_m/L_r*(x2*x5-x3*x4)-R_m/D_m*x1-tau...
f := [n_p*L_sr/D_m/L_r*(x2*x5-x3*x4)-R_m/D_m*x1-tau...

g := [[0, 0, 0, 1/(sigma*L_s), 0], [0, 0, 0, 0, 1/(...

xx := [x1, x2, x3, x4, x5]

> libname := "C:/mylib/aisys", libname:
with(AIsys);
`AIsys/Global/simp` := proc(x)
simplify(x);
end:

[ElimSt, GLie, Interior_prod, InvClos, Invar, Invol...
[ElimSt, GLie, Interior_prod, InvClos, Invar, Invol...

Setup mathematical model for AIsys package

> f := M_AIsys(f);
g := M_AIsys(g);

f := TABLE([compts = vector([n_p*L_sr/D_m/L_r*(x2*x...
f := TABLE([compts = vector([n_p*L_sr/D_m/L_r*(x2*x...
f := TABLE([compts = vector([n_p*L_sr/D_m/L_r*(x2*x...
f := TABLE([compts = vector([n_p*L_sr/D_m/L_r*(x2*x...

g := [TABLE([compts = vector([0, 0, 0, 1/(sigma*L_s...

Check, if system is reachable

> RDistr(f,g,xx):

`Dimension of reachability distribution:`, 5

Obviously, this is the fact.

Check, if the system is exact input-state linearizable

> MLinPart(f,g,xx):

`No Exact Linearization via feedback possible.`

Since this is not possible, we calculate the maximum possible relative degree

> MaxRelDeg(f,g,xx);

[2, 2]

The maximum possible relative degree is r = [2,2]. This can be obtained by the following choice H for the output.

> H := [x1,x2^2+x3^2];

H := [x1, x2^2+x3^2]

> MRelDeg(f,g,H,xx);

`Decoupling matrix seems to be regular.`

TABLE([degvec = [2, 2], G_matrix = _rtable[2984028]...

Check, if the system is observable, if x1 (angular velocity) and one transformed current (x4) is measurable.

> H := [x1,x4];

H := [x1, x4]

> ODistr(f,g,H,xx):

`Dimension of observability codistribution:`, 5

The dimension of observability codistribution is 5 and hence the system is observable.

Let us assume, the output H is the induction motor torque and the rotor flux norm.

> H := [n_p*L_sr/L_r*(x2*x5-x3*x4),x2^2+x3^2];

H := [n_p*L_sr/L_r*(x2*x5-x3*x4), x2^2+x3^2]

> MRelDeg(f,g,H,xx);

`Decoupling matrix seems to be regular.`

TABLE([degvec = [1, 2], G_matrix = _rtable[17123524...

The corresponding control law reads as:

> MUTrans(f,g,H,xx,new_inp);

`Checking the rank of the decoupling matrix`

`Decoupling matrix seems to be regular.`

`Computing inverse matrix`

`Computing Transformation`

TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...
TABLE([degvec = [1, 2], utrans = vector([-x3/(x2^2+...

>

>

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.