Control Engineering I.pdf

UNIVERSIDAD CARLOS III DE MADRID Departamento de Ingeniería de Sistemas y Automática

NOTES ON CONTROL SYSTEMS I

L. Moreno, S. Garrido, C. Balaguer and M. A. Salichs

ii

CONTENTS 1. Introduction . . . . . . . . . . . . . . . . . . . 1.1 Introduction . . . . . . . . . . . . . . . . . . 1.2 Transform concept . . . . . . . . . . . . . . 1.3 Fourier Transform . . . . . . . . . . . . . . 1.3.1 Inverse Fourier Transform . . . . . . 1.4 Laplace Transform . . . . . . . . . . . . . . 1.4.1 Properties of the Laplace transform . 1.4.2 Laplace transforms of some functions 1.4.3 Inverse Laplace transform . . . . . . 1.4.4 Resolution of differential equations . 1.4.5 Derivative operator . . . . . . . . . .

. . . . . . . . . . .

. . . . . . . . . . .

. . . . . . . . . . .

. . . . . . . . . . .

. . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1 1 1 1 3 3 4 6 7 10 11

2. Classical techniques for system modeling 2.1 Introduction . . . . . . . . . . . . . . . . . 2.2 Input/Output models . . . . . . . . . . . . 2.3 Time models . . . . . . . . . . . . . . . . 2.3.1 Basic concepts . . . . . . . . . . . 2.3.2 Models in continuous time . . . . . 2.3.3 Rotational motion . . . . . . . . . 2.3.4 Electrical systems . . . . . . . . . 2.4 Linearization . . . . . . . . . . . . . . . . 2.5 Block diagrams . . . . . . . . . . . . . . . 2.5.1 Diagrams reduction . . . . . . . . 2.5.2 Transposition of sum points . . . .

. . . . . . . . . . . .

. . . . . . . . . . . .

. . . . . . . . . . . .

. . . . . . . . . . . .

. . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13 13 14 14 14 17 20 21 24 26 26 27

3. Time Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1 Introduction. Time Analysis of Systems . . . . . . . . . . . . . . 3.2 Temporal Response . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.1 Impulse response . . . . . . . . . . . . . . . . . . . . . . 3.2.2 Response of First Order systems . . . . . . . . . . . . . . 3.2.3 Response of the second order systems . . . . . . . . . . . 3.3 Parameters characterizing the response of a second order system or specifications . . . . . . . . . . . . . . . . . . . . . . . . . . . .

45 45 45 46 47 50

. . . . . . . . . . . .

59

4. Stability of dynamical systems . . . . . . . . . . . . . . . . . . . 4.1 Classical stability analysis methods . . . . . . . . . . . . . . . . . 4.1.1 Routh’s method . . . . . . . . . . . . . . . . . . . . . .

73 74 74

5. Errors in steady state . . . . . . . . . . . . . . . . . . . . . . . . 5.1 Introduction to errors . . . . . . . . . . . . . . . . . . . . . . . . 5.1.1 Gain position (static) . . . . . . . . . . . . . . . . . . . . 5.1.2 Gain speed . . . . . . . . . . . . . . . . . . . . . . . . . 5.1.3 Gain acceleration . . . . . . . . . . . . . . . . . . . . . . 5.1.4 Type system . . . . . . . . . . . . . . . . . . . . . . . . 5.1.5 Relationship between type and gain of a system . . . . . . 5.1.6 Relationship between type and gain a system . . . . . . . 5.1.7 Steady state error of a unit feedback system . . . . . . . . 5.1.8 Errors to standard input . . . . . . . . . . . . . . . . . . . 5.1.9 The relationship between errors and system type . . . . . 5.1.10 Error input for polynomial in t . . . . . . . . . . . . . . .

81 81 82 83 83 84 84 85 86 86 88 88

6. Root Locus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.1 Concept of Root Locus . . . . . . . . . . . . . . . . . . . . . . . 6.2 Magnitude and Angle Condition (or Module and Argument Criteria) 6.3 Tracing Rules of the Root Locus. . . . . . . . . . . . . . . . . . .

91 91 92 93

7. Classical Control Techniques . . . . . . . . . . . . . . . . . . . 7.1 Problem approach . . . . . . . . . . . . . . . . . . . . . . . . . . 7.2 Design of classical controllers . . . . . . . . . . . . . . . . . . . 7.3 Specifications . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.4 PID controllers . . . . . . . . . . . . . . . . . . . . . . . . . . . 7.4.1 Basic structure of a PID controller . . . . . . . . . . . . . 7.4.2 Ziegler-Nichols methods . . . . . . . . . . . . . . . . . . 7.5 Analytical methods for PID controllers design . . . . . . . . . . . 7.5.1 Pole placement method for adjusting PID . . . . . . . . . 7.5.2 Method based on the dominant poles for adjusting PID . .

99 99 100 101 102 103 104 107 107 110

8. FreQuency models. Bode diagrams . . . . . . . . . . . . . . . . 8.1 Asymptotic Bode diagram . . . . . . . . . . . . . . . . . . . . . 8.1.1 Asymptotic Bode diagram. Constants K . . . . . . . . . . 8.1.2 Asymptotic Bode diagram. Poles / zeros in the origin s1n . 8.1.3 Asymptotic Bode plot. Negative real Polo/ Zeros (1+T1 s)n 8.1.4 Asymptotic Bode plot. Negative complex Pole/Zeros . . . 8.1.5 Asymptotic Bode plot. Positive real Poles/Zeros. . . . . . 8.1.6 Asymptotic Bode plot. Positive complex Poles/Zeros. . . .

125 127 128 128 129 130 132 134

9. Nyquist Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135 vi

10. Design of regulators by frecuency methods . . . . . . . . . . . 10.1 Relationship between parameters of relative stability and transient response . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10.2 Lead Networks . . . . . . . . . . . . . . . . . . . . . . . . . . . 10.3 Lag networks . . . . . . . . . . . . . . . . . . . . . . . . . . . .

149 151 151 155

Some historical data . . . . . . . . . . . . . . . . . . . . . . . . . . . 157

vii

viii

LIST OF TABLES 7.1 7.2

Parameter values given by Ziegler-Nichols method for the step response . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105 Parameter values given by Ziegler-Nichols frequency method . . 107

x

LIST OF FIGURES 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 2.9 2.10 2.11 2.12 2.13

Impulse response and a sum of delayed pulses . . . . . . . . . . Mass-spring-damper system . . . . . . . . . . . . . . . . . . . . Mechanical translation system with two inputs and two outputs . . Example of mechanical system: 1/4 Car Model . . . . . . . . . . Elements of mechanical systems of translation and rotation . . . . Example of mechanical rotation system . . . . . . . . . . . . . . RLC circuit with voltage as input . . . . . . . . . . . . . . . . . System example . . . . . . . . . . . . . . . . . . . . . . . . . . Example of RLC network with two meshes . . . . . . . . . . . . Linearization of a function on two points (x0 , y0 ) y (x1 , y1 ) . . . Diagrams reduction . . . . . . . . . . . . . . . . . . . . . . . . . Transposition of sum points . . . . . . . . . . . . . . . . . . . . a) Starting system of the example, and b) first simplification: feedback loop . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.14 a) Transposition of the bifurcation simplifying and streamlining the feedback blocks in parallel. b) Final transfer function. . . . . . . . 2.15 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.16 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

16 19 19 20 21 21 23 23 24 24 27 28

3.1 3.2 3.3 3.4 3.5 3.6

46 47 48 48 49

Response of a system . . . . . . . . . . . . . . . . . . . . . . . . Impulse response of a first order system . . . . . . . . . . . . . . Impulse response of a first order system . . . . . . . . . . . . . . Unit step response of a first order system . . . . . . . . . . . . . . Unit step response of a first order system . . . . . . . . . . . . . . Determination of the transfer function from the unit step response of a first order system . . . . . . . . . . . . . . . . . . . . . . . . 3.7 Response to a ramp unit of a first order system . . . . . . . . . . . 3.8 Unit ramp response of a first order system . . . . . . . . . . . . . 3.9 Relation among the parameters of a second order system. . . . . 3.10 Classification of second-order systems . . . . . . . . . . . . . . . 3.11 Impulse response of the second order underdamped systems. . . . 3.12 The real exponential function eσt is increased when the real coefficient σ is positive and decreasing when it is negative. . . . . . .

28 28 32 34

49 50 50 51 52 53 55

3.13 The function eαt sin(βt + φ) is a growing sinusoid if α is positive and decreasing if negative. . . . . . . . . . . . . . . . . . . . . . 3.14 The unit step responses of a second order system according to the coefficient of friction: 1) ζ = 0, oscillator (poles on the imaginary axis). 2) 0 < ζ < 1, underdamped system (complex conjugate poles in the left half. 3) ζ = 1, critically damped system (real and equal poles). 4) ζ > 1, overdamped system (real and different poles). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.15 Parameters of the temporal response . . . . . . . . . . . . . . . .

56 60

4.1 4.2

Feedback loop . . . . . . . . . . . . . . . . . . . . . . . . . . . . System of example 5.1.4 . . . . . . . . . . . . . . . . . . . . . .

73 77

5.1 5.2 5.3 5.4 5.5 5.6 5.7 5.8 5.9

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Relationship between type and gain a system . . Relationship between type and gain a system . . Unit feedback system. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . The relationship between errors and system type.

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

. . . . . . . . .

81 82 84 85 86 87 87 88 88

6.1 6.2 6.3 6.4 6.5 6.6 6.7 6.8

Meaning of the magnitudes and angles. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Tracing rules for Root Locus. . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

93 94 94 95 95 95 96 98

7.1 7.2 7.3 7.4 7.5

Block diagram of a general control scheme . . . . . . . . . . . . Controller in the principal loop . . . . . . . . . . . . . . . . . . . Parameters of the time response . . . . . . . . . . . . . . . . . . a −sL Response to step input of G(s) = sL e . . . . . . . . . . . . . Obtaining the critical gain and critical period for the second method of Ziegler-Nichols . . . . . . . . . . . . . . . . . . . . . . . . . . Response system of Example 7.5.1: (a) open loop and (b) closed loop with the designed controller . . . . . . . . . . . . . . . . . System of the example 7.5.2 . . . . . . . . . . . . . . . . . . . . Response system of Example 7.5.2: (a) open loop and (b) closed loop with the controller designed . . . . . . . . . . . . . . . . . .

99 100 102 105

7.6 7.7 7.8

xii

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

. . . . . . . .

55

106 110 113 114

7.9 7.10 7.11 7.12 7.13 7.14 7.15 7.16 7.17 7.18 7.19 8.1 8.2 8.3 8.4 8.5 8.6 8.7 8.8 8.9

Effect of the addition of a pole on the Root Locus . . . . . . . . . 115 Effect of the addition of a zero on the Root Locus . . . . . . . . . 116 Root locus for the example system 7.5.3 . . . . . . . . . . . . . . 116 Application of the magnitude and angle conditions for Example 7.5.3118 Response system of Example 7.5.3: (a) open loop and (b) closed loop with the controller designed . . . . . . . . . . . . . . . . . . 119 Example 7.5.5 System . . . . . . . . . . . . . . . . . . . . . . . 119 Distribution of poles and zeros of example 7.5.4 . . . . . . . . . . 120 Example 7.5.5 System . . . . . . . . . . . . . . . . . . . . . . . 121 The Root Locus doesn’t contains the desired poles positions. . . . 121 It is necessary to add a zero to displace the Root Locus to the left. 122 Use of the angle condition. . . . . . . . . . . . . . . . . . . . . . 122 Bode plot of a second order system . . . . . . . . . . . . . . . . . Asymptotic Bode Plot of Gain Term . . . . . . . . . . . . . . . . Asymptotic Bode plot. Poles in origin s1n . . . . . . . . . . . . . Asymptotic Bode plot. Negative real Poles/Zeros (1+T1 s)n . . . . Asymptotic Bode plot. Negative complex Poles/Zeros. . . . . . . Bode plot. Negative complex Poles/Zeros . . . . . . . . . . . . . Asymptotic Bode diagram of 1/s − 1. (Real positive poles). . . . Asymptotic Bode diagram ofs − 1. (Real positive zero). . . . . . Bode plot with complex poles with: a) negative real part, b) positive real part. . . . . . . . . . . . . . . . . . . . . . . . . . . . .

127 128 129 130 131 132 133 133

Closed path ρ and its image F (s) . . . . . . . . . . . . . . . . . . Sistema en bucle cerrado . . . . . . . . . . . . . . . . . . . . . . Metodo de Nyquist . . . . . . . . . . . . . . . . . . . . . . . . . Modified Nyquist method . . . . . . . . . . . . . . . . . . . . . . System of the example 5.1.6 . . . . . . . . . . . . . . . . . . . . Nyquist path for the example 5.1.8 . . . . . . . . . . . . . . . . . Nyquist diagram of example 5.1.6 . . . . . . . . . . . . . . . . . Nyquist path for a system with a pole at the origin . . . . . . . . . Nyquist path for a system with poles on the imaginary axis . . . . Stable system . . . . . . . . . . . . . . . . . . . . . . . . . . . . Unstable system . . . . . . . . . . . . . . . . . . . . . . . . . . . The system a) is more stable than the b) and c systems) . . . . . . Calculation of gain margin . . . . . . . . . . . . . . . . . . . . . Nyquist plot, Root locus and closed loop step response according 1 system gain increases G(s) = s(s+a)(s+b) con a, b > 0. . . . . . . 9.15 The two systems have the same profit margin, but a) is more stable than the b) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9.16 Phase Margin . . . . . . . . . . . . . . . . . . . . . . . . . . . .

135 136 137 137 138 138 139 139 140 141 141 142 143

9.1 9.2 9.3 9.4 9.5 9.6 9.7 9.8 9.9 9.10 9.11 9.12 9.13 9.14

xiii

134

144 145 146

9.17 Gain margin and phase margin in Bode diagram . . . . . . . . . . 146 10.1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10.2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10.3 RLC lead network and its corresponding diagram of poles and zeros. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10.4 Bode and Nyquist plot for a lead network according to different parameter values α. . . . . . . . . . . . . . . . . . . . . . . . . . 10.5 Bode diagram of the system and lead network. . . . . . . . . . . . 10.6 Comparison of Bode plots and unit step responses. . . . . . . . . 10.7 Adding a lead network to the system changes the magnitude and phase, so it is difficult to predict the new crossing point. . . . . . . 10.8 Designing a lead network. . . . . . . . . . . . . . . . . . . . . . 10.9 Lag network and pole-zero diagram. . . . . . . . . . . . . . . . . 10.10Frequency response of lag networks. . . . . . . . . . . . . . . . . 10.11Comparison of the response of the systems. . . . . . . . . . . . . 10.12Designing a lag network. . . . . . . . . . . . . . . . . . . . . . . 10.13The effect of the change of T. . . . . . . . . . . . . . . . . . . . .

xiv

149 150 152 152 152 153 153 154 155 155 155 156 156

1. INTRODUCTION 1.1

Introduction

One of the problems that have to face those working in the control or identification of processes lies in the existence of an variety of mathematical tools to study, analyze and design such systems. For example, if we analyze and control a continuous system, we can model the system in terms of its temporal response by using the transfer function. The system can be controlled in continuous time which will require the use of Laplace transform, in the case of continuous time. Sometimes, we are interested in the study of the frequency response of the system for this reason we will study the response in steady state (based on the transfer function). In the study and analysis of a system we will find different mathematical tools that allow us to characterize, analyze and design a control system so that we consider the different effects in the system. The basic objective of this chapter is to introduce the tools mathematics that are typically used to analyze and synthesize systems of deterministic control in the time domain, in continuous time as well as in frequency.

1.2

Transform concept

When analyzing and designing control systems for a wide variety of processes, we find the need to solve differential or integro-differential equations in order to know the response the process will give us to certain inputs or control actions. These differential equations can be solved directly or can be converted to solve other equations of easier resolution. There are a number of operational methods that allow the resolution of problems using transform equations in an easier way. Among the most frequently used transformations for solution of ordinary differential equations of linear transformations are an integral, and among these a widely used in control is called Laplace transform.

1.3

Fourier Transform Given a real function f (t), the Fourier transform of the function f (t) is defined

as:

2

1. Introduction

Z

+∞

F (ω) = F[f (t)] =

f (t)e−jωt dt

(1.3.1)

−∞

where ω a variable frequency, −∞ < ω < +∞, expressed in (rad/s). The function f (t) is said to have Fourier transform if the integral 1.3.1 exists for all values of ω. Sufficient conditions for the integral existence is that the function f (t) has only a finite number of discontinuities, maxima and minima over a finite time interval, and that function f (t) is absolutely integrable, ie Z

+∞

|f (t)|dt < ∞

(1.3.2)

−∞

The complex exponential term of the transformed 1.3.1 can be expressed by applying Euler’s formula, as e−jωt = cos(ωt) − j sin(ωt)

(1.3.3)

The complex exponential term of the transformed 1.3.1 can be expressed by applying Euler’s formula, as Z

+∞

f (t)[cos(ωt) − j sin(ωt)]dt

F (ω) = −∞ Z +∞

Z

+∞

f (t) cos(ωt)dt − j

= −∞

f (t) sin(ωt)dt

(1.3.4)

−∞

In this expression 1.3.4 of the Fourier transform you can see the existence of two terms, one real and one imaginary, which we call R(ω) and I(ω), respectively, so F (ω) is given by F (ω) = R(ω) + jI(ω)

(1.3.5)

In the equation 1.3.5 can be seen that the Fourier transform F (ω) is a complex function of the variable ω. For every value of the variable ω, will have a value of F (ω) which has a module, |F (ω)|, and an argument, ∠F (ω). p R2 (ω) + I 2 (ω)   I(ω) −1 ∠F (ω) = tan R(ω) |F (ω)| =

(1.3.6)

If we represent the values on a plot |F (ω)| for every ω would obtain what is called the spectrum signal magnitudes f (t), and if we did the same for argument ∠F (ω) would obtain what is known as phase spectrum signal f (t).

1.4. Laplace Transform

1.3.1

3

Inverse Fourier Transform

The inverse transform is an operation that allows us to obtain the original function f (t) from the Fourier transform F (ω). This operation is given by the following expression f (t) = F −1 [F (ω)] =

1.4

+∞

Z

1 2π

F (ω)ejωt dω

(1.3.7)

−∞

Laplace Transform

Given a real function f (t) such that f (t) = 0 < t < ∞, the unilateral Laplace transform of the function f (t) is defined as: ∞

Z F (s) = L[f (t)] =

f (t)e−st dt

(1.4.1)

0

There is also the bilateral Laplace transform, in which integration limits range from −∞ to +∞, but it is not typically used in control engineering, and we are not going to study it. Unilateral Laplace transform, which in and for simplicity we will call hereinafter as Laplace transform, exists always if the integral 1.4.1 converges, and it is a necessary and sufficient condition that: Z

∞

|f (t)e−σt |dt < ∞

(1.4.2)

0

for σ finite and real. Given the Laplace transform of a function f (t), the set of all complex numbers s for which the 1.4.1 exists is called region of convergence. Given a function f (t), if the Laplace transform L(f (s)) is particularized for those complex numbers in which s = jω we obtain that Z F (s)|s=jω = F (jω) =

∞

f (t)e−jωt dt

(1.4.3)

0

if in addition the function f (t) = 0 for t < 0, the equation 1.4.3 is equal to the Fourier transform F (ω) of the function f (t) (compare with the expression that we saw in 1.3.1). Namely F (ω) = F (s)|s=jω (1.4.4) This tells us that we can obtain the Fourier transform, F (ω) of a function f (t) (where the value of this function for t < 0 is zero) from its transform Laplace, F (s), replacing in the expression of the Laplace transform the variable s for jω.

4

1. Introduction

1.4.1

Properties of the Laplace transform

The Laplace transform has a number of very useful properties which means it can be considered as an operator that converts a temporal function f (t) in the time domain in a function F (s) of complex variable s. These properties are: 1. Linearity A very important property of the Laplace transform is that it is a linear operator. Given two functions f (t) and g(t) and a constant k is verified that: L[f (t) + g(t)] = L[f (t)] + L[g(t)]

(1.4.5)

L[kf (t)] = kL[f (t)]

(1.4.6)

and

2. Differentiation The Laplace transform of the differential of a function f (t) can be expressed as follows:  d f (t) = sF (s) − lim f (t) = sF (s) − f (0+ ) L t→0 dt 

(1.4.7)

where F (s) the transform of the function f (t) and f (0+ ) is the limit of the function f (t) when t tends to zero from the right. In general, for higher derivatives we have the following expression:  dn dn−1 n n−1 + n−2 d + L f (t) = s F (s) − s f (0 ) − s f (0 ) + . . . + f (0+ ) dtn dt dtn−1 (1.4.8) 

3. Integration The Laplace transform of the integral of a function f (t) with respect to time is the Laplace transform of the function divided by s, ie: Z L

t

 f (τ )dτ =

0

F (s) s

(1.4.9)

and if generalized to higher orders integrations Z

t1

Z

L

t2

Z ...

0

0

tn

 f (τ )dτ dτ1 . . . dτn =

0

F (s) sn

(1.4.10)

1.4. Laplace Transform

5

4. Time shift A time shift of a units of the function f (t), in the complex domain is equivalent to multiply the Laplace transform of the function F (s) for e−s , that is: L[f (t − a)u0 (t − a)] = e−sa F (s) (1.4.11) where u0 (t − a) represents the unit step function shifted in time a units. 5. Initial value theorem If the Laplace transform of the function f (t) is F (s), then lim f (t) = lim sF (s)

t→0

s→∞

(1.4.12)

if the limit when time t → 0 exists. 6. Final value theorem If the Laplace transform of the function f (t) is F (s), and if sF (s) is analytic on the right side of the plane-s and on the imaginary axis, then lim f (t) = lim sF (s)

t→∞

s→0

(1.4.13)

This theorem is useful in the analysis and design of control systems because it tells us the final value of a function in the time domain knowing the value of its Laplace transform for s = 0. 7. Complex shift The Laplace transform of a function f (t) multiplied by e∓αt , where α is a constant, is equivalent to replace the Laplace transform of the function, F (s), s with s ± a, ie L[e∓αt f (t]] = F (s ± a)

(1.4.14)

8. Real convolution Let F1 (s) and F2 (s) the Laplace transforms of the functions f1 (t) and f2 (t), respectively, and f1 (t) = 0, and f2 (t) = 0, for t < 0 then

F1 (s)F2 (s) = L[f1 (t) ∗ f2 (t)] Z t = L[ f1 (τ )f2 (t − τ )dτ ] 0 Z t = L[ f2 (τ )f1 (t − τ )dτ ] 0

(1.4.15) (1.4.16) (1.4.17)

6

1. Introduction

in this expression the symbol ∗ denotes the convolution operation in the time domain.

1.4.2

Laplace transforms of some functions

Example 1.4.1: The Laplace transform of the unit step function u0 (t) defined as  f (t) = u0 (t) =

1, t ≥ 0 0, t < 0

(1.4.18)

is obtained as Z

∞

F (s) = L[u0 (t)] =

−st

u0 (t)e 0

1 −st ∞ 1 dt = − e = s s 0

(1.4.19)

Example 1.4.2: The Laplace transform of the exponential function f (t) = e−αt for t ≥ 0 is obtained as

F (s) = L[e−αt ] =

Z

∞

e−αt e−st dt = −

0

1 −(s+α)t ∞ 1 e = s+α s+α 0

(1.4.20)

In the table you can see the transforms of some of the more usual functions.

f (t)

L(s)

δ(t)

1

u0 (t)

1 s 1 s2 1 s+a 1 (s+a)2

t −at

e

te−at tn , t ≥ 0 tn e−at , t ≥ 0 sin(ωt) cos(ωt) e−at sin ωt e−at cos ωt

n! , n = 1, 2, . . . sn+1 n! , n = 1, 2, . . . (s+a)n+1 ω (s2 +ω 2 ) s (s2 +ω 2 ) ω (s+a)2 +ω 2 s+a (s+a)2 +ω 2

1.4. Laplace Transform

1.4.3

7

Inverse Laplace transform

The operation to obtain the function f (t) from the Laplace transform F (s) is called the inverse Laplace transform, and it is written as f (t) = L−1 [F (s)] The inverse of the Laplace transformation can be calculated as Z c−j∞ 1 F (s)est ds f (t) = 2πj c+j∞

(1.4.21)

(1.4.22)

where c is a real constant greater than the real part of all singularities of F (s). This integral is evaluated in the complex s-plane. From the standpoint of process control, to obtain the inverse transformation is extremely costly and impractical. So, tables of Laplace transforms are used on the partial fraction expansion. Partial fraction expansion In most control applications, it is not necessary to use the integral inversion given by 1.4.22 for the inverse Laplace transform. Usually both the differential equations and its solutions in control problems take the form of rational functions in s, F (s) =

N (s) D(s)

(1.4.23)

where N (s) and D(s) are polynomials in s. It is assumed that the degree of the polynomial N (s) in the numerator s is less or equal than the degree of the polynomial D(s) in the denominator. This expression can be written as follows F (s) =

sn

+ a1

sn−1

N (s) + . . . + an−1 s + an

(1.4.24)

where the coefficients ai are real coefficients. The roots of the polynomial N (s) are called zeros and the roots of the denominator D(s) are called poles of the function F (s). The basic idea of the method of partial fraction expansion is to decompose the polynomial F (s) into a sum of simple fractions such that these are transformed into simple and known functions such that for each fraction becomes very easy to calculate the inverse transform. As the Laplace transform is a linear operator, after calculating the inverse transforms of the simple fractions, to obtain the inverse transform of the function F (s) simply we need to add the inverse transforms of the simple fractions. We consider the different cases.

8

1. Introduction

Case of simple and real poles If all poles of D(s) are simple and real, then the equation 1.4.24 takes the form F (s) =

N (s) (s + s1 )(s + s2 ) . . . (s + sn )

(1.4.25)

and it ∀i 6= j, Si 6= sj . Decomposing the expression 1.4.25 we obtain the partial fraction F (s) =

A1 A2 An + + ... + (s + s1 ) (s + s2 ) (s + sn )

where the coefficients A1 , . . . , An are obtained as   N (s) Ai = (s + si ) D(s) s=−si

(1.4.26)

(1.4.27)

Note that each of these simple fractions are Laplace transforms of exponential functions with exponent −si and coefficient Ai . Case of multiple real poles If some of the poles of D(s) are multiples, then the equation 1.4.24 takes the form F (s) =

N (s) (s + s1 )(s + s2 ) . . . (s + sk−1 )(s + sk )n−k+1

(1.4.28)

and it has a −sk root with multiplicity n − k + 1. Decomposing the expression in fractions 1.4.28 obtain the partial expansion F (s) =

A1 A2 Ak−1 + + ... + + (s + s1 ) (s + s2 ) (s + sk−1 ) Akn−k Ak 1 Ak2 + + + ... + 2 (s + sk ) (s + sk ) (s + sk )n−k+1

(1.4.29)

where the coefficients A1 , . . . , Ak−1 corresponding to the simple poles are obtained as in the previous case, and for the multiple pole the coefficients are obtained as  Akn−k

=

(s + sk )

 D(s)

n−k+1 N (s)

(1.4.30)

s=−sk

Akn−k+1

= ...

Ak1

=

  d n−k+1 N (s) (1.4.31) (s + sk ) ds D(s) s=−sk   1 dn−k n−k+1 N (s) (s + s ) (1.4.32) k (n − k)! dsn−k D(s) s=−sk

1.4. Laplace Transform

9

Case of complex conjugate poles The partial fraction expansion seen in the case of poles simple and real, it is also valid for the case of poles complex conjugates of the form s = −σ +jw and s = −σ −jw. In this case the coefficients A−σ+jw and A−σ−jw are obtained as follows  N (s)  (s + σ − jw) D(s) s=−σ+jw  N (s)  = (s + σ + jw) D(s) s=−σ−jw

A−σ+jw =

(1.4.33)

A−σ−jw

(1.4.34)

Example 1.4.3: Suppose the following function F (s) =

2s + 3 s(s + 1)(s + 3).

(1.4.35)

If the partial fraction is expanded as follows A1 A2 A3 + + s (s + 1) (s + 3)

F (s) =

(1.4.36)

the coefficients will be determined by 2s + 3 =1 = (s + 1)(s + 3) s=0 s=0   2s + 3 1 = (s + 1)F (s) = =− s(s + 3) s=−1 2 s=−1   2s + 3 1 = (s + 3)F (s) = =− s(s + 1) s=−3 2 s=−3

A1 = A2 A3



 sF (s)

thus F (s) =

− 12 − 12 1 + + s (s + 1) (s + 3)

(1.4.37)

and its inverse transform will be the function 1 1 f (t) = 1 − e−t − e−3t 2 2

(1.4.38)

Example 1.4.4: Suppose the following function F (s) =

s2 + 3s (s + 1)3 .

If the partial fraction will be expanded as follows

(1.4.39)

10

1. Introduction

F (s) =

A1 A3 A2 + + 2 (s + 1) (s + 1) (s + 1)3

(1.4.40)

The coefficients are determined by A3 = A2 = A1 =



 (s + 1)3 F (s)

s=−1

= s2 + 3s

s=−1

= −2

 d (s + 1)3 F (s) = 2s + 3 =1 ds s=−1 s=−1  1 d2  2 3 =1 (s + 1) F (s) = 2! ds2 2! s=−1 s=−1

thus F (s) =

1 1 −2 + + 2 (s + 1) (s + 1) (s + 1)3

(1.4.41)

and its inverse transform will be the function f (t) = e−t + te−t − t2 e−t

1.4.4

(1.4.42)

Resolution of differential equations

One of the most useful properties of the Laplace transform is the simplicity with which you can solve linear differential equations. Here’s an example. Example 1.4.5: Suppose a system whose differential equation is as follows d2 x dx +2 +δ =0 2 dt dt

(1.4.43)

besides, we know that for this system the values of x and x˙ for t = 0, are, respectively, x(0) = 1 and x(0) ˙ = 0, and we want to obtain the output x(t) that give the system 1.4.43. Applying the property of the Laplace transform of the derivative will have the following [s2 X(s) − sx(0) − x(0)] ˙ + 2[sX(s) − x(0)] + 1 = 0 Now if we replace x(0) and x(0) ˙ by its value, we it will be [s2 X(s) − s] + 2[sX(s) − 1] + 1 = 0 s2 X(s) + 2sX(s) = s + 1 s+1 X(s) = s(s + 2)

1.4. Laplace Transform

11

Expanding this expression into partial fractions X(s) =

A1 A2 + s (s + 2)

(1.4.44)

The coefficients are determined by   1 (s + 1) = A1 = sX(s) = (s + 2) s=0 2 s=0   1 s + 1 = A2 = (s + 2)X(s) = s s=−2 2 s=−2 thus X(s) =

1 2

+

1 2

s (s + 2) and its inverse transform is the solution of the differential equation 1.4.43 for the initial conditions are assumed, i.e. 1 x(t) = (1 + e−2t ). 2

1.4.5

(1.4.45)

Derivative operator

It is frequently used the expression derivative operator, which is usually represented by the symbol p or D, although here the first one will be used, d y(t) dt Using this symbol, a differential equation of order n py(t) =

dn dn−1 d y(t)+a y(t) + . . . + an−1 y(t) + an y(t) = 1 n n−1 dt dt dt dn−1 d = b1 n−1 u(t) + . . . + bn−1 u(t) + bn u(t) dt dt can be expressed as a polynomial

(1.4.46)

(1.4.47)

pn y(t)+a1 pn−1 y(t)+. . .+an−1 py(t)+an y(t) = b1 pn−1 u(t)+. . .+bn−1 pu(t)+bn u(t) (1.4.48) thus [pn + a1 pn−1 + . . . + an−1 p + an ]y(t) = [b1 pn−1 + . . . + bn−1 p + bn ]u(t) A(p)y(t) = B(p)u(t) B(p) y(t) = u(t) A(p)

(1.4.49)

This operator has the advantage of working directly in the time domain, which is sometimes useful. This operator is frequently used in the classical methods of resolution of ordinary differential equations.

12

1. Introduction

Bibliografia [1] L. Ljung, System Identification: Theory for the User, Prentice Hall, 1987. [2] L. Ljung and T. Söderstrom, Theory and Practice of Recursive Identification, The MIT Press, 1983. [3] T. Söderstrom and P. Stoica, System Identification, Prentice Hall, 1989. [4] J. P. Norton, An Introduction to Identification, Academic Press, 1986. [5] G.R. Cooper and C.D McGuillem, Probabilistic Methods of Signal and System Analysis, Oxford University Press, 1999. [6] A. Papoulis, Probabilidad, variables aleatorias y procesos estocasticos, Ed. Universitaria de Barcelona, 1980. [7] B.C. Kuo, Automatic Control Systems, Ed. Prentice Hall, 1991. [8] K. Ogata, Discrete-Time Control Systems, Prentice Hall International Edition, 1995. [9] K. Ogata, Ingenieria de control moderna, Prentice Hall, 1998. [10] E. Andres Puente, Regulacion automatica, Servicio de publicaciones ETSIIMUPM, 1993.

2. CLASSICAL TECHNIQUES FOR SYSTEM MODELING 2.1

Introduction

A basic problem in control engineering and in general in many other sciences is to be able to predict what effect it will have some action on a physical system. Because what you want is to predict future response of the system, it is necessary to use some kind of model that allows us to make this prediction. There are several types of models, some of them such as scale models of what is sought is to reproduce the actual physical system by some model or more manageable system on which to experiment and extract results. These, though interesting, are not the models more interested in control, but those where predictions are made on a mathematical model of the system. Within the mathematical models that can be used to analyze the effects that different actions will have on the system can distinguish two groups of models: • Input/Output models. In this class of models seeks a mathematical description that expresses the relationship between the system input and the output thereof. These models do not describe the inner workings of the system, but merely the relation between input and output. We may find different systems but having the same ratio input/output so lead to the same mathematical model. These types of models are what we might call classic models. • State Space model based. Such models seek a deeper description of the system, since not only the relationship between the input and output features, but besides the behavior of a number of variables or internal variables of the system features. There are many systems in which the internal variables that can be used to characterize it can be chosen in different ways, so that the same system may be described by different mathematical models. These types of techniques have given rise to what is called modern control theory. This chapter will address the study of what are the techniques classical modeling systems

14

2.2

2. Classical techniques for system modeling

Input/Output models

These types of models are what we might call classic models , and within them can distinguish several types : temporal patterns , and stochastic frecuenciales . The temporal patterns of input / output feature allow the study of the temporal response of the system to a any input . Knowing the temporal response of the system allow us to analyze transient effects of short duration in the system, such as response peaks , rise and settling times . These models can be temporal models based on differential continuous-time or discrete-time models described by differential equations equations . The frecuenciales called models are based on the characterization of the relationship / out of a steady-state system to sinusoidal inputs . Finally the most common models are stochastic or deterministic models combine temporal and frequential noise or disturbance signals that add a stochastic behavior to the system and that needs to be modeled when possible .

2.3 2.3.1

Time models Basic concepts

Linear systems invariant in time The model input/output linear systems corresponding to the following general form, dn dn−1 y(t) + a y(t) + . . . + a0 y(t) = n−1 dtn dtn−1 dm dm−1 bm m u(t) + bm−1 n−1 u(t) + . . . + b0 u(t) dt dt

(2.3.1) (2.3.2)

where the coefficients ai , i = 0, 1, . . . , n − 1 and bj , j = 0, 1, . . . , m are real numbers, u(t) is the input to system and y(t) is the output. If these coefficients ai and bj do not depend on time, the system is said to be time invariant. These linear models have several advantages, including and perhaps the most important is the linearity in the response. That is, if α and β are two arbitrary constants and u1 (t) and u2 (t) are two inputs to the system, it is verified that if u(t) = αu1 (t) + βu2 (t)

(2.3.3)

then the output of the system is given by y(t) = αy1 (t) + βy2 (t)

(2.3.4)

2.3. Time models

15

Impulse response Since an invariant linear system in time, which for a input u(t) gives an answer y(t), is defined as bf impulse response of the system output g(t) that would give the entry system at the same is a unit impulse δ(t). This impulse response g(t) contains all the information needed on the system, and you can get the output of same, at any input u(t), no more to perform convolution in the time domain with the signal, ie y(t) = g(t) ∗ u(t) Z t g(τ )u(t − τ )dτ = 0 Z t = u(τ )g(t − τ )dτ

(2.3.5) (2.3.6) (2.3.7)

0

In general terms the impulse response g(t) is not very practical when modeling a system and that it is much easier to work in the field of Laplace, to solve integral expressions can be complex. Figure 2.1 you can see that impulse δ(t), if an input the output of the system is g(t), then, by linearity, to an input δ(t), the output shall be subject of ag(t). Since the system is time invariant, momentum shifted to a δ(t − t0 ) shall be subject of a displaced g(t − t0 ) response. To a composite input a superposition of scaled pulse corresponds to one output will be the sum of the responses. Transfer function It is defined as transfer function of a linear and time invariant system the Laplace transform of the impulse response of the system, assumed zero initial conditions. G(s) = L[g(t)]

(2.3.8)

The concept of transfer function is very useful when determining the output y of the system to a given input u, since this is given by the expression Y (s) = G(s)U (s)

(2.3.9)

from whose values ??in the complex field, developing in partial fractions and making the inverse transforms of each of these partial fractions was easy to obtain the analytical response that give the system in the time domain as seen in the previous chapter. This concept of transfer function is widely used in the future and it is convenient to make some clarifications on the matter:

16

2. Classical techniques for system modeling

y(t)

u(t) δ(t)

g(t) t

t u(t)

g(t)

y(t)

y(t)

u(t)

ag(t)

Escalado

aδ(t)

t

t

y(t)

u(t) Desplazamiento

δ(t-to)

t

t

to

u(t) 2δ(t) δ(t-to) to

g(t-to)

y(t) 1.5δ(t-t1)

t1

2g(t)+g(t-to)+1.5g(t-t1)

Superposición t

to

t1

t

Fig. 2.1: Impulse response and a sum of delayed pulses • The transfer function is defined for linear time-invariant system over time, so it is not defined in the case of nonlinear systems. • For the transfer function are assumed zero initial conditions. • The transfer function expresses the relationship between an input signal and an output signal, so one can express completely the dynamics of a system with one input and one output. If the system had more of an input or output, various transfer functions would be required to express all the possible relationships between each input and each output. In this case we would have multivariable systems: MIMO (Multiple Input - Multiple Output).

2.3. Time models

17

• The transfer function is independent of the input to the system. • The transfer function of a system is only rational function of s or z, and also with a denominator greater than or equal to the numerator degree, to be physically realizable.

2.3.2

Models in continuous time

A very important part of the systems we are interested in corresponds to control physical systems whose dynamics is described by the following type of ordinary differential equations dn dn−1 y(t) + a y(t) + . . . + a0 y(t) = n−1 dtn dtn−1 (2.3.10) dm−1 dm bm m u(t) + bm−1 n−1 u(t) + . . . + b0 u(t) dt dt where u(t) is the input to the system and y(t) the output. Assuming zero initial conditions and take Laplace transforms on both sides of the equation 2.3.10, the result is (sn +an−1 sn−1 +. . . +a0 )Y (s) = (bm sm +bm−1 sm−1 +. . . +b0 )U (s) (2.3.11) Solving Y (s) in Eq 2.3.11 we obtain that Y (s) =

sn + an−1 sn−1 + . . . + a0 U (s) bm sm + bm−1 sm−1 + . . . + b0

(2.3.12)

and here the general form of the transfer function G(s) for a physical system described by an ordinary differential equation of the form 2.3.10 is a rational function of the type s G(s) =

sn + an−1 sn−1 + . . . + a0 bm sm + bm−1 sm−1 + . . . + b0

(2.3.13)

Example 2.3.1: A car responds to an applied force f (t) from the engine power and a force due to friction ρv(t), proportional to the velocity v(t) of the car. The differential equation of the system is mv(t) ˙ = f (t) − ρv(t)

(2.3.14)

Taking Laplace transforms is obtained msV (s) = F (s) − ρV (s)

(2.3.15)

18

2. Classical techniques for system modeling

i.e. (ms + ρ)V (s) = F (s)

(2.3.16)

and the transfer function of the relative speed of the force is V (s) 1 = F (s) ms + ρ

(2.3.17)

Example 2.3.2: Consider the mass-spring-damper system of Fig. 2.2. Newton’s second law states that the sum of all forces acting on a rigid body is equal to its mass times its acceleration. X Fi = m · y¨ (2.3.18) where the minus sign represents the force Fm opposes the change of position. A spring is an element that stores potential energy and when stretched or compressed exert a force that is expressed by Hooke’s law: Fm = −K · y

(2.3.19)

where the minus sign represents the force Fm opposes the change of position. A damper is a device that converts energy into heat and represents a viscous friction. The force exerted on their ends can be expressed as Fa = −b · y˙

(2.3.20)

That is, the damper exerts a force that opposes the speed changing. The mass-spring-damper system shown in Fig. 2.2 can be modelled using these three expressions in the form of differential equation m · y¨ = u − K · y − b · y˙ where u it is an external force, which is the system input. Rearranging m · y¨ + b · y˙ + K · y = u

(2.3.21)

(2.3.22)

By taking Laplace transforms ms2 Y (s) + bsY (s) + KY (s) = U (s)

(2.3.23)

The transfer function is the output input ratio (in the domain s) 1 Y (s) = U (s) ms2 + bs + K

(2.3.24)

2.3. Time models

19

Fig. 2.2: Mass-spring-damper system Example 2.3.3: In the multivariate case of Figure 2.3, you must write two differential equations, one for each mass. It is useful to draw the corresponding diagram of forces (free body diagram). In this case the differential equations are m1 y¨1 = u1 + b1(y˙ 2 − y˙ 1 ) − k1 y1

(2.3.25)

m2 y¨2 = u2 − b1 (y˙ 2 − y˙ 1 ) − k2 y2

(2.3.26)

Fig. 2.3: Mechanical translation system with two inputs and two outputs

Example 2.3.4:

20

2. Classical techniques for system modeling

Fig. 2.4: Example of mechanical system: 1/4 Car Model

2.3.3

Rotational motion

In this case considered movements about a fixed axis. It is rather like translation movements, but the variables of interest are angular displacement, angular velocity and angular acceleration and instead of mass and forces, consider moments of inertia and torques. Newton’s second law for rotational motion says that the sum of all torques is equal to the moment of inertia by the angular acceleration. X ¨ Ti (t) = J · θ(t) (2.3.27) In addition to Newton’s second law can be considered the torsion springs, which Hooke’s law Tm (t) = −K · θ(t) (2.3.28) and the corresponding elements to the viscous dampers are friction rotation, whose expression is ˙ Ta (t) = −b · θ(t) (2.3.29) Example 2.3.5: Consider the mechanical system of rotation of the figure 2.6. If it starts from rest and a small torque is applied to the flywheel J1 , in the early stages of the movement, θ1 > θ2 > θ3 . This means that the torsion spring K is θ1 − θ2 and angle in the damper D2 is θ2 − θ3 . Therefore equations are J1 θ¨1 = T (t) − D1 θ˙1 − K(θ1 − θ2 ) J2 θ¨2 = K(θ1 − θ2 ) − D2 (θ˙2 − θ˙3 )

(2.3.30)

J3 θ¨3 = D2 (θ˙2 − θ˙3 ) − D3 θ˙3

(2.3.32)

(2.3.31)

2.3. Time models

21

Fig. 2.5: Elements of mechanical systems of translation and rotation

Fig. 2.6: Example of mechanical rotation system

2.3.4

Electrical systems

The three basic elements of electrical circuits are resistors, capacitors and inductors. Resistance is an element that dissipates energy by converting it into heat and its equation is given by Ohm’s law v(t) = R · i(t) A capacitor is an element that stores energy and its equation is Z 1 v(t) = i(t)dt C

(2.3.33)

(2.3.34)

Finally, the coil or inductance has the equation u(t) = L

di(t) dt

In addition, the basic laws are the laws of Kirchoff:

(2.3.35)

22

2. Classical techniques for system modeling

• At any node (junction) in an electrical circuit, the sum of currents flowing into that node is equal to the sum of currents flowing out of that node or equivalently The algebraic sum of currents in a network of conductors meeting at a point is zero. • The directed sum of the electrical potential differences (voltage) around any closed network is zero, or: More simply, the sum of the emfs in any closed loop is equivalent to the sum of the potential drops in that loop, or: The algebraic sum of the products of the resistances of the conductors and the currents in them in a closed loop is equal to the total emf available in that loop. In addition, to calculate the transfer function is preferable working in the frequency domain, ie, with the Laplace transforms of the equations: • For resistors, the expression v(t) = Ri(t) becomes V (s) = RI(s). You can also say that its complex impedance is ZR (s) = • For capacitors, the term v(t) =

1 C

R

V (s) I(s)

=R

1 Cs I(s). V (s) 1 I(s) = Cs

i(t)dt becomes V (s) =

You can also say that its complex impedance is ZC (s) =

• For the coils or inductors, the expression v(t) = L di(t) dt becomes V (s) = LSI(s). In this case the complex impedance is ZL (s) =

V (s) I(s)

= Ls.

Example 2.3.6: Consider the RLC network in Fig. 2.7. At the right appears the result when using complex impedances. The input loop, according to the 2nd Law of Kirchoff is 1 V (s) = LSI(s) + RI(s) + I(s) Cs and the output voltage is 1 VC (s) = I(s) Cs . The transfer function is the ratio 1 VC (s) Cs I(s) = V (s) LsI(s) + RI(s) +

1 Cs I(s)

=

1 Cs

Ls + R +

1 Cs

=

LCs2

1 + RCs + 1

2.3. Time models

23

Fig. 2.7: RLC circuit with voltage as input

Fig. 2.8: System example Example 2.3.7:

24

2. Classical techniques for system modeling

Fig. 2.9: Example of RLC network with two meshes

2.4

Linearization

Fig. 2.10: Linearization of a function on two points (x0 , y0 ) y (x1 , y1 ) The actual physical systems are nonlinear. In the above examples, the nonlinearities are small and can be considered negligible. But in many cases the system is clearly nonlinear. In these cases, you can get a simplified model, which is a linear approximation in a working point, for small increases in the variable around that point (x0 , y0 ). As you can see in Figure ref Picture9, linearization elsewhere (x1 , y1 ), will generally be different. If we consider a small variation around the working or equilibrium point , ex-

2.4. Linearization

25

panding in Taylor series 1 d2 y dy 1 dm y 2 (x − x0 ) + y = y0 + (x − x0 ) + . . . + (x − x0 )m + . . . dx 0 2! dx2 0 m! dxm 0 (2.4.1) The best linear approximation is the tangent line, corresponding to the expansion for m = 1 dy y = y0 + (x − x0 ) (2.4.2) dx 0 For multivariate functions, the linear approximation is ∂y ∂y ∂y (x1 − x10 ) + (x2 − x20 ) + · · · + (xn − xn0 ) (2.4.3) y = y0 + ∂x1 0 ∂x2 0 ∂xn 0 and corresponds to a plane or a hyperplane. Example 2.4.1: Obtain the transfer function VL (s)/V (s) of the nonlinear electrical network figure for small signs v(t), containing a non-linear resistor whose

i(t)

1H

v (t) L

+ v(t) -

Resistencia no lineal 20V

voltage-current relationship is defined byir = 2e0.1vr , where ir y vr are the current and voltage of the resistor, respectively. If we consider the system at rest or in equilibrium, the small signal v(t) = 0 and all derivatives are zero, i.e. VL (t) = L di(t) dt = 0. Then, at the equilibrium point, I0 = 2e0.1·20 = 14.8. The linearized equation of resistance is ∂2e0.1vr ir = ir0 + (vr − vr0 ) = 14.8 + 2 · 0.1e0.1vr0 (vr − vr0 ) ∂vr 0 = 14.8 + 0.2e2 (v − 20) = 14.8 + 1.48(vr − 20) Considering increments around the equilibrium point ∆(ir (t)) = 1.48∆(vr (t))

26

2. Classical techniques for system modeling

and the Laplace transform is ir (s) = 1.48Vr (s). Once linearized the nonlinear element, we can write the equations of the circuit ∆v(t) = L∆

di(t) + R∆i(t) dt

where R = 1.48 is the value obtained in the linearization. The Laplace transform is V (s) = VL (s) + Vr (s) = LsI(s) + RI(s) The output is VL (s) = LsI(s). Dividing the output from the input, we obtain the transfer function LsI(s) Ls s VL (s) = = = V (s) LsI(s) + RI(s) Ls + R s + 1.48

2.5

Block diagrams

The block diagrams are a way to represent systems that are composed of several parts which are represented by blocks and which are linked together by arrows depicting the flow of signals. The aim is to represent as blocks the transfer functions of each part, because it is easier, and then simplify the diagram to obtain the transfer function of the complete system.

2.5.1

Diagrams reduction

a) The output of a system is the product of the transfer function from the input (this is true only in the frequency domain). b) The transfer function of two blocks in series is the product of the transfer functions of the blocks. c) The transfer function of two blocks in parallel is the sum of the transfer functions of the blocks. d) The transfer function of a negative feedback loop is the transfer function of the direct line divided by one plus the transfer function of the direct line and the feedback.

2.5. Block diagrams

27

Fig. 2.11: Diagrams reduction

2.5.2

Transposition of sum points

Example 2.5.1: Consider the block diagram of Figure ref BL3a. To simplify the diagram the steps include: 1. Simplification of the negative feedback loop consisting of H1 (s) and H2 (s) for the diagram in Figure ref BL3b. 2. Transposition of the bifurcation of the above line , you will H6 (s), for the diagram in Figure ref BL4a. 3. Replacing the negative feedback loop for the TF and the corresponding sum of the signals in parallel on the right side as seen in the diagram of Fig ref BL4a. 4. Finally, the product is both TF and it is simplified as shown in Figure ref BL4b.

28

2. Classical techniques for system modeling

Fig. 2.12: Transposition of sum points

Fig. 2.13: a) Starting system of the example, and b) first simplification: feedback loop

2.5. Block diagrams

29

Example 2.5.2: For each of the systems given in the sketches below, find the transfer functions.

Solution: a) Using the 2nd Newton’s Law J θ¨ = τ − mg · sin θ · l where J = ml2 and τ is the external applied torque. The linearization is trivial because sin θ ≈ θ, when θ ≈ 0. Applying Laplace transform ml2 s2 Θ(s) + mglΘ(s) = τ (s) The transfer function is

Θ(s) 1 = 2 2 τ (s) ml s + mgl

b) Aplying the 2nd Newton’s Law J θ¨ = τ (t) − kθ − bθ˙

30

2. Classical techniques for system modeling

Applying Laplace transform Js2 Θ(s) = τ (s) − KΘ(s) − bsΘ(s) The transfer function is Θ(s) 1 = 2 τ (s) Js + bs + K c) The equations of the system are I(s)=i1 (s) + i2 (s) V (s)=Ri1 (s) 1 V (s)= i2 (s) Cs considering i1 the intensity trough the resistor and i2 the the intensity trough the capacitor. Solving for the intensities I(s)=i1 (s) + i2 (s) R i1 (s)= V (s) i2 (s)=CsV (s) Substituting I(s) =

1 V (s) + CsV (s) R

and the transfer function is V (s) = I(s)

1 R

1 R = RCs + 1 + Cs

d) Considering the position u(t) as input and y(t) as output, when the suspension system is in steady state and suddenly a small positive input is applied, the spring and the damper are compressed and it appear in the mass two forces: Fs = K(u − y) and Fb = b(u˙ − y). ˙ Using the 2nd Newton’s Law m¨ y = K(u − y) + b(u˙ − y) ˙ The Laplace transform is ms2 Y (s) + bsY (s) + KY (s) = bsU (s) + KU (s) and the transfer function Y (s) bs + K = 2 U (s) ms + bs + K

2.5. Block diagrams

31

Example 2.5.3: For raising and lowering a mass m = 1 kg, it is used a pulley with a radius r = 0.5 m and moment of inertia J = 0.1 kg ∆m2 powered by a DC motor whose characteristics are: • Internal resistance Ri = 70 Ω • Inductance Li = 10 H • Torque constant Ki = 2 N · m/A • Constant of electromotive force Ke = 6 V · s

The vertical displacement of the mass is proportional to the position of the pulley and is measured using the dynamic behavior which is given by the position sensor expression: r 2θ (t) us (t) = −5 cos(θ(t)) π Measuring voltage us (t) is compared with the reference voltage uref (t) with which the system is governed, being the gain of the amplifier Ku = 1. Find: • 1. Linearize the system around the equilibrium point defined by x0 = π4 . • 2. Draw the block diagram of the above system having as signal input Uref (t) and as output the weight displacement x(t). • 3. Calculate the transfer function of the system

X(s) Uref (s)

• 4. The displacement changes of the weight if the voltage value becomes uref = 180.5 V . • 5. If you want to convert the weight that drives the motor to a variable value (Such as what occurs in an elevator, when you do not know the accurate number people entering), where the disturbance signal has to be introduced the block diagram?

32

2. Classical techniques for system modeling

Fig. 2.15:

Example 2.5.4: In the figure it is shown a system that controls the height of liquid in a tank h(t) depending on the reference signal, r(t).

The float sensor is a negligible dimensions weight attached to a negligible weight rod whose end is the cursor of a potentiometer of 1/4 in circumference. Tension rod end, b(t), is compared to the reference, r(t), and the difference is amplified with a gain A = 20. This tension regulates the inlet valve whose equation is: Q1 (t) = Kv · [50v(t)] with 0.5Kv =

m3 s·V

2.5. Block diagrams

33

p √ 5 Data: K = 5 2 ms2 (coefficient of Torricelli ’s equation q2 (t) = K h(t), drain of the tank). a) Write the equations that represent the physical model of the system. b) Linearize around the equilibrium point given by H0 = 0.5m. c) Draw a representative block diagram of the system. d) Calculate the transfer function H(s)/R(s). e) Calculate the final value of h(t) when the input undergoes a change abrupt value of 1. f) Does the liquid overflow the tank in this case? Resolution: a) Write the equations that represent the physical model of the system. • The amplifier output is: V (t) = A[r(t) − b(t)] • The valve regulates the inflow according to the relationship: Q1 (t) = Kv [50 − v(t)] • The equations of deposit in the figure with section S = 1 m3 are: q1 (t) − q2 (t) = S q2 (t) = K

dh(t) dt

p h(t)

√ 5 with K = 5 2 ms2 . • The potentiometer equation is: b(t) = Vmax

θ(t) θmax

with Vmax = 20 V y θmax = π/2. • The angle is related to the height of liquid: hmax − h(t) = l cos θ(t) with hmax = 1m y l = 1m.

34

2. Classical techniques for system modeling

Fig. 2.16: • Finally we write the system equations: v(t) = 20[r(t) − b(t)] 1 q1 (t) = [50 − v(t)] 2 dh(t) q1 (t) − q2 (t) = √ p dt q2 (t) = 5 2 h(t) cos θ(t) = 1 − h(t) θ(t) b(t) = 40 π

b) Linearize around the given equilibrium point H0 = 0.5 m. The equilibrium condition is:

h0 = 1 rad θ0 = arccos(1 − h0 ) =

π rad 3

40 θ0 = 13.3 V πp = 5 2h0 = 5 m3 /s

b0 = q2,0

q1,0 = q2,0 v0 = 50 − 2q1,0 = 40 V 1 r0 = v0 + b0 = 15.3 V 20 The linearized system around the equilibrium condition is:

2.5. Block diagrams

35

∆v = A(∆r − ∆b) = 20(∆r − ∆b) 1 ∆q1 = −Kv ∆v = − ∆v 2 ∆q1 − ∆q2 = S∆h˙ = ∆h˙ K ∆q2 = √ ∆h = 5∆h 2 h0 √ 3 ∆h = l sin θ0 ∆θ = ∆θ 2 ∆θ 40 ∆b = Vmax = ∆θ θmax π c) Draw a representative block diagram of the system. To draw the block diagram, we must apply the Laplace Transform to the linearized equations: V (s) = A(R(s) − B(s)) Q1 (s) = −Kv V (s) Q1 (s) − Q2 (s) = SsH(s) K Q2 (s) = √ H(s) 2 h0 H(s) = l sin θ0 Θ(s) Vm ax B(s) = Θ(s) θm ax(s) The block diagram is:

d) Calculate the transfer function H(s)/R(s). H(s) C · G(s) = R(s) 1 + C · G(s) · Z

36

2. Classical techniques for system modeling

where C = −Kv A = −10 G(s) = Z=

m3 s·V

H(s) = Q1 (s) 1+

1 sS 1 √ K sS 2 h0

=

1 s

1+

5 s

=

1 s+5

B(s) 80 B(s) Θ(s) Vmax 1 = √ = = H(s) Θ(s) H(s) θmax sin θ0 π 3

Finally: 10 − s+5 −10 H(s) = 80 10 = √ √ R(s) 1 − π 3 s+5 s + 5 − π800 3

e) Calculate the final value of h (t) when the input suffers a sudden change of value 1. Applying the final value theorem, it is obtained: ∆h∞ = lim ∆h(t) = lim sH(s) = lim s t→∞

s→0

= lim s s→0

s→0

−10 R(s) = √ s + 5 − π800 3

−10 1 = −0.0704 m 800 s + 5 − π √3 s

f) It is overflowing the fluid from the reservoir in this case? The height reached by the liquid is: h∞ = h0 + ∆h∞ = 0.5 − 0.0704 m = 0.5704 m < hmax = 1 m So the liquid does not overflow. Example 2.5.5: It a furnace in which the heat, generated by an internal resistance qg (t), is proportional (K1 = 10) to the square of the current flowing through it (q(t) = K1 i2 (t)), while the lost heat qp (t) is proportional (K2 = 2.5) to the difference between the temperature inside the oven θi (t) and the external temperature θe (t). Knowing that the variation of the internal furnace temperature is proportional (K3 = 1) to the difference between the heat generated and wasted heat, and the system control is performed by an amplifier which generates a current i(t) = θref (t) − θi (t) + 5. Find: a) Equations of the system. b) Linearize the system around the equilibrium condition θi0 = θref 0 = 100o C y θe0 = 0o C

2.5. Block diagrams

37

c) Draw the block diagram of the device. Inputs: ∆θref (s), ∆θe (s). Output: ∆θi (s). Resolution: a) System equations. q(t) = K1 i2 (t) qp (t) = K2 [θi (t) − θe (t)] dθi (t) = K3 [qg (t) − qp (t)] dt i(t) = θref (t) − θi (t) + 5 with K1 = 10, K2 = 2.5, and K3 = 1. b) Linearize the system around the equilibrium condition θi0 = θref 0 = 1000 C y θe0 = 00 C. Initial conditions: qg0 = K1 i20 qp0 = K2 [θi0 − θe0 ] 0 = K3 [qg0 − qp0 ] i0 = θref 0 − θi0 + 5 namely i0 = 5 qg0 = 250 qp0 = 250

38

2. Classical techniques for system modeling

The linearized system around the equilibrium point is ∆qg (t) = 2i0 K1 ∆i(t) ∆qp (t) = K2 [∆θi (t) − ∆θe (t)] ∆θ˙i (t) = K3 [∆qg (t) − ∆qp (t)] ∆i(t) = ∆θref (t) − ∆θi (t) c) Draw the block diagram of the device. Inputs: ∆θref (s), ∆θe (s). Output: ∆θi (s) Laplace transforms are: ∆Qg (s) = 2i0 K1 ∆I(s) ∆Qp (s) = K2 [∆Θi (s) − ∆Θe (s)] s∆Θi (s) = K3 [∆Qg (s) − ∆Qp (s)] ∆I(s) = ∆Θref (s) − ∆Θi (s) The block diagram is:

Example 2.5.6: Is a physical system consisting of a power amplifier, a motor, a gear train and an encoder as shown in the figure:

2.5. Block diagrams

39

The voltage that reaches the power section is proportional to the difference between the reference voltage Ur and the output voltage of the encoder Ur , where the constant of proportionality Kp = 2. The power section produces a current i from the input voltage i, according to the following differential equation: p d i(t) + i(t) − 2 e(t) + 1.75 = 0 dt In the engine, the torque τ1 is proportional to the intensity i , being the proportionality constant Km = 4 N m/A. The motor shaft is connected to a reducer, being N1 = 16 teeth, N2 = 32 teeth. On the second axis we have coupled a flywheel with a moment of inertia I = 10 Kg m2 , a friction B = 5 N ms/rad, showing the axis an elastic behavior with K = 2N m/rad. Finally, the angle rotated by the second axis, θ2 , is measured with an encoder, that measures the output voltage Us iand it is proportional to the rotated angle θ2 , being the proportionality constant Ke = 0.5 V /rad. Find: a) System equations. b) Linearized around equilibrium point given by τ20 = 2N m. c) Blocks diagram. Us (s) d) Find the relation U . r (s) e) Find the output voltage of the encoder Us when the reference voltage experiences a sharp jump of two units. I(t)

Solution: a) System equations. The equations of the system are as follows Amplifier: e(t) = Kp [Ur (t) − Us (t)], con Kp = 2 Power section: i(t)

p di(t) + i(t) − 2 e(t) + 1.75 = 0 dt

Motor: τ1 (t) = Km i(t),

con Km = 4N m/A

Reductor: τ1 (t) N1 = , τ2 (t) N2

con

N1 = 0.5 N2

Flywheel I

d2 θ2 (t) dθ2 (t) = τ2 (t) − Kθ2 (t) − B 2 dt dt

40

2. Classical techniques for system modeling

con I = 10Kgm2 , B = 5N ms/rad y K = 2N m/rad Encoder Us (t) = Ke θ2 (t),

con Ke = 0.5V /rad

b) Linearize around equilibrium point given by τ = 220N m. We consider the equilibrium point: τ20 = 2N m 1 τ10 = τ20 N N2 = 1 N m τ10 i0 = Km = 0.25 A  2 e0 = i0 +1.75 = 1V 2 The only nonlinear equation is the one that characterizes the power section. So: 1 i0 ∆i˙ + ∆i − √ ∆e = 0 e0 The equations of changes are: ∆e = Kp (∆Ur − ∆Us ) i0 ∆i˙ + ∆i − √1e0 ∆e = 0 ∆τ1 = Km ∆i N2 ∆τ ∆τ2 = N 1 I∆θ¨2 = ∆τ2 − K∆θ2 − B∆θ˙2 ∆Us = Ke ∆θ2 Applying the Laplace transform E(s) = Kp [Ur (s) − Us (s)] i0 sI(s) + I(s) − √1e0 E(s) = 0 τ1 (s) = Km I(s) 2 τ2 (s) = N N 1 τ1 (s) 2 Is Θ(s) = τ2 (s) − KΘ(s) − BsΘ(s) Us (s) = Ke Θ2 (s) The blocks diagram is:

2.5. Block diagrams

41

c) Find the relation

Us (s) Ur (s)

Us (s) G(s) = Ur (s) 1 + G(s) where G(s) =

Us (s) 1 1 N2 1 16 = Kp √ Km Ke = U (s) e0 i0 s + 1 N 1 Is2 + Bs + K 2.5s3 + 11.25s2 + 5.5s + 2

So Us (s) 1 = Ur (s) 2.5s3 + 11.25s2 + 5.5s + 18 d) Find the output voltage of the encoder Us when the reference voltage experiences a sharp jump of two units. Applying the final value theorem is obtained:

∆Us,∞ = lim ∆Us (t) = lim sU s(s) = lim s t→∞

s→0

= lim s s→0

2.5s3

s→0

2.5s3

16 Ur (s) = + 11.25s2 + 5.5s + 18

2 16 = 1.778 V 2 + 11.25s + 5.5s + 18 s Us,∞ = Us,0 + ∆Us,∞

Us0 = Ke θ20 =

Ke τ20 = 0.25V K

So Us,∞ = Us,0 + ∆Us,∞ = 2.028 V Example 2.5.7: A thermistor has a response to the temperature represented by R = R0 e−0.1T , where R = resistance and T = temperature in Celsius. Find the linear model for the thermistor operating around T = 20o C with a small range of temperature variation. Solution: In the balance point or point of operation T = 20o C Req = R0 e−0.1·20 = 0.135R0 Linearizing around this point of equilibrium ∂R0 e−0.1T ∆R = ∆T = R0 (−0.1)e−0.1T T =20 ∆T = −0.0135∆T ∂T T =20

42

2. Classical techniques for system modeling

By taking Laplace transforms R(s) = −0.0135T (s) Example 2.5.8: Determine the transfer function H(s)/Qe (s) of the system consisting of a water reservoir, linearizing around the point of equilibrium given by qe0 = 1m3 /s. if the section of the outlet pipe is A = 1m and the area of the base is B = 1m2 . √ Data: Consider the Torricelli’s equation qs (t) = K · A h.

Solution: The inflow minus outflow is the variation in volume qe (t) − qs (t) = B

dh(t) dh(t) = dt dt

According to Torricelli’s equation √ qs (t) =

h

In the equilibrium point, qe0 = 1m3 /s qe0 =qs0 ⇒ qs0 = 1m3 /s √ 1 qs0 =K h0 ⇒ h0 = m K

(2.5.1)

Linearizing ∆qe − ∆qs =∆h˙ K K 3/2 ∆qs = √ ∆h = ∆h 2 2 h0

(2.5.2)

2.5. Block diagrams

43

By taking Laplace transforms Qe (s) − Qs (s)=sH(s) Qs (s)=

K 3/2 H(s) 2

Sustituting Qe (s) − The transfer function is

K 3/2 H(s) = sH(s) 2

1 H(s) = 3/2 Qe (s) s + K2

(2.5.3)

44

2. Classical techniques for system modeling

Bibliografia [1] E.A. Puente, Regulacion Automatica I, Seccion de Publicaciones ETSIIM, 1980. [2] K. Ogata, Modern Control Engineering, Prentice Hall International Edition, 1995. [3] B. Kuo, Automatic Control Systems, Prentice Hall, 1996.

3. TIME ANALYSIS 3.1

Introduction. Time Analysis of Systems

The previous chapters have studied the mathematical modeling on the basis of their relationship between input and output of the systems. This leads us to the concept of transfer function. In this chapter we analyze how the systems behave using their transfer functions.

3.2

Temporal Response

When the input signal u(t) is introduced into a system transfer function G(s), its output is Y (s) = G(s)U (s) ⇒ y(t) = g(t) ∗ u(t) = L−1 (G(s)U (s)) This time response has two parts: the transient response and steady state response. Transient response is part of the temporal response that tends to zero as time goes to infinity and the steady-state response is that which remains over time. y(t) = yrt (t) + yrp (t) lim yrt (t) = 0 t→∞

yrp (t) = lim y(t) t→∞

The dynamic behavior of the system is primarily related to transient response and the speed and stability can be studied. In the steady-state response the accuracy of the response can be studied, ie the stationary system error. To analyze the response, different test signals are introduced. The most important are: The signal impulse or Dirac delta u(t) = δ(t)⇒U (s) = 1 The unit step signal  u(t) = u0 (t) =

0,t < 0 ⇒U (s) = 1,t > 0

1 s

46

3. Time Analysis

Fig. 3.1: Response of a system and the ramp signal unit u(t) = tu0 (t)⇒U (s) =

3.2.1

1 s2

Impulse response

The main reason for the classical control theory uses the frequency domain rather than the time domain is given by the convolution theorem. In the time domain, the output is the convolution of the input and the impulse response of the system Z

t

Z

t

g(τ )u(t − τ )dτ

g(t − τ )u(τ )dτ =

y(t) = g(t) ∗ u(t) =

0

0

while in the frequency domain, the output is a product Y (s) = G(s)U (s). In the event that the input is an impulse unit u(t) = δ(t), its Laplace transform is the corresponding signal in the frequency domain, U (s) = 1. Z y(t) = g(t) ∗ δ(t) =

t

g(t − τ )δ(τ )dτ = g(t) 0

That is, when the input signal is an impulse, the output is y(t) = g(t), where g(t) = L−1 (G(s)) The function g(t) is called the impulse response of the system and is a mathematical model of the system, as can be the differential equations or transfer function G(s).

3.2. Temporal Response

3.2.2

47

Response of First Order systems

The First Order systems are the simplest and are given by the differential equation T y(t) ˙ + y(t) = Ku(t)

(3.2.1)

its Laplace transform with zero initial conditions is T · sY (s) + Y (s) = KU (s) and its transfer function can be written as Y (s) K = U (s) 1 + Ts The constant K is the static gain, T is the time constant and σ = of attenuation. It is interesting to note that the pole is located at −1/T .

1 T

is the constant

Impulse response Consider a first order system transfer function G(s) = unit impulse u(t) = δ(t) ⇒ U (s) = 1, the output is Y (s) = G(s) · 1 =

K 1+T s .

If the input is a

K 1 + Ts

and from the time domain is y(t) = g(t) = L−1



K 1 + Ts

 =

k −t/ e T T

t≥0

Fig. 3.2: Impulse response of a first order system As shown in Fig. 3.3, the constant T is the time it takes the system to reach the 37% of the initial value, g(0) = K/T .

48

3. Time Analysis

Fig. 3.3: Impulse response of a first order system Step response If the input signal to the system G(s) is a unit step u(t) = u0 (t) ⇒ U (s) = 1/s

Fig. 3.4: Unit step response of a first order system then output in the time domain is      k −1 −1 1 − y(t) = L−1 G(s) = L = kL s s s(1+T s) t = k(1 − e− /T ) ,

T 1+T s



=

when t ≥ 0.

In figure Fig. 3.5 the output of a first-order system to a step input unit is shown. The time constant T is the time it takes for the system to reach the 63% of the final value K and measures the speed of system response. The time it takes the system to reach the 95% of the final value is called settling time and can be calculated as ts ≈ 3T ≈ πσ . Example 3.2.1: It is also interesting to note that from the response unit step is possible to determine the parameters K (final value) and T (time it takes for the system to reach the 63% of the final value) and therefore , you can determine the transfer function. For example if the unit step response is shown in Fig. 3.6, the final value K = 0.72 (which is also the gain because the input is a unit step) and the time

3.2. Temporal Response

49

Fig. 3.5: Unit step response of a first order system

Fig. 3.6: Determination of the transfer function from the unit step response of a first order system constant is corresponding to the time at which 0.63 · 0.72 , which is T = 0.13, and the transfer function is K 0.72 G(s) = = 1 + Ts 1 + 0.13s Response to ramp If the system input G(s) is the unit ramp u(t) = t · u0 (t) ⇒ U (s) = 1/s2 the output y(t) can be expressed as

50

3. Time Analysis

Fig. 3.7: Response to a ramp unit of a first order system

y(t) = L−1



G(s) s2



= L−1



k s2 (1+T s)



 = KL−1 − Ts +

t t = K(−T + t + T e− /T ) = K(t − T ) + kT e− /T ,

1 s2

+

T2 1+T s



=

para t ≥ 0

In this expression, T is the delay between the ramp input and the response.

Fig. 3.8: Unit ramp response of a first order system

3.2.3

Response of the second order systems

The standard equation of a second order system is G(s) =

Kωn 2 s2 + 2ζωn s + ωn 2

where K = Static gain. ωn = Undamped natural frequency. ζ = Damping ratio. σ=ζωn = p Attenuation. ωd = ωn 1 − ζ 2 = Damped frequency.

(3.2.2)

3.2. Temporal Response

51

The system behavior is determined by K, ζ y ωn . The poles can be determined from the characteristic equation (denominator of G(s) equated to zero). s2 + 2ζωn s + ωn 2 = 0

(3.2.3)

p p s = −ζωn ±p ζ 2 ωn 2 − ωn2 = −ζωn ± ωn ζ 2 − 1 = −ζωn ± jωn 1 − ζ 2

(3.2.4)

s = −σ ± jωd

(3.2.5)

The poles are

Fig. 3.9: Relation among the parameters of a second order system.

52

3. Time Analysis

Fig. 3.10: Classification of second-order systems Impulse response Consider a second order system G(s) =

kKωn 2 s2 + 2ζωn s + ωn 2

(3.2.6)

and assume that your input is a unit impulse u(t) = δ(t) ⇒ U (s) = 1 The most interesting case is when the system is underdamped 0 < ζ < 1. In

3.2. Temporal Response

53

this case the poles are complex conjugates with negative real part. p −ζωn ± jωn 1 − ζ 2 = −σ ± jωd The response is Kωn −σt e sin(ωd t + ϑ), y(t) = p 1 − ζ2

para t ≥ 0

That is a damped sine with envelope e−σt . The decay factor σ = ζωn indicates the rate at which decreases the amplitude of the oscillation.

Fig. 3.11: Impulse response of the second order underdamped systems.

Step response Consider a system with transfer function G(s) and step input unit, G(s) s

Y (s) = U (s)G(s) = y(t) = L

−1



G(s) s



Z =

t

g(τ )dτ 0

If the denominator is decomposed into simple factors, which may correspond to simple roots or pairs of complex conjugate roots. The output to a unit step is G(s) = Q  Y (s) =

N (s) 2

(s − αi ) + βi 2

Q

(s − σi )

G(s) N (s)  = = Q Q s s (s − σi ) (s − αi )2 + βi 2

54

3. Time Analysis

X A X Bi Ci s + D i + + s s − σi (s − αi )2 + βi 2 X X y(t) = L−1 (Y (s)) = A + Bi eσi t + Ei eαi t sen(βi t + ϕi ) =

y(t) = G(0) +

X

Bi eσi t +

X

Ei eαi t sen(βi t + ϕi )

t≥0

t≥0

(3.2.7)

The output consists of three parts: the constant A = G(0) is due to the input and it is the forced response and the stationary system response gain and the second part X X Bi eσti + Ei eαi t sin(βi + tϕi ) which is the natural response of the system, or the part that it is due to the transfer function. In it, the expression X Bi eσti corresponds to the real roots of the denominator (real poles) and X Ei eαt i sin(βi + tϕi ) corresponds to pairs of complex conjugate roots of the denominator (complex conjugate poles). It is interesting to note that if any of the coefficients σi or αi of the real exponential is positive, the corresponding exponential is growing and the response y(t) tends to infinity, ie, the system is unstable. A dynamic system is stable if and only if all poles are located in the left half semi plane. If the poles are real, ie the coefficients Ei are zero, and the response X y(t) = A + Bi eσi t consists of two parts: a constant and an exponential that is increasing if σi is positive, and it is decreasing if σi is negative, and one if σi = 0. Adding the constant A and the part of the exponential, we have the responses shown in Fig. 3.14. The constant A is the steady-state response yrp (t) = lim y(t) = A t→∞

and according to the final value theorem coincides with G(0): lim y(t) = lim s

t→∞

s→0

G(s) = G(0) s

3.2. Temporal Response

55

Fig. 3.12: The real exponential function eσt is increased when the real coefficient σ is positive and decreasing when it is negative.

Fig. 3.13: The function eαt sin(βt + φ) is a growing sinusoid if α is positive and decreasing if negative.

56

3. Time Analysis

Fig. 3.14: The unit step responses of a second order system according to the coefficient of friction: 1) ζ = 0, oscillator (poles on the imaginary axis). 2) 0 < ζ < 1, underdamped system (complex conjugate poles in the left half. 3) ζ = 1, critically damped system (real and equal poles). 4) ζ > 1, overdamped system (real and different poles).

3.2. Temporal Response

57

Example 3.2.2: Consider the system of Fig.3.2.3, but with the parameters: m=1, b=3, K=2. a) Find the transfer function G(s) = X(s)/F (s). b) Find the response to unit step of the system. The differential equation system is m¨ x + bx˙ + kx = f (t) By taking Laplace transforms (ms2 + bs + k)X(s) = F (s) The transfer function is X(s) 1 1 = = 2 2 F (s) ms + bs + k s + 3s + 2 b) If the input is a unit step F (s) = X(s) =

1 s

and output is

1 1 s s2 + 3s + 2

Decomposing into partial fractions X(s) =

1 1 A B C = + + 2 s s + 3s + 2 s s+1 s+2

where A = lim

s→0 s2

1 = 0.5 + 3s + 2

1 1 = −1 s→−1 s s + 2

B = lim C = lim

s→−2

1 1 = 0.5 ss+1

58

3. Time Analysis

namely 0.5 −1 0.5 + + s s+1 s+2 Taking inverse Laplace transforms, we have X(s) =

x(t) = 0.5 − e−t + 0.5e−2t para t > 0 Example 3.2.3: 1.- For the system of Fig. 3.2.3 a) Find the transfer function G(s) = X(s)/F (s). b) Find the unit step response of the system. a) The differential equation of the system is m¨ x + bx˙ + kx = f (t) By taking Laplace transforms (ms2 + bs + k)X(s) = F (s) The transfer function is X(s) 1 1 = = 2 2 F (s) ms + bs + k 3s + 15s + 33 b) If the input is a unit step F (s) = X(s) =

1 s

and output is

1 1 1 1 = 2 2 s 3s + 15s + 33 3 s(s + 5s + 11)

Decomposing into partial fractions X(s) =

1 A Ms + N = 3 s s2 + 5s + 11

 1/11 −(1/11)s − 5/11 + s s2 + 5s + 11   1 1 s+5 X(s) = − 33 s (s2 + 2.5s)2 + 2.182 1 = 3



This last expression can be written as the sum of the transformed damped sine and cosine   1 1 s + 4.75 0.25 − − X(s) = 33 s (s2 + 2.5s)2 + 2.182 (s2 + 2.5s)2 + 2.182

3.3. Parameters characterizing the response of a second order system or specifications 59

1 33

X(s) =



1 0.25 s + 2.182 − − 2 s (s + 2.5s)2 + 2.182 (s2 + 2.5s)2 + 2.182



The temporal response is the inverse transform x(t) =

3.3


1 1 − e−2.5t cos(2.18t) − e−2.5t sin(2.18t) for t > 0 33

Parameters characterizing the response of a second order system or specifications

The specifications express the characteristics to be achieved in the system when a controller is introduced into it. They can be expressed in different ways by using the frequency domain or the time domain. In our case we will refer to specifications in the time domain. Specifications regarding the answer in the time domain of a system are usually defined with respect to the temporary response to a step input of a second order system. To do so a reference system is taken, with a transfer function: G(s) =

Kωn2 s2 + 2ζωn s + ωn2

(3.3.1)

whose unit step response is p 1 y(t) = 1 − p e−ζωn t sin(ωn t 1 − ζ 2 + φ) 1 − ζ2

(3.3.2)

with p 1 − ζ2 φ = arctan ζ This system has its poles located in p1,2 = −ζωn ± jωn

p 1 − ζ2

(3.3.3)

(3.3.4)

where ζ is called damping ratio damping ratio and ωn is the undamped natural frequency of the system. For a system with damping rate 0 < ζ ≤ 1 (underdamped system) the system response curve has the form shown in Fig.7.3 Among the parameters that can be used to specify the response, we have the following: • Peak time, Tp . Tp =

π ωd

where ωd is the damped frequency ωd = ωn

(3.3.5) p 1 − ζ 2.

60

3. Time Analysis

y

To Mp

d·Mp

ep

Tp

t Ts

Fig. 3.15: Parameters of the temporal response • Overshoot, Mp . −π √

Mp = e • Settling time, Ts .

Ts '

√

ζ 1−ζ 2

π ζωn

=

d

(3.3.6)

(3.3.7)

Finally, the steady state position error, Ep , is defined as the difference between the desired reference value and the output value steady state of the controlled variable of the system (we are assuming for simplicity that the feedback signal is unitary). Example 3.3.1: Given the following circuit, find the maximum voltage on the capacitor when closing the switch:

Solution: The system equations are calculated, considering as input the battery-switch voltage (0 to 10 V, depending on the switch is open or closed) and as output the voltage across the capacitor: Z di(t) 1 + i(t)dt vi = vR + vL + vC = Ri(t) + L dt C

3.3. Parameters characterizing the response of a second order system or specifications 61

vo =

Z

1 C

i(t)dt

Applying the Laplace transform Vi (s) = RI(s) + LsI(s) + Vo (s) =

1 I(s) Cs

1 I(s) Cs

Dividing the transfer function is obtained 1 Vo (s) Cs I(s) = Vi (s) RI(s) + LsI(s) +

=

1 Cs I(s)

1 Cs

R + Ls +

1 Cs

=

LCs2

1 + RCs + 1

With the values of the statement you have 108 Vo (s) = 2 Vi (s) s + 12 · 103 + 108 As seen is a second order system. Compared with the standard expression Kωn2 s2 + 2ζωn s + ωn2 it is obtained ωn2 =108 ⇒ ωn = 104 Kωn2 =108 ⇒ K = 1 2ζωn =12 · 103 ⇒ ζ =

12 · 103 = 0.6 ⇒ Underdampened. 2 · 10

(3.3.8)

The maximum voltage is reached when the overshoot is maximized Mp = e

− √ πζ

1−ζ 2

= 0.0948 = 9.48%

When the input is a step of 10V , then Vo,max = 10 + 10 ∗ 0.0948 = 10.948V Example 3.3.2: Compare these systems in terms of overshoot and settling time, drawing responses from the closest possible way: Ga = Gc =

10 s2 +2s+5 5 (s2 +2s+5)(s+0.5)

Gb = Gd =

20 (s2 +2s+5)(s+2) 5(s+2) (s2 +2s+5)

62

3. Time Analysis

Solution: a) First system

Ga =

10 s2 + 2s + 5

3.3. Parameters characterizing the response of a second order system or specifications 63

b) Comparison among the systems:

64

3. Time Analysis

• The Gb (s) system equals Ga (s) with a further (farther from the imaginary axis) pole. • The Gc (s) system equals Ga (s) with an additional influential pole on the case of Gb (s): the system will be much slower. • The Gd (s) equals Ga (s) system with an extra zero: the system will be faster. Looking at the graph can be concluded that: • Adding a zero to a system makes the peak occurs earlier and is greater. • The closer the zero will be to the imaginary axis, the faster the system is. • Adding a real pole system makes the peak occurs later and is minor. • The closer the pole to the imaginary axis is, the slower the system is. Example 3.3.3: A system consists of the following block diagram figure:

The response of each subsystem is shown in the following figure:

3.3. Parameters characterizing the response of a second order system or specifications 65

Calculate the final value when a 2 units step is applied. Solution: The subsystem response to a unit step input is:

It can be seen that: • That’s the typical response of an underdamped second order system. • It stabilizes around the value 1, thus the static gain of the system is one. • Reaches a maximum at 1.08 therefore Mp =

1.08−1 1

= 0.08.

• The stabilization time is ts = 3.14 With these data we will calculate the transfer function of the system. Compared with the standard expression

66

3. Time Analysis

Kωn2 s2 + 2ζωn s + ωn2 It can be obtained: • Static gain: K = 1. π σ

• Stabilization time: ts =

→σ=

π ts

= 1 → ζωn = 1.

• Damping ratio: Mp = e

π − tan θ



 −π → θ = atan = 0.8936rad. ln(Mp )

To calculate the damping ratio: ζ = cos(θ) = 0.6266 < 1 To calculate the undamped natural frequency: ζωn = 1 → ωn =

1 = 1.596rad/s ζ

The system is: GA =

2.5 s2 + 2s + 2.5

You can verify in this plot:

The impulse response of the system B is:

3.3. Parameters characterizing the response of a second order system or specifications 67

It can be seen that this is the typical response to a system of first order. GB (s) =

K 1 + Ts

The generic response to a first order system facing an input impulse is:

We note that y(0) = 2 →

K =2 T

T =1 With the data T = 1 and K = 2 we will calculate the transfer function of the system: 2 GB = 1+s

68

3. Time Analysis

Finally, the response of the system C to a ramp with slope 0.5 is:

That is the answer to a ramp of a system of first order:

Gc (s) =

K 1 + Ts

The response to a ramp with slope 0.5 has a value:

y(t) =

The generic response is

K K (t − T ) + T e−t/T 2 2

3.3. Parameters characterizing the response of a second order system or specifications 69

We note that T = 1, the output delay with respect to the line r(t) = output signal. The steady slope of the output signal is 1:

K 2 (t − T )

once stabilized

K =1→K=2 2 We can also verify that using the final value theorem:

The transfer function is: Gc (s) =

2 1+s

The complete system will be:

Example 3.3.4: If Vi (t) is a voltage step in the network shown in the figure, find the resistance value so that the overshoot is 20% in the voltage across the capacitor, if C = 10−6 F y L = 1H.

70

3. Time Analysis

Solution: The system equations in the frequency domain are 1 Vi (s)=LsI(s) + RI(s) + I(s) Cs 1 Vo (s)= I(s) Cs Dividing the transfer function is obtained 1 Vo Cs = Vi Ls + R +

1 Cs

=

s2 +

1 LC R 1 L s + LC

Compared to the standard form of a second order system Kωn2 s2 + 2ζωn s + ωn2 q q 1 , y ζ = R2 C Se obtiene que K = 1, ωn = LC L. Furthermore, the overshoot is required to be the 20% π

Mp = e− tanθ = 0.20 ⇒ θ = 62.370 ⇒ ζ = cosθ = 0.456 Considering C = 10−6 F y L = 1H, r ωn =

1 = 103 LC

R R = = 2ζωn = 2 · 0.456 · 103 = 912Ω ⇒ R = 912Ω L 1 Example 3.3.5: Consider the following RC-circuit implementation of a lowpass filter.

3.3. Parameters characterizing the response of a second order system or specifications 71

(a) The transfer function for this circuit is given below. Employing this model, find the output response of the circuit e(∞) for a step input of magnitude E: Eo (s) 1 = Ei (s) R1 C1 + 1 (b) If R1 = 5Ω and C1 = 2F , determine the settling time ts = σ4 for the filter, that is, the time it takes for the circuit to reach and stay within 2% of its final value. (c) An alternative implementation of a low-pass filter is given below. Find the transfer function for this circuit for an input of ei (t) and an output of eo (t) .

(d) You wish to design the filter in (c) to achieve similar performance to the filter in (a). Specifically, if L2 = 5H, then find the resistance R2 and capacitance C2 so that: i. The two filters have the same settling time ts (corresponds to the filter’s pass band). ii. The two filters reach the same steady-state value for a constant input. iii. The second filter has a maximum percent overshoot equal to 5%. Solution: a) e(∞) = E. 4 b) ts = σ4 = 0.1 = 40s. c) 1

Eo (s) C2 s = Ei (s) L2 s + R2 +

1 C2 s

=

s2 +

1 LC R 1 L s + LC

=

Kωn2 s2 + 2ζωn s + ωn2

d)  π  π Mp = exp − = 0.05 ⇒ tan θ = − = 1.0487 tan θ ln(0.05) θ = 46.36o = 0.81rad ⇒ ζ = 0.7 σ 0.1 ωn = = = 0.14 ζ 0.7 Using the value L = 5 R = 2ζωn ⇒ R = 2 · 0.7 · 0.14 · 5 = 0.98Ω L

72

3. Time Analysis

1 1 = ωn2 ⇒ C = = 10F LC 0.144 · 5

4. STABILITY OF DYNAMICAL SYSTEMS The first property they must have control systems is that they are stable. Moreover, in these systems, the closed loop operation is the risk that they are unstable in open loop systems were stable. The simplest concept of stability is the stability of input-output that is bounded input bounded output will correspond, i.e. the gain is always finite. If we consider the system of Figure 4.1, where G(s) and H(s) are linear or nonlinear, the error e is e = r + H(s)y = r + H(s)G(s)e. The norm error is given

r(t)

+

e(t)

G(s)

y(t)

-

H(s) Fig. 4.1: Feedback loop by the expression kek =

kuk kuk ≤ k1 − G(s)H(s)k 1 − kG(s)kkH(s)k

for any norm we take. If kG(s)kkH(s)k < 1 the system gain is finite and therefore the system is input-output stable. Theorem 4.0.1: Small gain theorem Consider the system of Figure 4.1. A sufficient condition for the system is input-output stable is that kG(s)kkH(s)k < 1.

(4.0.1)

In the case that the system is linear, it can require a little less demanding condition, it is sufficient that kG(s)H(s)k < 1. Of the numerous methods to study the stability of dynamical systems, they are going to be seen in this chapter some of the best known, both for linear and nonlinear systems.

74

4. Stability of dynamical systems

4.1

Classical stability analysis methods

In linear systems, the stability of the system is determined by the position of the closed-loop poles. Consider a system with an input and an output whose transfer function in a closed loop Qm (s + zi ) Qi=1 M (s) = K Qn p 2 2 (s + p ) j j=1 k=1 (s + 2ζk ωk s + ωk ) The impulse response is y(t) =

k X

Aj e

−pj t

+

j=1

p X

 Bk

k=1

1 ωk



e−ζk t sin(ωk t)

Therefore, for a bounded input bounded output to a need for both the real poles pj , as the real parts of the complex poles ζk are negative. The sufficient condition for a linear SISO system (Single Input-Single Output) in continuous time, is stable is that All of the poles are located in the left half-plane (negative real part). The problem with this approach is that it is very difficult to apply when parameters. In this case we must apply other criteria such as Routh.

4.1.1

Routh’s method

Consider a continuous-time system with characteristic equation

p(s) = an sn + an−1 sn−1 + an−2 sn−2 + . . . + a1 s + a0 = 0

(4.1.1)

As prerequisites we require that: 1. All coefficients have the same sign 2. No coefficient is zero. The Routh criterion is based on the construction of a table called Routh table and which is constructed in the following way: The first row is formed by the coefficients an , an−2 , an−4 . . . and the second row of the remaining an−1 , an−3 , an−5 , . . . . Namely sn sn−1

an an−1

an−2 an−3

The remaining rows are constructed as

an−4 an−5

... ...

4.1. Classical stability analysis methods

sn sn−1 sn−2 sn−3 ... s1 s0

an an−1 bn−1 cn−1 ... u1 v1

75

an−2 an−3 bn−3 cn−3 ... u2

an−4 an−5 bn−5 cn−5 ...

... ... ... ... ...

(4.1.2)

wherein the third row is

bn−1

bn−3

an an−2 an−1 an−3 1 =− = (an−1 an−2 − an an−3 ) an−1 an−1 an an−4 an−1 an−5 1 = (an−1 an−4 − an an−5 ) =− an−1 an−1

Analogously, the fourth row is

cn−1

cn−3

an−1 an−3 bn−1 bn−3 1 =− = (bn−1 an−3 − an−1 bn−3 ) bn−1 bn−1 an−1 an−5 bn−1 bn−5 1 =− = (bn−1 an−5 − an−1 bn−5 ) bn−1 bn−1

and so on down to zero degree. Note 4.1.1: In developing the table, we can multiply or divide a row by a positive constant without the result is changed (this is useful to simplify the table). The Routh criterion states that the number of roots in the right half-plane coincides with the number of sign changes of the coefficients of the first column. That is, the system will be stable if all the signs in the first column are equal. There are a number of remarkable special cases: There is a zero in the first column. The value 0 is replaced by a very small positive number ε and the development of the table continues. To count the number of sign changes in the first column is necessary to consider that ε goes to zero.

76

4. Stability of dynamical systems

There is a row of zeros. ( Poles symmetrically placed about the origin of the complex plane). The coefficients of the row above the row of zeros are considered. With these coefficients the auxiliary equation is constructed, which will be the grade that indicates the power of s respectively. The row of zeros is replaced by the coefficients of the derivative and the method continues. Example 4.1.1: Normal case. It is considered the system with characteristic equation a0 s3 + a1 s2 + a2 s + a3 = 0 with a0 , a1 , a2 y a3 positives. The corresponding Routh table is s3 s2 s1 s0

a0 a1

a1 a2 −a0 a3 a1

a2 a3

a3

For all coefficients in the first column be positive, shall be verified that a1 a2 − a0 a3 >0 a1 Therefore, the system will be stable when a1 a2 > a0 a3 . Example 4.1.2: Case with a zero in the first column. Consider the system with characteristic equation p(s) = s4 + 2s3 + 4s2 + 8s + 5 The corresponding Routh’s table is s4 s3 s2 s1 s0

1 2 0≈ε

4 8 5

5

8ε−10 ε

5

where the zero value of the first column by the small positive number is replaced ε.

When ε → 0, sign 8ε−10 = sign −10 = −1. Since there are two sign ε ε changes in the first column, there are two poles in the right half, and therefore the system is unstable.

4.1. Classical stability analysis methods

77

Example 4.1.3: Case with a zeros row. Consider a system with characteristic equation p(s) = s4 + 2s3 + 11s2 + 18s + 18 = 0. Solution: This system has its roots located symmetrically respect to the origin of the plane s: ±σ ± jω. Routh’s table is s4 s3 s2 s1 s0

1 2 2 0→4 18

11 18 18 0→0

18

The auxiliary equation associated with the third row is A(s) = 2s2 + 18 and its derivative is A(s) ds = 4s + 0 to replace the old values third row. Since all the coefficients of the first column are positive, the system is stable. Example 4.1.4: Consider the system of Fig. 4.2 Is desired to find the stability of r(t) +

-

e(t)

1 2 s +2s+5

K

y(t)

1 s+3 Fig. 4.2: System of example 5.1.4 the system as a function of the parameter K. Solution: The transfer function in closed loop is K

K(s + 3) KG(s) s2 +2s+5 Gc (s) = = = K 1 2 + 2s + 5) + K 1 + KG(s)H(s) (s + 3)(s 1 + s2 +2s+5 s+3 Equating the denominator to zero the characteristic equation is obtained, whose roots are the poles of the closed loop system: p(s) = (s + 3)(s2 + 2s + 5) + K = s3 + 5s2 + 11s + 15 + K = 0

78

4. Stability of dynamical systems

The precondition requires that all the coefficients are of the same sign: 15 + K > 0 ⇒ K > −15 Routh’s table is s3 s2 s1 s0

1 5

11 15 + K

40−K 5

15 + K

For the system to be stable it is necessary that all the coefficients of the first column be greater than zero: 40 − K > 0 ⇒ 40 − K > 0 ⇒ K < 40 5 Putting together the two conditions, we find that the system is stable for −15 < K < 40.

4.1. Classical stability analysis methods

79

Bibliografia [1] E. A. Puente, Regulacion Automatica Seccion de Publicaciones de la ETSIIM, 1993. [2] E. Umez-Eronini, Dinamica de Sistemas y Control, Thomson Learning, 2001. [3] S. Barnett y R.G. Cameron, Introduction to Mathematical Control Theory, Oxford University Press, 1992. [4] T. Glad y L. Ljung, Control Theory, Taylor and Francis, 2000. [5] JJ.E. Slotine y W.Li, Applied Nonlinear Control, Prentice Hall, 1991.

80

4. Stability of dynamical systems

5. ERRORS IN STEADY STATE 5.1

Introduction to errors

Consider the system errors to the standard inputs. The behavior of different systems with the same input signal can be very different. As an example, consider the behavior of an engine when you want to control it in position and speed: Control of the axis position of a DC motor

Fig. 5.1: By varying the reference signal: • varies the output of the adder (error) • is given to the motor voltage • The motor rotates • When the motor position coincides with the reference output the adder is zero • The engine is off power • (possible oscillations) • The motor is stopped at the position indicated The steady system output matches the reference signal

82

5. Errors in steady state

Speed control of a DC motor

Fig. 5.2: By varying the reference signal: • varies the output of the sum block (error) • voltage to the motor is given • The motor rotates • When the engine speed matches the reference output the sum is zero • The engine power is off • The engine brakes • varies the output of the sum • ... • The system will stabilize at a speed different from the reference (If the error were zero, it would stop) The steady system output does not match the reference signal These two examples show the different behavior of the two systems to the step input. While in the case of motor position control, the error in the steady may be zero, in the case of speed control, the error must be nonzero.

5.1.1

Gain position (static)

Constant position error is the value of the steady state output when the input is a unit step, i.e. Kp = lim y(t). t→∞ It is usually calculated as Kp = lim G(s) s→0

Is only meaningful in stable systems. Interpretation Kp = lim y(t) = lim s t→∞

s→0

G(s) = lim G(s) s→0 s

5.1. Introduction to errors

5.1.2

83

Gain speed

Error constant acceleration in steady state is the value of the first derivative of the output (slope) when time tends to infinity, i.e. Kv = lim y(t). ˙ It is usually t→∞ calculated as Kv = lim sG(s) s→0

Interpretation

Kv = lim y 0 (t) = lim s s t→∞

5.1.3

s→0

G(s) = lim sG(s) s→0 s

Gain acceleration

Error constant acceleration in steady state is the value of the second derivative of the output (curvature) when time tends to infinity, i.e. Kv = lim y¨(t). t→∞

It is usually calculated as Ka = lim s2 G(s) s→0

Ka = lim y 00 (t) = lim s s2 t→∞

s→0

G(s) = lim s2 G(s) s→0 s

84

5. Errors in steady state

5.1.4

Type system

Is named type of a system the number of integrators of the system (poles at the origin). That is, if the transfer function is factored as Q (s − zi ) G(s) = K r Q s (s − pi ) the system type is r. Examples G(s) = G(s) = G(s) = G(s) =

5.1.5

3 Type 0 (order 1) s+2 s+1 Type 1 (order 2) s(s+3) 4 Type 2 (order s2 (s+2)(s+1) 5 Type 3 (order 3) s3

4)

Relationship between type and gain of a system

The relationship between gain and type of system is:  Kp Type = 0 Kp = lim G(s) = ∞ Type > 0 s→0  0 Type = 0 Kv = lim sG(s) = Kv Type = 1 s→0  ∞ Type > 1  0 Type < 2 Ka = lim s2 G(s) = Ka Type = 2 s→0  ∞ Type > 2

Fig. 5.3: Relationship between type and gain a system

5.1. Introduction to errors

5.1.6

Relationship between type and gain a system

Fig. 5.4: Relationship between type and gain a system

85

86

5. Errors in steady state

5.1.7

Steady state error of a unit feedback system

Consider a unitary feedback system as shown on Fig. 5.5

Fig. 5.5: Unit feedback system. The error signal is e(t) = r(t) − y(t) in the frequency domain is E(s) = R(s) − Y (s) = R(s) −

E(s) =

G(s) R(s) R(s) = 1 + G(s) 1 + G(s)

R(s) 1 + G(s)

and the steady-state error erp = lim e(t) = lim sE(s) = lim s t→∞

s→0

erp = lim s s→0

5.1.8

s→0

R(s) 1 + G(s)

R(s) 1 + G(s)

(5.1.1)

Errors to standard input

To find the steady-state response of the systems the value of the steady-state error to the standard input step, ramp and parabola is calculated. • Unit step input The Laplace transform of the unit step is R(s) =

1 s

Substituting in the expression of the steady-state error 5.1.1 erp = ep = lim s s→0

R(s) 1 1 1 = lim = = 1 + G(s) s→0 1 + G(s) 1 + lim G(s) 1 + Kp s→0

(5.1.2) This error is called position error.

5.1. Introduction to errors

87

Fig. 5.6: • Unit ramp input The Laplace transform of the unit ramp is 1 s2 Substituting in the expression of steady-state error 5.1.1 R(s) =

erp = ev = lim s s→0

1 1 1 R(s) = lim = = 1 + G(s) s→0 s + sG(s) lim sG(s) Kv s→0

(5.1.3) This error is called velocity error.

Fig. 5.7: • Unit parabola input The Laplace transform of 1/2 unit parabola is   1 2 r(t) = t2 ⇒ R(s) = 3 R(s) = 3 s s Substituting in the expression of steady-state error 5.1.1 erp = ea = lim s s→0

R(s) 1 1 1 = lim 2 = = 2 2 s→0 1 + G(s) s + s G(s) lim s G(s) Ka s→0

(5.1.4) This error is called acceleration error.

88

5. Errors in steady state

Fig. 5.8:

5.1.9

The relationship between errors and system type

The relationship between the constant error and system type is represented in the table 5.9

Fig. 5.9: The relationship between errors and system type. and the relationship between errors and the types in the following table

ep Type 0 Type 1 Type 2 5.1.10

1 1+kp

0 0

ev ea ∞ ∞ 1 ∞ kv 0 k1a

Error input for polynomial in t

A lot of input signals can be approximated by a polynomial entry. For example, if the input is r(t) = 1 + 3t + 5t2

5.1. Introduction to errors

89

then the error in steady state is a composition of the errors of position, speed and acceleration 3 10 1 + + erp = ep + 3ev + 10ea = 1 + kp kv ka

90

5. Errors in steady state

6. ROOT LOCUS In control theory and stability theory, root locus analysis is a graphical method for examining how the roots of a system change with variation of a certain system parameter, commonly a gain within a feedback system. This is a technique used as a stability criterion in the field of control systems developed by Walter R. Evans which can determine stability of the system. The root locus plots the poles of the closed loop transfer function in the complex S plane as a function of a gain parameter. The dynamical behavior of a system is determined by the position of its poles. The place of the roots to determine the position of the poles in closed loop varies according constant K the system.

6.1

Concept of Root Locus Consider the closed loop system of the figure

Its transfer function is closed loop M (s) =

G(s) 1 + G(s)H(s)

Factoring the transfer function in open loop Q (s − zi ) G(s)H(s) = K Q (s − pi ) If we change the value of K, the poles of M (s) vary. (zeros remain unchanged). The Root Locus of M (s) is the place followed by its poles when K vary from zero to infinity.

92

6. Root Locus

The poles of F (s) are the roots of 1 + G(s)H(s) = 0, nad this is called characteristic equation. Q (s − zi ) 1+KQ =0 (s − pi )

(6.1.1)

where K is proportional to the static gain.

6.2

Magnitude and Angle Condition (or Module and Argument Criteria) From factoring in open chain Q (s − zi ) G(s)H(s) = K Q (s − pi ) The characteristic equation is Q (s − zi ) Q =0 1+K (s − pi )

namely, Q (s − zi ) Q K = −1 (s − pi ) As this equation is complex, it can be expressed as two real equations, first magnitudes (modules) and secondly the angles (arguments). The equation of the magnitudes Q |s − zi | Q K =1 |s − pi | It is usually written as Q |s − pi | product of the distances to the poles K= Q ⇒K≡ |s − zi | product of the distances to the zeros

(6.2.1)

This equation is the Magnitude Condition. Furthermore, the equation of the angles is X X arg(s − zi ) − arg(s − pi ) = (2q + 1)π that it is the Angle Condition and is usually written in the form X

arg(s − pi ) −

X

arg(s − zi ) = (2q + 1)π

(6.2.2)

The meaning of the magnitudes and the angles of these formulas can be in Fig. 6.1.

6.3. Tracing Rules of the Root Locus.

93

Fig. 6.1: Meaning of the magnitudes and angles.

6.3

Tracing Rules of the Root Locus.

• Number of branches If n = number of poles of G(s)H(s) m = number of zeroes of G(s)H(s) then Number of branches = max number of poles and number of zeroes of M(s) M (s) = max(m, n) This is because the characteristic equation Q (s − zi ) =0 1+KQ (s − pi ) This is because the characteristic equation Y Y (s − pi )+K (s − zi ) = 0 is a real polynomial coefficients whose degree is equal to the number of roots (Fundamental Theorem of Algebra). Each of these roots, to vary K, produces a branch. • Rule 2. Initial and Final points For K = 0, the characteristic equation is Y

(s − pi ) = 0

94

6. Root Locus

That is the Root Locus, starts in the closed-loop poles. On the other hand, if K = ∞, the characteristic equation can be written as Y 1 Y (s − pi )+ (s − zi ) = 0 K and as

1 K

=0 K=∞⇒

Y

(s − zi ) = 0

meaning that the place of roots ends in open loop zeros. The difference between poles and zeros corresponds to points in the infinite. • Rule 3. Real axis. If we apply the angle condition to a point of the real axis we get: X X arg(s − pi ) − arg(s − zi ) = (2q + 1)π The contribution of any of these terms if the real root is on the right is π.

Fig. 6.2: The contribution of any real root on the right is 0.

Fig. 6.3: Finally, the contribution of any pair of complex conjugate roots is es π −π = 0.

6.3. Tracing Rules of the Root Locus.

95

Fig. 6.4: As a consequence, to satisfy the angle condition we must have an odd number of roots on the right side of the point s.

Fig. 6.5: A point on the real axis will belong to the root locus if it has an odd number of real roots on the right. • Rule 4. Symmetry of the Root Locus The root locus is symmetric w.r.t. to the real axis. • Rule 5. Asymptotes

Fig. 6.6:

96

6. Root Locus

If we consider a point that is far, as is shown in Figure Fig. 6.6 point at infinity, then all angles of view of the argument X X arg(s − pi ) − arg(s − zi ) = (2q + 1)π They can be considered equal, and therefore X

arg(s − pi ) −

X

arg(s − zi ) = nθa − mθa = (2q + 1)π

namely, θa =

2q + 1 π n−m

(6.3.1)

• Rule 6. Centroid The cutting point of the asymptotes. This point is called ‘centroid’. P σ0 =

P p i − zi n−m

(6.3.2)

• Rule 7. Angles of departure and arrival Departure angles: If we take a point infinitely close to a pole in a departure branch.

Fig. 6.7:

X

arg(s − pi ) −

X

arg(s − zi ) =

X

αi + θ −

X

βi = (2q + 1)π

6.3. Tracing Rules of the Root Locus.

97

θ = (2q + 1)π −

X

αi +

X

βi

(6.3.3)

Arrival angles: The calculation is the same with a point in an arrival branch. • Rule 8. Points of dispersion and confluence They coincide with local maxima (dispersions) and minima (confluence) of K on the real axis dK(s) =0 ds or also X

(6.3.4)

X 1 1 = σ − pi σ − zi

• Rule 9. Intersection with the imaginary axis It is calculated by using the Routh method. • Rule 10. Value of K in a point of the Root Locus It is calculated by using the magnitude condition. Q |s − pi | K= Q |s − zi | • Rule 11. Sum of the roots If the characteristic eqn. is sn + an−1 sn−1 + an−2 sn−2 + ... + a1 s + a0 = 0 We have X

raices = − an−1

(6.3.5)

98

6. Root Locus

Fig. 6.8: Tracing rules for Root Locus.

7. CLASSICAL CONTROL TECHNIQUES 7.1

Problem approach

The objective of any control system is that given a certain system, which provides a certain amount of information based on a series of measures of the input / output signals, try to determine a certain control inputs feasible so that the system variable to be control continues to track as accurate as possible a certain reference signal, although the influence of the possible disturbances, measurement errors and variations in load system. A very general process control scheme Gp is displayed on the figure7.1. r Gf

u

+

y Gp

-

Gc

Fig. 7.1: Block diagram of a general control scheme In this scheme, the control signal is U (s) = Gf (s)R(s) − Gc (s)Y (s)

(7.1.1)

We can see a part formed by a controller Gc located in the feedback and a filter Gf at the input of system. This control scheme has as joint transfer function of the entire system the following expression Gp (s) Y (s) = Gf (s) . R(s) 1 + Gp (s)Gc (s)

(7.1.2)

In this expression it can be seen that the term filtering Gf directly affects the steady-state gain of the feedback system as well as to the poles and zeros of the same to be located in the main chain loop before the feedback. However, does not change the positions of the poles and zeros of the feedback loop , except that pole-zero cancellations occur between Gf and the closed loop transfer function.

100

7. Classical Control Techniques

Moreover the term Gc to be located in the feedback loop significantly affect the positions of the poles of the closed loop system, altering system dynamics . Simplifying a bit , we may think that while the design Gf adjust or improve steady-state response of the system, is sought Gc adjust the transient response. A particular case of this control structure is when Gf = Gc . In this case, the freedom in design is lost as it only has a transfer function . Thus, the control scheme is simplified as can seen in the diagram (a) of Figure ?? , being able to be expressed by the equivalent circuit diagram ( b ) of this Figure. r Gc

u

+

y Gp

(a) Gc

r

+ -

u

e Gc

y Gp (b)

Fig. 7.2: Controller in the principal loop This equivalent circuit (b) is the we could call classic control scheme, at least academically. Within this group we can include controllers PID control, which is what we will study in this chapter.

7.2

Design of classical controllers

Among the drivers that match the control scheme Classic can distinguish two basic approaches when address controller design Gc . These two approaches are fundamentally different ways of defining the function transfer controller . begin enumerate item controller with fixed structure. This group of techniques are characterized by using a transfer function defined a priori. This group includes the so-called PID control. In these controllers the transfer function of the controller is fixed and in them are present: proportional control to error ( P action ), proportional control to the integral of error ( action I) and proportional control to the derivative of the error ( action D ). The weight given to each of the control actions can be determined from theoretical form ( in the case of the known transfer function of the process) or the empirical form if such transfer function is not known. Clearly it can happen that

7.3. Specifications

101

one or more of the parameters of control override some actions if they do not will be necessary. item controller with variable structure. This second controllers group has not a predefined structure of its transfer function, but this result is obtained using the desired specifications for the system and the transfer function of the process to be controlled. This is the case of the design of controllers by direct synthesis technique. end enumerate

7.3

Specifications

The specifications express the characteristics that are desired to achieve in the system when a controller is introduced therein. Can be expressed in different ways using for this the frequency domain or the time domain. In our case we will refer to domain specifications time. Specifications regarding the domain response of time of a system are often defined with respect to the response time to a step input of a second system order. For this purpose, taking as reference system a system with theoretical transfer function: G(s) =

Kωn2 s2 + 2ζωn s + ωn2

(7.3.1)

Whose response to a step is p 1 y(t) = 1 − p e−ζωn t sin(ωn t 1 − ζ 2 + φ) 1 − ζ2

(7.3.2)

p 1 − ζ2 φ = arctan ζ

(7.3.3)

con

This system has its poles located at p1,2 = −ζωn ± jωn

p 1 − ζ2

(7.3.4)

where ζ is called damping ratio and damping ratio ωn is the undamped natural frequency index undamped natural frequency of the system. For a system with damping ratio 0 < ζ ≤ 1 (underdamped system) the system response curve has the form shown in Figure7.3 Among the parameters that can be used to specify the response can be included the following: • Peak time, Tp =

π ωd

where ωd is the damped frequency ωd = ωn

(7.3.5) p 1 − ζ 2.

102

7. Classical Control Techniques

y

To d·Mp

Mp

ep

Tp

t Ts

Fig. 7.3: Parameters of the time response • Overshoot, Mp . −π √

Mp = e • Settling time, Ts .

Ts '

√

ζ 1−ζ 2

π ζωn

=

d

(7.3.6)

(7.3.7)

Finally, the position error in steady-state Ep is defined as the difference between the desired reference value and the output value of the controlled variable of the system (for simplicity we are assuming that the feedback signal is unitary). It is relatively easy to obtain, from a certain specifications, the positions which must have the poles in the closed loop system that meet these specifications.

7.4

PID controllers

PID controllers are by far the most frequently used for industrial process control , where more than 95% of the control loops used are PID controllers. Their simplicity, versatility and ability to solve the basic problems that arise in the process favorable dynamics and modest performance requirements make them an essential tool in controlling processes. Over the past decades the PID controllers have survived various technological changes ranging from first controllers developed from pneumatic elements to the microprocessors developed going by electronic valves, the analog devices, transistors and integrated circuits . The incorporation of microprocessors has had a big impact because it has allowed the PID controllers to be enriched in other functions without losing any of their properties, so they have added possibility of automatic

7.4. PID controllers

103

adjustment of the parameters, possibilities of execution of logical rules or sequential automatisms. Therefore, and although there are more sophisticated control techniques, in the lowest level are still present or even is the majority of control processes. In the next few pages we will refer to them.

7.4.1

Basic structure of a PID controller

The academic version of this driver is characterized by a control law of the form:   Z 1 t de(t) u(t) = K e(t) + e(τ )dτ + Td (7.4.1) Ti 0 dt where u(t) is the control variable and e(t) is the error expressed as e(t) = r(t) − y(t) (Figure 7.2 b). The expression 7.4.1 shows the three types of control actions: the proportional action with weight K, the integral action with weight TKi and the derivative action with weight KTd . To set the behavior of the controller has three parameters: the proportional gain K, the integral time Ti and the derivative time Td . It can be seen that if the time integrated becomes infinite, the integral action disappears and if the time derivative becomes zero the derivative action disappears. The transfer function of the PID controller is as follows U (s) 1 K(1 + sTi + s2 Ti Td ) = K[1 + + Td s] = E(s) Ti s sTi

(7.4.2)

In this controller, the setting seeks to determine the positions of the two zeros of the transfer function and the static gain of the same, so that the desired design specifications are met the best possible in the system . Equation 7.4.2 is not physically realizable , so it is called theoretical PID controller. For a PID to be physically realizable there will be added to the derivative action a pole with limited influence. A real PID controller will have the form Gc (s) = K

1 +2 Ti Td sTi + s (1 + sTi )(1 + sTd )

(7.4.3)

The adjustment of the parameters of the controller can be done in two ways : 1. Empirically . For what it is adjusted experimental values ??to the desired response. This system may be too slow in many systems. Especially if the response time of these is great. To avoid this problem, the Ziegler-Nichols methods are used. These methods are based on some very simple measures observed in the response of the system and give some theoretical values ??of the parameters of the controller.

104

7. Classical Control Techniques

Values ??are taken as reference and then they require a finer adjustment. These methods and others derived from them are considerably used in industry, both for purely manual adjustment and industrial PID controller with auto-tuning systems. 2. Theoretically. To determine what are the analytical values ??of the regulator. The explicit knowledge of the transfer function of the process is required.

7.4.2

Ziegler-Nichols methods

The Ziegler -Nichols methods are two classic methods for the empirical adjustment of the parameters of a PID controller. They were presented by these authors in 1942. These methods are still widely used, either in its original form or with improved versions. Both methods are based on the determination of some characteristics of the process response, (time or frequency), to determine from these features and by means of simple mathematical relationships the controller parameters that givie a response rate of decay of a fourth between the values ??of the first and the second overshoot . Step response Tuning PID method Step response Tuning PID method This first method is based on the observation of the open loop response of the system to a step input. Noting this response, two parameters that are obtained are determined as follows : 1. The point of maximum slope is determined by the response curve to a step, and at that point the tangent is drawn to the response curve at that point. 2. The intersection of the tangent obtained above with the coordinate axes is determined and the distances a and L. These two parameters correspond to the theoretical response of a mathematical model of the form G(s) =

a −sL e sL

(7.4.4)

corresponding to a time delay integrator. This system can be characterized by two parameters a and L as can be seen in the figure 7.4. 3. Once you have determined the parameters a and L in response, Ziegler and Nichols proposed as the PID controller parameters indicated in Table 7.1 obtained directly as a function of the parameters a and L measured on the system response.

7.4. PID controllers

105

y

L

t

a

Fig. 7.4: Response to step input of G(s) = Controller P PI PID

K 1/a 0.9/a 1.2/a

Ti

Td

3L 2L

L/2

a −sL sL e

Tab. 7.1: Parameter values given by Ziegler-Nichols method for the step response This controller is designed to give a decay rate of a quarter ( d = 0.25) between the magnitudes of the first and second overshoot, so it usually has a high overshoot. It has the advantage that from these values ??is easy to make a finer adjustment to fit the response you want without the need of a long process of trial and error . subsubsection Method of the frequency response index method to adjust PID’s using the frequency response This second method developed by Ziegler and Nichols is based on determining the point on the Nyquist curve for the transfer function G(s), in which that curve intersects the negative real axis ( point of relative stability ) . Once obtained this point is characterized by two parameters named : Gain critical Ku and critical period Tu . This method is based on the observation that many systems become unstable when proportional control has high gain. The method of determination of the parameters of the PID controller has the following steps: 1. the feedback loop is closed by a proportional regulator (Figure 7.5a). It necessary to increase the gain of the proportional controller until you get the critical gain Ku that correspond to an oscillatory response as indicated in the Figure 7.5. For this gain Ku , is determined the Tu critical period of the oscillatory output signal. 2. As in the first method, are determined the values of the parameters of the con-

106

7. Classical Control Techniques

r

+ -

e

u K

y G(s)

a)

y

K
b)

K=Ku, Tu

c)

K>Ku

d)

t y

t

y

t

Fig. 7.5: Obtaining the critical gain and critical period for the second method of Ziegler-Nichols troller from from Ku and Tu . The values proposed by Ziegler and Nichols is shown in Table 7.2. As in the previous method the design criterion is to regulators decay rate of a quarter. These controllers give good response to changes in system load, but produce frequently a high overshoot, as already discussed above. In any case, it’s easy to make a further adjustment of the regulatory on the basis of these parameters.

7.5. Analytical methods for PID controllers design

Controlador P PI PID

K 0.5Ku 0.4Ku 0.6Ku

107

Ti

Td

0.8Tu 0.5Tu

0.125Tu

Tab. 7.2: Parameter values given by Ziegler-Nichols frequency method

7.5

Analytical methods for PID controllers design

These methods seek to determine the most appropriate form analytical values ??of the parameters of a PID controller so that certain specifications are met. See two analytical techniques for calculating these parameters: the first of these techniques (pole assignment) is applicable to systems of second order, and the second (design based on the dominant poles) is applicable to complex systems. Both techniques require prior knowledge of the transfer function of the process to be controlled.

7.5.1

Pole placement method for adjusting PID

Pole placement method for adjusting PID The objective of the technique of assignment of poles index pole assignment is to find the parameters of a controller such that the transfer function in closed loop has its poles in the predetermined positions. Closed loop poles are determined so that a certain specifications are obtained in the system response. In general terms a regulator of type P to have a single parameter alone would allow us to assign only a pole, a PI controller or a PD having two parameters would allow us to assign two poles and a PID controller has three parameters that would allow us to assign three poles. The use of either type of controller depends on the transfer function of the system. Suppose first a process whose transfer function is first order, with the following parameterization Kp Gp (s) = (7.5.1) 1 + sT and suppose we want to control the system through a PI controller type, which responds to the following expression   1 Gc (s) = K 1 + (7.5.2) 1 + sTi the closed loop system has the following transfer function G(s) =

Gp (s)Gc (s) 1 + Gp (s)Gc (s)

(7.5.3)

108

7. Classical Control Techniques

The characteristic equation of the closed loop system defined by the equation 7.5.3 is the following 1 + Gp (s)Gc (s) = 0 (7.5.4) and substituting in 7.5.4 Gp (s) y Gc (s) by the expessions 7.5.1 y 7.5.2 obtain that   Kp 1 1+K 1+ =0 sTi (1 + sT )   Kp K 1 + Kp K + =0 s2 + s T T Ti

(7.5.5)

Suppose we want the response of closed loop system is that of a second order system with a damping rate ζ and frequency ωn . That is, you want p to place the poles in the closed chain system at positions p1,2 = −ζωn ± ωn 1 − ζ 2 . So the closed loop characteristic equation will have the expression s2 + 2ζωn s + ωn2 = 0

(7.5.6)

Igualando las expresiones 7.5.5 y 7.5.6 obtenemos que 1 + Kp K T K K p ωn2 = T Ti

2ζωn =

(7.5.7) (7.5.8)

y despejando en 7.5.7 y 7.5.8 tenemos que 2ζωn T − 1 Kp 2ζωn T − 1 Ti = ωn2 T

K=

(7.5.9) (7.5.10)

This technique requires the introduction of a controller which have two possible adjustment parameters, a PI controller, when the open-chain system is first order since it is necessary to fix arbitrarily two pole positions in a closed loop. Similarly, when the system is second order, it is usually convenient to use a PID controller that has three degrees of freedom to set arbitrarily the three poles that will have the closed-loop system. Let the transfer function of the second order system of the form Kp (7.5.11) (1 + sT1 )(1 + sT2 ) and suppose we want to control the system via a PID controller, which responds to the following expression Gp (s) =

Gc (s) = K

(1 + sTi + s2 Ti Td ) sTi

(7.5.12)

7.5. Analytical methods for PID controllers design

109

so that the characteristic equation of the closed loop system 1 + Gp (s)Gc (s) = 0, will be     Kp KTd Kp K Kp K 1 1 1 3 2 + + +s + + = 0 (7.5.13) s +s Ti T2 T1 T2 T1 T2 T1 T2 T1 T2 If you want the system to have the following closed-loop characteristic equation, (s + α)(s2 + 2ζωn s + ωn2 ) = 0

(7.5.14)

equating coefficients of the polynomial expressions 7.5.13 y 7.5.14 we have the following expressions Kp K ωn2 α= T T 1 2  Kp K 1 ωn2 + 2ζωn α= + T1 T2 T1 T2   Kp KTd 1 1 2ζωn + α= + + Ti T2 T1 T2

(7.5.15) (7.5.16) (7.5.17)

and solving this system of equations we obtain T1 T2 (2ζωn + α) + 1 Kp Kp K Ti = T1 T2 ωn2 α T1 T2 ωn3 α(ωn + 2ζα) Td = T1 T2 (2ζωn + α) − 1 + T1 ωn2 α(ωn2 α + 1) K=

(7.5.18) (7.5.19) (7.5.20)

Example 7.5.1: Suppose you want to control a system whose transfer function 1 so that has a closed loop damping rate ζ = 0.5 and a frequency is GP (s) = s(s+2) ωn = 3.14 rad s . To control this system we need is a PD type controller with transfer function of the form Gc (s) = K(1 + Td s) (7.5.21) because the system is second order but has a pole at the origin in the transfer function, so that the system does not require adding an integral effect on the controller to override the position error in steady-state. By not including regulatory effect on the whole system will have a closed loop characteristic equation of order 2. The closed loop system has a characteristic equation of the form

110

7. Classical Control Techniques

1 + K(1 + Td s)

1 =0 s(s + 2)

(7.5.22)

and you want to have the characteristic equation ?? closed loop. To do this we equate 7.5.6 to 7.5.22, s2 + s(2 + KTd ) + K = s2 + 2ζωn s + ωn2 = 0

(7.5.23)

and equating coefficients of two polynomials have to K=ωn2 2ζωn − 2 Td = ωn2

(7.5.24) (7.5.25)

and substituting the values ??of the specifications we have indicated that K = 9.86 and Td = 0.115. So that the the resulting controller has the expression: Gc (s) = 9.86(1 + 0.115s)

(7.5.26)

2.5 y(t) 2 (a) 1.5

1

(b)

0.5

0

0

0.5

1

1.5

2

2.5

3

3.5

4

4.5

t

5

Fig. 7.6: Response system of Example 7.5.1: (a) open loop and (b) closed loop with the designed controller

7.5.2

Method based on the dominant poles for adjusting PID

The pole assignment technique is to locate the positions of all the poles in the closed loop system. This complicates considerably the method for the case of

7.5. Analytical methods for PID controllers design

111

complex systems, the problem may even have no feasible solution. Thus, if you want to control with a PID controller a system through this technique, the system should be restricted to first and second order to assures the existence of a solution. If the system is complex, or we are looking for an approximation of the model to a first or second order, this is the technique commonly used design based on the dominant poles . This technique assigns only a few poles of the system. This allows to design simple controllers for controlling complex systems, although the controller behaviour do not assure the given specifications since the design ignores the effect of some poles. However, if the choice of positions is appropriate, the behaviour will be very close to what is desired . This technique also requires that the transfer function of the system to be known, and the idea of assigning only a few poles of the closed loop system has been widely used in the early design of PID control, especially when the power of computers did not allow other possibilities. Even today , this technique has considerable interest since many industrial controllers implement a fixed PID structure and are used to be applied to processes of variable complexity. This allocation based on the dominant poles can be performed in two ways: either by determining the values of the controller parameters algebraically which place the dominant poles in the desired positions or by using the Root Locus. Both techniques allow us to obtain analytically the parameters of the PID controllers.

Method based on algebraic determination for adjusting PID The problem is similar to the previous section. If you enter P controller, we need to put a pole in order to determine the controller parameter. If you want to introduce a PI or PD controller is necessary to assign two poles to calculate the parameters controller and a PID controller if it is necessary to set three poles. Suppose, for illustrating the operation of this technique, it is desired to design a PI controller for controlling a system. This controller has two parameters by which to determine the value of them need to locate the position of both poles of the closed loop system. Assume the controller equation is of the form  Gc (s) = K 1 +

1 1 + sTi

 (7.5.27)

the closed loop system has the following transfer function G(s) =

Gp (s)Gc (s) 1 + Gp (s)Gc (s)

(7.5.28)

112

7. Classical Control Techniques

so that the characteristic equation of the closed loop system is defined by 

1 1 + Gp (s)K 1 + 1 + sTi

 =0

(7.5.29)

We need to set two pole positions p1,2 , so we can assume that these are complex conjugates so that they have a certain damping ζ and a certain frequency ωn , ie p p1,2 = ωn (−ζ ± j 1 − ζ 2 ) (7.5.30) In the event that the driver you need is type of PD are in a similar situation, because we have to determine two parameters of the controller, whose transfer function is Gc (s) = K(1 + sTd )

(7.5.31)

If what you want is to design a PID controller to control such a system. This controller has three parameters by which to determine the value of them need to locate the position of three poles of the closed loop system. Suppose that the controller equation is of the form  Gc (s) = K 1 +

1 + sTd 1 + sTi

 (7.5.32)

the closed loop system has the following transfer function G(s) =

Gp (s)Gc (s) 1 + Gp (s)Gc (s)

(7.5.33)

so that the characteristic equation of the closed loop system is defined by  1 + Gp (s)K 1 +

1 + sTd 1 + sTi

 =0

(7.5.34)

We need to set three pole positions p1,2,3 , so we can assume that two of these are complex conjugates so that they have a certain damping ζ and a certain frequency ωn , and the other is a pole located at −α ie p p1,2 =ωn (−ζ ± j 1 − ζ 2 ) (7.5.35) p3 =−α

(7.5.36)

The determination of the values of the parameters of the PID controller would be analogous to that seen for controllers with two parameters.

7.5. Analytical methods for PID controllers design

113

Fig. 7.7: System of the example 7.5.2 Example 7.5.2: Consider the system shown in Figure 7.7, whose transfer function is 1 Gp (s) = (7.5.37) s(s + 2)(s + 8) and suppose we want to introduce a controller in the system so that the system has a damping rate ζ = 0.5 and ωn = 3.14. If we fix only the dominant poles, these specifications would have to said poles should be p1,2 = ωn (−ζ ± j

p 1 − ζ 2 ) = −1.57 ± j2.72

(7.5.38)

Since the system is of type 1, which includes an integrator in open chain, introduce a PD type controller in the K(s + a) system. First obtain the characteristic equation of the closed loop system 1 + K(s + a)

1 = s3 + 10s2 + (16 + K)s + Ka = 0 (7.5.39) s(s + 2)(s + 8)

Moreover, we fix the positions of two of the poles at positions 7.5.38. As the thirdorder system is the third pole will be in a position (s + b) so that the closed-loop characteristic equation will have the following appearance (s+b)(s2 +2ζωn s+ωn2 ) = s3 +(2ζωn +b)s2 +(2ζωn b+ωn2 )s+ωn2 b = 0 (7.5.40) By matching the characteristics equations 7.5.39,7.5.40 we obtain the following equations 10=2ζωn + b 16 + K=2ζωn b + ωn2

(7.5.41)

Ka=ωn2 b and solving for in 7.5.41 we have: b=10 − 2ζωn = 6.86 K=2ζωn b + ωn2 − 16 = 15.4 ω2 b a= n = 4.392 K

(7.5.42)

114

7. Classical Control Techniques

so that the PID controller that puts the dominant poles at desired positions in closed chain is Gc (s) = 15.4(s + 4.392)

(7.5.43)

Note that the position of the third pole in b = 6.86 is far away from the dominant poles, so their influence is very small. 1.4

y(t)

1.2 (b) 1 0.8 0.6 0.4 (a) 0.2 0

0

0.5

1

1.5

2

2.5

3

3.5

4

4.5 t (s)

5

Fig. 7.8: Response system of Example 7.5.2: (a) open loop and (b) closed loop with the controller designed

Method based on the Root Locus for adjusting PID Method The method of Root Locus , introduced by Evans in 1954, is a graphical method to determine the positions of the poles in closed loop from the positions of the poles and zeros of the system in open open as any of the parameters (typically the gain) varies from zero to infinity . In practice, the drawing of the Root Locus of some system can tell us if it will be sufficient to vary gain system to meet the specifications that have given us. That is, in case that the location of the roots pass through the desired locations for the poles closed chain with enough to adjust to the proper value the system gain to meet the required transient response specifications. If this not, it will be required to the system controller to alter the place of the required roots in the right way. If we remember the effect that had the addition of poles and zeros system is relatively simple to determine what type of structure must have the driver that we introduce in the system.

7.5. Analytical methods for PID controllers design

115

• Effect of addition of poles. The addition of a pole to the transfer function of an open-loop system has the effect of adding a branch more to the Root Locus and branches as a result tend to approach or even enter the positive real half-plane. This tends to decrease the relative stability of the system and that the settling time of the response is longer. jω

x

(a)

0

jω

σ

x

x

0

(b)

jω

σ

x x

x

0

σ

(c)

Fig. 7.9: Effect of the addition of a pole on the Root Locus This effect can be seen in Figure 7.9, where poles are added as the number of branches and these are increasingly closer to the half-plane with positive real part. Coming in case (c), the introduction of the branches on the positive real part, for high values ??of gain and to become the system unstable. • Effect of the addition of zeros. The addition of a zero to the transfer function in open loop decreases the number of branches that go to infinity so increases the relative stability of the system to move the branches to the negative real part of the complex plane. The system as a faster response by introducing a derivative effect that anticipates the evolution of the system. This effect can be seen in Figure 7.10, where adding a zero to a system with three poles as shown in 7.9c only two branches go to infinity and none of them comes in part positive real, so that the system becomes unstable for any value of the gain. Example 7.5.3: For the system of Figure 7.7 can be seen that if a proportional controller is introduced and the value of gain K is varied in 1 + KGp (s) = 0

(7.5.44)

the poles of the characteristic equation of the system move through the branches of the Root Locus, which for this system would have an aspect as indicated in fig-

116

7. Classical Control Techniques

jω

ο

x

x

x

0

(a)

jω

σ

x οx

x

0

jω

xο x

x

σ

(b)

0

σ

(c)

Fig. 7.10: Effect of the addition of a zero on the Root Locus

ure 7.11. Looking at the positions of the dominant poles of the system PD can be seen that the Root Locus does not pass through them, so it is not enough with a proportional controller. If you want the Root Locus pass through these points and that these become the dominant poles of the closed chain system, we just need to introduce a PD type controller. This regulator introduce a zero and if the position of it is adequately chosen, the Root Locus can pass trough the location of the desired positions for the roots.

Pd ωn x -8

x -2

θ

x 0

P'd

Fig. 7.11: Root locus for the example system 7.5.3

7.5. Analytical methods for PID controllers design

117

In this case, it is necessary to determine the positionof the zero and the value of the gain necessary for the poles to be located exactly at the desired points. Recalling that the Root Locus was the locus of the points that verify the characteristic equation of the system in closed loop 1+G(s) = 0, this equation led to the module and angle conditions. |G(s)|=1

(7.5.45)

∠G(s)=(2q + 1)π

q = 0, ±1, ±2, . . . (K > 0)

(7.5.46)

if the transfer function G(s) has the form G(s) = K

(s + z1 )(s + z2 ) . . . (s + zm ) (s + p1 )(s + p2 ) . . . (s + pn )

where n ≥ m the expressions 7.5.45 and 7.5.46 will be as follows Qm |s + zi | Q =1 |G(s)| = |K| ni=1 |s + pj | j=1 or what is the same Qn |K| =

|s + pj | Qj=1 m i=1 |s + zi |

Qn

j=1 = Qm

distances from pj to Pd

i=1 distances

desde zi to Pd

(7.5.47)

(7.5.48)

(7.5.49)

and

∠G(s) = ∠K +

m X

∠(s + zi ) −

i=1

n X

∠(s + pj ) = (2q + 1)π

(7.5.50)

j=1

or what is the same ∠K +

m n X −−→ X −−→ ∠zi Pd − ∠pj Pd = (2q + 1)π i=1

(7.5.51)

j=1

In our case, given the positions of the poles of the process and controller, we have to if we want the Root Locus pass through the positions p p1,2 = ωn (−ζ ± j 1 − ζ 2 ) = −1.57 ± j2.72 (7.5.52) must meet the magnitude 7.5.48 and angle conditions 7.5.50. The arguments for the poles whose positions are known (fig 7.12): β1 =120o 2.72 β2 =arctan = 81.01o 2 − 1.57 2.72 β3 =arctan = 22.93o 8 − 1.57

(7.5.53) (7.5.54) (7.5.55)

118

7. Classical Control Techniques

as the angle condition must be met 180o = α1 − (β1 + β2 + β3 )

(7.5.56)

and solving for α1 = 43.93◦ . Knowing this angle, we can determine the position where it should be located the zero to form that angle with the point where the dominant pole is located −1.57 + 2.72j. The result is that the zero is located at a = 4.39. To determine the gain K we apply the magnitude condition K=

n1 n2 n3 m1

(7.5.57)

and substituting the values of n1 , n2 , n3 , we obtain that m1 = 15.4K. So that the PD controller that puts the dominant poles at desired positions in closed chain is Gc (s) = 15.4(s + 4.39)

n3 x -8

(7.5.58)

Pd j2.72 ωn=n1 m1

β3 -a

α1

n2 x -2

β2 θ

x

β1 0

Fig. 7.12: Application of the magnitude and angle conditions for Example 7.5.3 As can be seen the values obtained by both methods are similar (Fig.7.13).

7.5. Analytical methods for PID controllers design

119

1.4

y(t)

1.2 (b) 1 0.8 0.6 0.4 (a) 0.2 0

0

0.5

1

1.5

2

2.5

3

3.5

4

4.5 t (s)

5

Fig. 7.13: Response system of Example 7.5.3: (a) open loop and (b) closed loop with the controller designed Example 7.5.4: Find a controller for an continuous system with transfer function 1 , with temporal specifications damping coefficient ζ = 0.5 and G(s) = s(s+2) settling time ts = 2 sec. r + -

e

Gc(s)

u

1 s(s+2)

y

Fig. 7.14: Example 7.5.5 System The position of the dominant pole of the system is obtained from the temporal specifications ζ = cos(θ) → θ = 60o ts ≈

π → σ = 1.57 σ

p1,2 = −1.57 ± j2.72 As shown in Fig. 7.15, to place the dominant poles in the desired positions, it is necessary to add a regulator P D to the the system with transfer function Gc (s) =

120

7. Classical Control Techniques

n3 x -b

n2

Pd j2.72 ωn=n1

β1 β2 θ β3 x x 0 -a=-2

m1=n2 α1=β2

Fig. 7.15: Distribution of poles and zeros of example 7.5.4 K s+a s+b . To set the pole position in −b, use the method of cancellation of the zero in −a with a regulator with a pole in −2. Using the angle condition we have X

α−

X

β = (2q + 1)π

α1 = 81o , β1 = 120o , β2 = 81o

α1 − β1 − β2 − β3 = (2q + 1)π β3 = −(2q + 1)π − 120o = 60o Thus the position of pole of the controller b is −3.14. For the determination of the gain we use the magnitude condition K=

n1 n2 n3 3.14 × 2.75 × 3.14 = = 9.86 m1 2.75

Thus the continuous controller that meets specifications is shaped as Gc (s) = 9.86

s+2 s + 3.14

7.5. Analytical methods for PID controllers design

121

Fig. 7.16: Example 7.5.5 System Example 7.5.5: Find a controller for an continuous system with transfer function 1 G(s) = (s+1)(s+2) , with temporal specifications: overshoot Mp = 5% and settling time ts = 1 sec. The position of the dominant pole of the system is obtained from the temporal specifications −π Mp = e tan θ = 0.05 → θ = 45o ts ≈

π = 1 → σ = 3.14 σ

ωd = σ tan θ = 3.14 p1,2 = −3.14 ± 3.14 ± j2.72 The Root Locus doesn’t contains the desired poles positions, as shown in Fig 7.17.

Fig. 7.17: The Root Locus doesn’t contains the desired poles positions.

122

7. Classical Control Techniques

It is necessary to add a zero on the left to displace the Root Locus to the left, as shown in Fig 7.18.

Fig. 7.18: It is necessary to add a zero to displace the Root Locus to the left. Using the angle condition, as shown in Fig 7.19, we have

Fig. 7.19: Use of the angle condition.

X i

α−b = atan

αi −

X

βj = (2q + 1)180o

j

3.14 3.14 3.14 , β−1 = 180o − atan , β−2 = 180o − atan b − 3.14 1.14 2.14 α−b − β−1 − β−2 = (2q + 1)180o

7.5. Analytical methods for PID controllers design

α−b = −(2q + 1)180o − 125.8o = 54.2o Thus the position of pole of the controller b is 5.4. For the determination of the gain we use the magnitude condition √ √ d−1 d−2 3.142 + 1.142 × 3.142 + 2.142 p K= = 3.28 = d−b 3.142 + (5.4 − 314)2 Thus the continuous controller that meets specifications is shaped as Gc (s) = 3.28(s + 5.4)

123

124

7. Classical Control Techniques

Bibliografia [1] R. Aracil y A. Jimenez, Sistemas discretos de control, Servicio de publicaciones ETSIIM-UPM, 1980. [2] K. Aström and T. Hägglund, PID Controllers:Theory, Design, and Tuning, ISA - The Instrumentation, Systems, and Automation Society, 1994. [3] G.F. Franklin, J.D. Powell and M. Workman, Digital Control of Dynamics Systems, Addison Wesley, 1997. [4] R. Isermann, Digital Control Systems: Volume I Fundamentals, Deterministic Control , Springer-Verlag, 1989. [5] K. Ogata, Discrete-Time Control Systems, Prentice Hall International Edition, 1995. [6] K. Ogata, Ingeniería de Control Moderna, Prentice Hall, 1998. [7] E. Andres Puente, Regulación Automática, Servicio de publicaciones ETSIIMUPM, 1993.

8. FREQUENCY MODELS. BODE DIAGRAMS When a sinusoidal signal of unit amplitude and variable frequency as u(t) = 1 sin ωt is introduced into a stable linear system, the system responds by outputting a signal whose Laplace transform is given by w Y (s) = G(s) 2 . (8.0.1) s + ω2 If this signal is expanded in partial fractions, it appears a set of terms corresponding to the poles of the transfer function G(s) and the terms corresponding to the sinusoidal input. Namely c c∗ + . (8.0.2) s + jω s − jω If all poles of the system give a stable response, your response will decrease over time so the the steady-state response will retain only the last two terms of expansion 8.0.2. So the steady-state response in the time domain will take the form Y (s) = . . . +

y(t) = 2|c| sin(ωt + φ) = A sin(ωt + φ)

(8.0.3)

A = |G(jω)|

(8.0.4)

where

Im[G(jω)] = ∠G(jω) (8.0.5) Re[G(jω)] The frequency response of a system is usually represented by two curves, of which, the first is the magnitude |G(jω)| and in the other is the angle ∠G(jω) for each value of ω. The the stability of the closed loop system can be determined using the frequency response open loop system. Besides these response curves may be obtained experimentally without have knowledge of the mathematical model of the system, no more than excite the system with a variable frequency sinusoidal signal and measuring the amplitude and the phase of the sinusoidal output signal. φ = arctan

126

8. FreQuency models. Bode diagrams

When in the diagrams represent the amplitude in decibels 20 log|G(jω)| and the frequency represent logarithmically, the resulting diagram is called Bode diagram. This diagram is formed by two response curves: one it is the magnitude response of the system for different frequencies of the input signal and another one in which the phase delay of the response is represented for different frequencies of the input signal. This way of characterizing a system has some properties that make it very interesting, the first is the ease with which they can be obtained experimentally these curves for a considerable number of (electrical, electronic, mechanical, etc.) systems. The second is interpretative in that it allows us to see or treat any physical system as a filter it were (normally low pass) that tend to amplify or dampen sinusoidal signals of certain frequencies to be presented at the entrance system (remember that any periodic signal could be decomposed into series of Fourier). Example 8.0.1: Consider a system whose transfer function is 5 (8.0.6) +s+5 In order to represent its frequency response, we use the command bf bode () MATLAB, as shown below G(s) =

s2

num=[5]; den=[1 1 5]; sys=tf(num,den); bode(sys); and the result can be seen in Fig. 8.1. In the reasoning process that was followed to determine the frequency response of a system we started the mathematical model in s = σ+jω, and has been replaced by jω considering that the value is studied steady state. It is clear therefore that this model G(jω) mathematically derived from the temporal model.

8.1. Asymptotic Bode diagram

127

Bode Diagrams From: U(1) 10 5 0 −5

Phase (deg); Magnitude (dB)

−10 −15 −20 −25 −30 0

To: Y(1)

−50

−100

−150

−200 −1 10

0

1

10

10

Frequency (rad/sec)

Fig. 8.1: Bode plot of a second order system

8.1

Asymptotic Bode diagram

The Asymptotic Bode diagram can be obtained as sum of the Asymptotic Bode diagram of each of the factors of the transfer function. Consider a system with transfer function factored as Q (s − zi ) G(s) = K Q (s − pi ) The frequency response is obtained by replacing s with jω. Q (jω − zi ) G(jω) = K Q ⇒ (jω − pi ) The magnitude and the angle are Q |jω − zi | |G(jω)| = |k| Q |jω − pi | arg (G(jω)) = arg(k) +

X

arg(jω − zi ) −

By taking logarithms products become sums

X

arg(jω − pi )

128

8. FreQuency models. Bode diagrams

 ⇒

P P log |G(jω)| = log |k| + Plog |jω − zi | − P log |jω − pi | arg (G(jω)) = arg(k) + arg(jω − zi ) − arg(jω − pi )

(8.1.1)

For this reason, the modules are taken in decibels (dB), which are 20 log10 of the given quantity. Also on the horizontal axis, the frequency is taken on a logarithmic scale. The asymptotic Bode diagram is an approximate diagram of the transfer function from the sum of the representation of each term of the expression 8.1.1.

8.1.1

Asymptotic Bode diagram. Constants K

The constant K is a real number. Therefore, the magnitude is 20log10 |K| dB and the phase is 0 if K > 0 and −1800 if K < 0.  |K| = 20 log |K| dB | | =  00 k > 0 K= ϕ = 1800 k < 0

Fig. 8.2: Asymptotic Bode Plot of Gain Term

8.1.2

Asymptotic Bode diagram. Poles / zeros in the origin

1 sn

The poles at the origin s1n of order n have as module −20n log(ω) as shown in the semilog plot of magnitude as a straight line with a negative slope of −20n and the phase as a horizontal line at the position −900 n. 1 = (jω)n

1 = (jω)n



| | = ω1n = −20n log(ω) dB ϕ = −n 900

( 1 1 −n = −20n log(ω) dB 20 log (jω) n = 20 log n = 20 log ω ω ϕ = −n900

8.1. Asymptotic Bode diagram

129

Fig. 8.3: Asymptotic Bode plot. Poles in origin

1 sn

The magnitude diagram is a straight line that intersects the horizontal line of 0 dB in ω = 1, because −20n log ω = 0 =⇒ ω = 1 . For the zeros at the origin of order n module is 20n log(ω) as it is shown in the semi-log plot as a straight line with slope 20n positive argument and as a horizontal line in the 900 n position.

8.1.3

Asymptotic Bode plot. Negative real Polo/ Zeros

1 (1+T s)n

The negative real poles /zeros can be written as (1+T1 s)n , with n > 0 for the poles y n < 0 for the zeros. For negative real poles, the amplitude is a horizontal line 0 dB for low frequencies and a line of slope −20N dB/dec to high frequencies. The phase is a horizontal line 00 for low frequencies and −900 n for high frequencies. 1 = (1 + T jω)n   | | = 0dB   ω → 0 1 = ϕ = 00  = 1 | | = (T ω) n = −20n log(T ω) dB  1  ω → ∞ (T jω)n = ϕ = −n900 For negative real zeros, the amplitude is a horizontal line 0 dB/dec for low frequencies and a line of slope 20n dB/dec for high frequencies. The phase is a horizontal line of 00 for low frequencies and for high frequencies 0 90 n.

130

8. FreQuency models. Bode diagrams

(1  + T jω)n =  | | = 0dB   ω → 0 1 = ϕ = 00 = | | = (T ω)n = 20n log(T ω) dB   ω → ∞(T jω)n = ϕ = +n900 Cutoff Frequency Cut of the two asymptotes 0 db = −20 log(ωc T ) log(ωc T ) = 0 ⇒ ωc T = 1 ωc = T1

Fig. 8.4: Asymptotic Bode plot. Negative real Poles/Zeros

8.1.4

1 (1+T s)n

Asymptotic Bode plot. Negative complex Pole/Zeros

The complex negative Poles/Zeros can be represented as ωn 2 = s2 + 2ζωn s + ωn 2 1+

1 2ζ ωn s

+



s ωn

2  n

with n > 0 for the poles and n < 0 for the zeros. If the equation s2 + 2ζωn s + ωn2 = 0 has complex roots with negative real part, then 0 < ζ < 1. The asymptotic approximation is

8.1. Asymptotic Bode diagram

1   2 n 1+ ω2ζ jω+ ωjω n

131

=

n

  | | = 0dB    ω → 0 1 = ϕ = 00 (  −2n =  −2n  | | = ωωn = −40n log( ωωn ) dB jω  = ω → ∞   ωn ϕ = −2n900 = −n1800 Cutoff frequency Cut of the two asymptotes 0 db = −40n log( ωωn ) dB log( ωωn ) = 0 ⇒ ωωn = 1 ωc = ωn

Fig. 8.5: Asymptotic Bode plot. Negative complex Poles/Zeros. Resonance Frequency ωr ωr = ωn

p 1 − 2ξ 2

(0 < ξ < 0.707)

The Resonance Frequency ωr exists when: 1 − 2ξ 2 > 0

⇒

1 ξ < √ = 0.707 2

Resonance peak according to ζ: |G(jω)|max = |G(jωr )| = Mr =

1 p (0 < ξ < 0, 707) 2ξ 1 − ξ 2

132

8. FreQuency models. Bode diagrams

Fig. 8.6: Bode plot. Negative complex Poles/Zeros

8.1.5

Asymptotic Bode plot. Positive real Poles/Zeros.

1 Positive real Poles/Zeros can be expressed as (−1+T , with n > 0 for the s)n poles and n < 0 for the zeros. For positive real poles, the amplitude is a horizontal line 0 dB for low frequencies and a line of slope −20n dB / dec to high frequencies.

The phase is a horizontal line of −1800 n for low frequencies and −900 n for high frequencies.

1 n (−1+T  jω)

=

=

 1  ω → ∞ (T jω) n



| | = 0dB ϕ = −1800 n  1 | | = (T ω) n = −20n log(T ω) dB = 0 ϕ = −90 n

 n  ω → 0 (−1) =

For positive real zeros, the amplitude is a horizontal line 0 dB for low frequencies and a line of slope 20n dB/dec to high frequencies. The phase is a horizontal line of −1800 n for low frequencies and 900 n for high

8.1. Asymptotic Bode diagram

133

frequencies. n (−1  + T jω) = | | = 0dB  n  ω → 0 (−1) = ϕ = 1800 n  = | | = (T ω)n = 20n log(T ω) dB   ω → ∞(T jω)n = ϕ = +900 n

Fig. 8.7: Asymptotic Bode diagram of 1/s − 1. (Real positive poles).

Fig. 8.8: Asymptotic Bode diagram ofs − 1. (Real positive zero). In short, in the Bode plots of positive real poles and zeros, the magnitude is the same as in the negatives and the phase is reversed.

134

8. FreQuency models. Bode diagrams

8.1.6

Asymptotic Bode plot. Positive complex Poles/Zeros.

Negative complex Poles/Zeros can be expressed as ωn 2 = s2 + 2ζωn s + ωn 2 1+

1 2ζ ωn s

+



s ωn

2  n

with n > 0 for the poles and n < 0 for the zeros. The asymptotic approximation for the poles is 1   2 n 1+ ω2ζ jω+ ωjω n

=

n

   ω → 0 1 = | | = 0dB   ϕ = 00 (  −2n =  −2n ω  = −40n log( ωωn ) dB | | = jω  ωn = ω → ∞   ωn ϕ = 2n900 = n1800

Fig. 8.9: Bode plot with complex poles with: a) negative real part, b) positive real part.

9. NYQUIST METHOD Nyquist stability criterion The method determines Nyquist stability through a graphic-analytical method, the stability of a closed-loop system based on the frequency response G(jω)H(jω) in open loop. In addition to determining the absolute stability of a system, we also used to determine the degree of stability-instability of a system of regulation. Provides information about how to improve the system performance in both steady-state and transient. It is based on the argument principle Cauchy for functions of a complex variable. Concepts and principles of the complex variable theory If we choose a closed path ρ in the complex plane, being F (s) analytic at all points of the road. Then the curve image of ρ in the plane F (s) will be also closed (function analytical means that they exist all its derivatives).

ρ F(s) F(ρ)

Fig. 9.1: Closed path ρ and its image F (s)

Theorem 9.0.1: Cauchy argument principle Let F (s) a complex function of a complex variable, and ρ a closed path in the complex plane which does not pass through singular points of F (s). Then F

136

9. Nyquist Method

(ρ) gives a number of turns around the origin equal to the difference between the number of zeros and poles of F (s) within ρ (Fig. 9.1). That is, if Z =number of zeros of F (s) within ρ P = number of poles of F (s) within ρ N = number of turns of F (ρ) around the origin. Then N = Z − P. Introduction to Nyquist method The Nyquist method allows us to study the stability of a feedback system from open loop frequency response. If we consider the system of Figure 9.2, the transfer function in closed loop will be: G(s) M (s) = 1 + G(s)H(s)

r(t) +

G(s)

y(t)

H(s)

Fig. 9.2: Sistema en bucle cerrado The system will be stable if 1 + G(s)H(S) has no zeros in the positive halfplane (the zeros 1 + G(s)H(s) are the poles of M (s)). Q (s + zi ) G(s)H(s) = k Q (s + pi ) Q Q Q (s + zi ) (s + pi ) + k (s + zi ) Q F (s) = 1 + G(s)H(s) = 1 + k Q = (s + pi ) (s + pi ) Applying the Cauchy argument principle to the function F (s), and taking a path that surrounds the positive half-plane (Nyquist path). If P is known (positive poles in the path ρ, in the unstable half plane) and N (turns to the origin of F (ρ)) it is possible to know Z (positive zeros of F (s), ie, unstable) making Z = N + P Since the zeros Z of F (s) are the closed loop poles system, they define whether the system is stable. If Z is nonzero, the system will be unstable.

137

F(s) = 1+G(s)H(s) ρ

F(ρ)

∞

Fig. 9.3: Metodo de Nyquist Modified method. One way to simplify the calculations in complex variable, maintaining the principle of the argument is to take the function F (s) = G(s)H(s), instead of 1 + G(s)H(s). Thus, the curve γ corresponding to F (ρ), remain displaced a unit on the real axis, and must now count the number N of turns around the point −1, instead from around the origin (Figure 9.4).

F(s) = G(s)H(s) ρ

F(ρ)

∞ -1

Fig. 9.4: Modified Nyquist method Example 9.0.1: Study the system stability of the Figure9.5 using the method of Nyquist. Solution: The transfer function in open loop is G(s)H(s) =

s+2 (s + 1)(s + 3)

The Nyquist path and the location of the zeros and poles of open loop is:

138

9. Nyquist Method

y(t)

s+2 s+1

e(t)

r(t) +

-

1 s+3 Fig. 9.5: System of the example 5.1.6

II

I x o

x

-3

-1

-2

III

Fig. 9.6: Nyquist path for the example 5.1.8 The image of section I is the polar diagram. Taking the limit for ω → 0 y ω → ∞ we have jω + 2 G(jω)H(jω) = (jω + 1)(jω + 3) lim |G(jω)| =

ω→0

2 3

lim arg(G(jω)) =

ω→0

lim |G(jω)| = 0

0◦

ω→∞

lim arg(G(jω)) = −90◦

ω→∞

The image of the section II will be the origin, since Num poles> Num zeroes. The image of the section III will be the symmetric section of I’ to the real axis. As the number of poles inside the Nyquist path is P = 0, and the image of the Nyquist path does not encircle the point -1 (N = 0), then Z = N + P = 0. That is, there are no zeros of the characteristic equation (poles of the closed loop transfer function) in the right half plane and therefore the system is stable (see Figure 9.6).

139

0.3 0.2 0.1 0 −0.1 −0.2 −0.3 −1

−0.8

−0.6

−0.4

−0.2

0

0.2

0.4

0.6

Fig. 9.7: Nyquist diagram of example 5.1.6 For systems with poles on the imaginary axis. If there are singular points in the Nyquist path, ie, if GH has poles at the origin or on the imaginary axis, it will be necessary to select a path that avoids these poles Nyquist. • Open loop system with pole at the origin. In this case you should add a new section (IV) around this pole. The section IV is defined by a semicircle of infinitesimal radius and centered at the origin:  jθ ε → 0
s = εe π θ ∈ −π 2 , 2

I

II IV

III

Fig. 9.8: Nyquist path for a system with a pole at the origin

140

9. Nyquist Method

• Open loop system with pure imaginary poles The new sections will now be: I s = jω ω ∈ (0, a) π jθ II s = aj + εe ε → 0 θ ∈ ( −π 4 , 4) III s = jω ω ∈ (a, ∞) IV s = Rejθ R → ∞ θ ∈ ( π4 , −π 4 ) V s = jω ω ∈ (−∞, −a) π V I s = −aj + εejθ ε → 0 θ ∈ ( −π 4 , 4) V IIs = jω ω ∈ (−a, 0)

Im III II

IV

I Re VII VI V

Fig. 9.9: Nyquist path for a system with poles on the imaginary axis

Open-loop stable systems (minimum phase systems) In the open loop stable systems poles are located in the left half of the plane s therefore P = 0. In this case, the Nyquist criterion is reduced to Z = N. Therefore, the system is stable if N = 0. To see if the system is stable just draw the image of the section I (polar diagram) and check that cuts on the right of the point -1 of the real axis. Most real systems have open loop transfer function of G(s)H(s) with all its poles and zeros in the left half plane. This kind of systems are called Minimum phase.

141

Then, since G(s)H(s) has no poles in the right half-plane, we have that P = 0, which will be fulfilled Z=N and therefore for the system to be stable must apply: N =0 In minimum phase systems the condition of stability means that the polar diagram of G(s)H(s) does not enclose the point -1

-1

Fig. 9.10: Stable system

-1

Fig. 9.11: Unstable system

Relative stability The Nyquist method allows the study of the relative stability of a closed loop system from the frequency response of the open loop system. This relative stability

142

9. Nyquist Method

is given by the proximity of the polar diagram at point s = -1 (see Figure ref estabrel1). In this figure the polar diagrams of three systems are shown. From the observation of them can be said for the system a) will be more stable than the b) and c) as its polar diagram is farthest from -1 than those of the other two systems (assuming all three are of minimum phase). To have a more quantitative and accurate picture of the relative stability, the concepts of gain margin and phase are used. Im(G(jω))

Im(G(jω))

-1

-1

Re(G(jω))

Re(G(jω))

b)

a) Im(G(jω))

-1

Re(G(jω))

c)

Fig. 9.12: The system a) is more stable than the b) and c systems)

Gain Margin M G. Provides a measure of the distance in magnitude between point -1 and the polar diagram of G(jω)H(jω). It is the inverse of the value of the polar diagram where this intersects the negative real axis (Figure 9.13). Kg =

1 A

The frequency at which this occurs is called emph phase crossover frequency, and is the frequency at which the phase is 180◦ . arg(GH(ωϕ )) = π

|GH(ωϕ )| = A

143

Im(G(jω))

ω=ωϕ

-1

A

Re(G(jω))

Fig. 9.13: Calculation of gain margin For a system to be stable is needed that Kg > 1. Kg > 1⇒sistema estable Kg = 1⇒sistema marginalmente estable Kg < 1⇒sistema inestable Expressed in dB; the gain margin M G is: M G = 20 log

1 = −20 log |G(ωϕ )H(ωϕ )| |G(ωϕ )H(ωϕ )|

For a stable system is verified |G(ωφ )H(ωφ )| < 1. Therefore, M G > 0. This value indicates how much gain can be increased without causing instability. For an unstable system is verified |G(ωφ )H(ωϕ )| > 1. Therefore, M G < 0, indicating there is to reduce the gain so that the system becomes stable. When the system is critically stable, |G(ωφ )H(ωφ )| = 1. Therefore, M G = 0, no gain can be increased without causing open loop instability. If the curve of the polar diagram does not intersect the negative real axis, as occurs in systems of first and second order, M G = ∞. Theoretically gain can be increased indefinitely without cause instability (FigureFor an unstable system is verified |G(ωφ )H(ωϕ )| > 1. Therefore, M G < 0, indicating there is to reduce the gain so that the system becomes stable. When the system is critically stable, |G(ωφ )H(ωφ )| = 1. Therefore, M G = 0, no gain can be increased without causing open loop instability. If the curve of the polar diagram does not intersect the negative real axis, as occurs in systems of first and second order, M G = inf ty. Theoretically gain can be increased indefinitely without cause instability (Figure9.14).

144

9. Nyquist Method

jω

Im(G(jω))

. ..

ω= ∞

-1

A

y(t)

K1

x

Re(G(jω))

K1

x

x

t

σ

K1

estable ω=0 Im(G(jω))

ω= ∞ -1 A

y(t)

jω

.

Re(G(jω))

K2

. .

K2

x

x

x

σ

t

K2

estable ω=0

jω

Im(G(jω))

K3

ω= ∞

.

K3

Re(G(jω))

A= -1

x

x

.

y(t)

x

σ

.

t

K3

ω=0

marginalmente estable

Im(G(jω))

A -1

.

jω ω= ∞

Re(G(jω))

.

x

x

x

σ

K4

inestable

y(t)

K4

t

.

K4

ω=0

Fig. 9.14: Nyquist plot, Root locus and closed loop step response according system 1 gain increases G(s) = s(s+a)(s+b) con a, b > 0. Phase margin γ. Figure 9.15 shows polar diagrams of two systems whose gain margin is the same. However, the System a) is more stable than the b). To quantify

145

Im(G(jω))

-1

Im(G(jω))

A

Re(G(jω))

-1

A

a)

Re(G(jω))

b)

Fig. 9.15: The two systems have the same profit margin, but a) is more stable than the b) this circumstance the phase margin is used. Gain crossover frequency (ωg ): the one for which |G(jω)H(jω)| = 1 (0 dB) The phase margin is the difference between the phase angle to ωg and 180◦ . M F = γ = arg(GH(ωg )) + π Graphically, γ is reflected in Figure 9.16. In stable and minimum phase open loop systems, the condition of stability of the closed loop system is that  M G|dB > 0 MF = γ > 0 Calculation of gain margin (Kg ) and phase margin (γ) Given the definition of gain margin, to calculate its value we will have to find the point of intersection with the real axis and find the value of phase margin as the constant by which multiply |G(jω)H(jω)| for the value 1, which is the module |(−1.0)|. In decibels (20log1 0 of the corresponding value) this product becomes sum. therefore arg(G(jωϕ )H(jωϕ )) = −180 |G(jωϕ )H(jωϕ )|dB + Kg dB = 0 Kg dB = − |G(jωϕ )H(jωϕ )|dB 1 Kg = |G(jωϕ )H(jωϕ )|

146

9. Nyquist Method

Im(G(jω))

Im(G(jω))

γ(−) ω=∞

ω=∞ Re(G(jω)) A

-1 A

ωg

-1

Re(G(jω)) ϕ

ϕ γ(+) ωγ

ω=0

ω=0

Margen de fase positivo

Margen de fase negativo

∞

Fig. 9.16: Phase Margin

ωg

0 dB

ω Kg

-180o

γ

ω ωϕ

Fig. 9.17: Gain margin and phase margin in Bode diagram

To find the phase margin must be calculated as the cut of the trace lef t|G(j omega)H(j omega) right|

147

with the unit circle and, then the frequency corresponding to this cut point. |G(jωg )H(jωg )| = 1 = 0 dB − arg(G(jωg )H(jωg )) − γ = −180 γ = 180 − arg(G(jωg )H(jωg )) The phase and gain margins are measured on the loop system opened. The gain/phase margin represents how you can increase the gain/phase shift of the open-loop system without the system closed loop becomes unstable.

148

9. Nyquist Method

10. DESIGN OF REGULATORS BY FRECUENCY METHODS The design of regulators for frecuency methods are based on the concepts of phase margin and gain margin in relation to Nyquist and Bode diagrams. The margins of phase and gain are measured in the open loop system. Phase margin and gain represent how much you can increase the gain or phase of the open loop system without turning the closed loop system unstable. The gain margin is the ratio of the distance from the origin to the point −1 + 0j and the distance to the cutoff of G(jω)H(jω) with the real axis.

Fig. 10.1: In this cut point of G(jω)H(jω) with the real axis, the argument is zero: arg(G(jωϕ )H(jωϕ )) = −180◦ and its magnitude

150

10. Design of regulators by frecuency methods

Fig. 10.2:

|G(jωϕ )H(jωϕ )|dB + Kg dB = 0 Therefore, Kg dB = −|G(jωϕ )H(jωϕ )|dB and natural scale (not decibels) Kg =

1 |G(jωϕ )H(jωϕ )|

(10.0.1)

On the other hand, if we consider cutting G(jω)H(jω) with the unit circle, at this point the module is one. |G(jωg )H(jωg )| = 1 = 0 dB The difference to the argument of −1 is the missing phase to reach the instability. arg(G(jωg )H(jωg )) − γ = −180◦ Therefore, the phase margin can be calculated as γ = 180◦ + arg(G(jωg )H(jωg ))

(10.0.2)

10.1. Relationship between parameters of relative stability and transient response

10.1

151

Relationship between parameters of relative stability and transient response

The approximate relationship between the parameters of relative stability and transient response γ ≈ 100ζ para 0 < ζ < 0, 7 (10.1.1) ts ≈

2π 1 = ωg tan(γ) fg tan(γ)

(10.1.2)

where • ζ It is the damping ratio • ts It is the settling time Conclusions: In the underdamped systems, the phase margin (γ) is directly related to the damping ratio (ζ). When the phase margin decreases, overshoot increases. When the cutoff frequency gain (ωg ) and the phase margin (γ) increase the system becomes faster (the time of settlement (ts ) decreases).

10.2

Lead Networks

The RC network shown in Fig. is equivalent to a phase lead controller We denote the impedances Z1 = R1 / (1 + R1 Cs)

y

Z2 = R2

The transfer function of the phase-lead network can be written as

152

10. Design of regulators by frecuency methods

Fig. 10.3: RLC lead network and its corresponding diagram of poles and zeros.

Fig. 10.4: Bode and Nyquist plot for a lead network according to different parameter values α.

Fig. 10.5: Bode diagram of the system and lead network.

10.2. Lead Networks

153

Fig. 10.6: Comparison of Bode plots and unit step responses.

Fig. 10.7: Adding a lead network to the system changes the magnitude and phase, so it is difficult to predict the new crossing point.

154

10. Design of regulators by frecuency methods

Fig. 10.8: Designing a lead network.

10.3. Lag networks

10.3

Lag networks

Fig. 10.9: Lag network and pole-zero diagram.

Fig. 10.10: Frequency response of lag networks.

Fig. 10.11: Comparison of the response of the systems.

155

156

10. Design of regulators by frecuency methods

Fig. 10.12: Designing a lag network.

Fig. 10.13: The effect of the change of T.

SOME HISTORICAL DATA

158

10. Design of regulators by frecuency methods

Jacopo Francesco Riccati (1676 Venice (Italia) – 1754 Treviso (Italia))

Taylor, Brooks (1685 Edmonton (Reino Unido) – 1731 Londres (Reino Unido))

Leonhard Euler (1707 Basilea (Suiza) – 1783 San Petersburgo (Rusia))

Watt, James (1736 Greenock, land) – 1736)

Laplace, Pierre Simon (1749 Normandia (Francia) – 1827 Paris (Francia))

Fourier, Jean Baptiste Joseph (1768 Borgoña (Francia) – 1830 Paris (Francia))

(Scot-

10.3. Lag networks

159

Cauchy, Agustin Louis (1789 Paris (Francia) – 1857 Paris (Francia))

Sir William Rowan Hamilton (1805 Dublin (Irlanda) 1865 Dublin (Irlanda))

Routh, Edward John (1831 Quebec (Canada) – 1907 Cambridge (Reino unido))

Jordan, Marie Ennemond Camille (1838 Lyon (Francia) – 1922 Paris (Francia))

Liapunov, Aleksandrer Mijailovich (1857 Yaroslavl (Rusia) – 1918 Odessa (Rusia))

Hurwitz, Adolf (1859 Hanover (Germany)1919 Zurich (Switzerland))

160

10. Design of regulators by frecuency methods

Evans, Griffith Conrad (1887 Boston (USA) – 1973 Berkeley (USA))

Nyquist, Harry (1889, Nilsby (Suecia) – 1976 (USA))

Norbert Wiener (1894, Columbia (USA) – 1964 Estocolmo (Suecia))

Bode, Henrick W. (1905 Madison (USA) – 1982 Maryland (USA))

Pontriagin, Lev Semionovich (1908 Moscú (Rusia) – 1988 Moscú (Rusia))

Shannon, Claude E. (1916 Michigan(USA) – 2001 Massachusetts (USA))

10.3. Lag networks

161

Yakov Zalmanovitch Tsypkin ( 1919 (Rusia) – 1997 (Rusia))

Richard Bellman (1920 Brooklyn USA) – actualidad)

Kalman, Rudolf E. (1930 Budapest (Hungria) – actualidad)

Truxal, John G.

Eliahu I. Jury

Karl J. Åström (1934 Östersund (Sweden) - actualidad)

INDEX , 60, 101, 104, 107

Small gain theorem, 73

Bode plot, 124

Time models, 14 Transfer function, 15

Cauchy argument principle, 133

undamped natural frequency, 60

Final value theorem, 5 Fourier transform, 1

Ziegler-Nichols, 104 Ziegler-Nichols methods, 104

Gain Margin, 140 impulse response, 15 Initial value theorem, 5 Input/Output models, 14 Inverse Fourier Transform , 3 Inverse Laplace transform, 7 Laplace Transform, 3 method based on algebraic determination for adjusting PID, 111 Method based on the dominant poles for adjusting PID, 110 Method based on the Root Locus for adjusting PID Method, 114 Nyquist, 133 Nyquist stability criterion, 133 Overshoot, 61, 102 Peak, 61 Peak time, 101 Phase margin, 142 PID controllers, 102 Routh’s method, 74 Settling time, 61, 102 162