Differential Equations and the Laplace Transform

A differential equation is an equation involving a function and its
derivatives. Eg.

dx
-- = 4
dt

is a differential equation. Thus a differential equation is an
equation whose solution is a function, or a family of functions. The
above equation can be simply solved by integration, giving a value for
x

where c is an arbitrary constant - ie the solution is a family of
functions parameterised byc.
An ordinary differential equation (ODE) is a differential equation
in which differentiation is with respect to one variable only.
A partial differential equation is a differential equation in
which differentiation with respect to several variables may occur.
Eg.

is a partial differential equation which describes many forms of wave motion.
Robotics, for the most part, makes use of Ordinary Differential Equations,
where differentiation takes place with respect to time.

Clearly, differential equations of the form

dx
-- = f(t)
dt

where f(t) is a known function of t can be solved by integration, which
is not itself an operation that can always be performed analytically -
sometimes integration gives rise to new functions. Thus solving differential
equations is at least as hard as integration, and may sometimes have
to be performed numerically. However for the control of robots (as
in any other problem of Control Theory) we are not so much concerned
to solve a differential equation as to choose a differential equation
whose solution has certain desirable properties.

Linear Differential Equations

A differential equation is linear is linear if the unknown function(s)
and its (their) derivatives occur linearly in the equation. Eg.

is linear, while

is not. Homogeneous
linear differential equations have the property of superposition
of solutions. Given an equation

??? an image missing here

an(t) \dxbydtn + a(t)n-1
\dxbydtnpred + ... +a0(t)x = 0

if x1(t) and x2(t)
are solutions, and if c1c2 are constants, then
c1 x1(t) + c2 x2(t) is also a solution. (This is immediate from the
fact that differentiation itself is a linear operator).
The same property holds for simultaneous linear differential equations.

Linear differential equations with constant coefficients (ie where
the an(t) are constants) are by far the most tractable kind of differential
equations, and there is a rich body of knowledge about them. Therefore
we try to cast our problems in a form in which these are applicable.
This can be approached through suitable robot design (for example reducing
friction and backlash). However, the differential equations for the
dynamics of any robot with revolute joints will be non-linear, but for
control purposes the technique of linearising the equations
can be effective since the purpose of a controller is to ensure that the
controlled system deviates from its intended behavior only by a small amount.

A simple example of linearising a differential equation can be seen
in the case of a pendulum. The differential equation of motion is

whose exact solution requires elliptic functions. However the
equation can be linearised, by approximating \sin(\theta) by \theta
(measured in radians of course). In this case we get

m\dbydtsq{\theta} = \frac{m}{l}\theta$$

which is linear, with constant coefficients.

The Thermostat, a Simple controller.

Consider a room, of heat capacity H, heated by a heater, of
capacity W watts, We will assume that the temperature in the room
has a uniform value T, and that heat is lost to the outside at a
rate that depends on the temperature difference T-T0, of
L watts per degree centigrade. If the heater is turned on
then the differential equation for T is:-

while if it is turned off, then the ODE for T is:-

??? missing image
$$ \dbydt{T} = \frac{ - L(T-T0)}{H} $$

Both of these are linear differential equations for T, with
constant coefficients. The first equation can be written
(with a positive denominator on the LHS):-

$$\frac{1}{a-T}\dbydt{T} = b$$

Integrating with respect to t we get:-

$$\int\frac{dT}{a-T} = \int b\, dt$$

So that:-

$$ -\ln(a-T) = bt+c$$

$$a - T = e^{-bt-c}$$

?????? picture missing here

$$T = a - de^{-bt}$$

where d=e^{-c}. That is T tends to a from below.

In the second case we have

$$\int\frac{dT}{T-a } = -\int b\, dt$$

with a positive denominator. Integrating, we get

$$ \ln(T-a) = -bt+c$$
$$ T = a + e^{-bt+c} = a + ke^{-bt} $$

So that T approaches a from above. Note that both the heating up
and cooling down equations have the same form if we allow a positive
or negative coefficient of e^{-bt}.

In the cooling case a=T0, so that the room temperature
approaches T0 in the limit. In the heating case, let
Tmax
= T0+W/L, then T tends to Tmax in the limit.

If now we install a thermostat, which turns off the heater when the
T>T2 and turns it on when T1, where T12, then
\begin{enumerate}
\item If T1 < T0 then
the room will stabilise to T0 with the heat turned off.
\item If T2 > Tmax then the room will stabilise to Tmax with the
heat turned on
\item Othewise the room will go into a limit cycle, in which the
temperature rises to T2, the thermostat turns off, the temperature falls
to T1, the thermostat turns on, and the cycle continues.
\end{enumerate}

It is important to remember that a controller can only accomplish
what is physically possible. A thermostat cannot by itself cool a room
down below the outside temperature- you need an air conditioner.

Laplace Transformations

We saw above that the solution of a what is a rather simple differential
equation is rather hard work, and might despair of dealing with the complicated
equations that can arise when we have multiple motors acting on
a robot through flexible drive trains. In the case of linear ODE's with
constant coefficients there is a very effective treatment through
the following:--
Definition

If f is a function, then the Laplace transform
is defined by:-

$$F(s) = \laptrans{f(t)}$$

Clearly, there will be functions for which the Laplace transform is not
defined because the integral does not converge, but such functions are
not of significance from the point of view of robotics.

The reason Laplace transforms are so useful is in view of the following
result, obtained by integration by parts:-

Now, provided e-stf(t) tends to zero as t increases
unboundedly for some s this will evaluate to

This will indeed happen for any functions of engineering interest.

Thus, if we arrange, by adding suitable offsets if required, that all
functions we apply the Laplace transform have the property that f(t)=0,
we see that we can map the operation of differentiating
a function into
multiplying its Laplace transform by s.

Since, for practical purposes, functions are uniquely determined by
their Laplace transform for the values of t>0, we can translate problems
in differential equations into manipulations of Laplace transforms, which
are often simple polynomials.

Laplace transforms can be determined by standard integration techniques.
Eg.

The Laplace Transforms of Some Important Functions.

In this section we tabulate the Laplace transforms of a number of functions.
For control, important functions are the unit step:-

and the unit ramp:-

The following table can be obtained by standard integration techniques
(the sin and cosine formulas require a double integration by parts).

Note that the Laplace transforms of all of these functions are
rational functions of s. Moreover we can see that:-

(consider the differentiation rule and notice that the integral is zero
if t=0).

So we see that for the whole vocabulary of elementary functions, and
operations of differentiation and integration on them, the
corresponding Laplace transforms are rational functions of s.
Now any rational function over the reals can be expressed as a sum of
partial fractions with linear or quadratic demoninators, or powers of the
same. (Apply the fundamental theorem of algebra, and group complex
conjugate pairs of roots together.) Thus we begin to see that it may be
possible to find a function whose Laplace transform is any given rational
function of s, by converting into a sum of partial fractions and
then finding functions whose Laplace transform is each partial fraction.

The Inverse Laplace Transform.

The Inverse Laplace Transform
allows us to translate from the s domain (or frequency domain) back to
the t domain, or time domain. While there is a definition of the inverse
transform as a complex integral, in practice it is more useful to
derive it from the table of forward transforms, together with
the completion of the square and the following two properties of
the forward transform in which, as usual
.

The shifting property.

{\cal L}(eatf(t)) = F(s-a)$$

This result is readily obtained by change of the variable of integration.

The tn property.

This result can be shown by differentiating under the integral sign.

With these results under our belt we can specify the following steps to
finding the inverse transform for a rational function:-

Where N(s) and D(s) are polynomials in s.
We may assume that the highest power of s in D(S)
is of the form sn (for if not we can make it so by dividing N(s) and
D(s) by the coefficient of sn in D(s)).
\begin{enumerate}
\item Factor D(s) into linear complex factors

In general this step will have to be performed numerically, and requires
D(s) to have numeric coefficients. Some algebraic manipulation programs
(eg. Macsyma) provide symbolic factorisation algorithms, but they cannot
always work (Galois). Note that finding eigenvectors of matrices also
requires the factorisation of polynomials, and it should therefore be no
surprise that a matrix formulation of the problem of solving Linear ODE's
with constant coefficients exists.
Thus D(s) will have factors of the form s-ai. If ai is complex,
then there will be a corresponding factor of the form s-\bar{ai},
where \bar{ai} is the complex conjugate of ai.
\item Combine the complex conjugate factors together to form real quadratic
factors. Thus

Where some of the factors may be repeated, and all of the a'i are real.
In practice these first two steps will usually be combined.
\item Express F(s) as partial fractions.

\item Form f(t), the inverse transform, as

f(t) = \sum fi(t) + \sum gi(t) $$

where the fi are derived by from the terms with linear denominators, the
gi from those with quadratic denominators. Then, if ni = 1 :-

fi(t) = pi e^{a'it}$$
or if ni>1 :-

fi(t) = \frac{pi}{(ni-1)!} t^{ni-1} e^{a'it} $$

Proceeding to the quadratic terms we first complete the square in the
denominator, obtaining

How to Treat Controllers with Laplace Transforms.

Let us now consider a simple controller used this
time to control a motor which actuates a joint of a 1-link robot.
Recall that the thermostat was
unable to control the room temperature very precisely. Perhaps if we used
a continuously variable amplifier instead of a discrete switch we might
get better control. This is indicated in the following diagram:-

Here the amplifier is supposed to be a current amplifier, ie one
whose output current is specified by the input (this simplifies our
description of the motor).
Now the motor produces a torque
related to the current i
by

where
cM is a constant of the motor. Let the inertia of the
link be J, and the angle of the joint be y. Then the motion of the
joint is described by the angular analog of Newton's first law:-

Now the current i is commanded by the controller as follows:-

i = g(x-y)$$

Where g is the gain. (we assume that the amplifier converts 1 unit
of output from the gain to 1 unit of current). Thus we have the
differential equation for y

For convenience we let \omega2=\frac{gcM}{J}.
This gives

\ddot{y} = \omega2(x-y)$$

If x=0 this will be recognised as the equation of simple harmonic motion -
so it seems that our controller will oscillate, which is not a desirable state
of affairs. Let us verify this by transforming the equation into the
s-domain.

Solving for Y(s) we get:-

Thus if we choose x(t) to be the unit step u(t), we get

Expanding as partial fractions we get:-

From which we see, by equating coefficients, that:-

Transforming back into the t-domain, we obtain:-

y(t) = u(t) - \cos\omega t$$

That is to say, y differs from x by the term \cos\omega t,
so that the link oscillates for ever about the goal position.
We can see how to improve the controller by considering a pendulum,
or mariner's compass, whose motion is damped by being placed
in fluid. The fluid exerts a viscous friction which slows down the
oscillations. For certain conditions it is
reasonable to treat this damping as being proportional to the
velocity. (You may know that sub-sonic aerodynamic form drag
is taken to be proportional to the square of the velocity).

Likewise we can damp out the oscillations of our controlled
link by subtracting off a term proportional to the velocity. The velocity
can either be found from a tachomoter, or by differentiating the position
(but this is a noisy estimate).

In the critically damped case, we have a perfect square down below, so
we can either apply the tn rule, or treat it as a limiting case
of underdamped behavior, remembering that \cos(x)\rightarrow 1 and
sin(x)/x \rightarrow 1 as x\rightarrow 1. Thus we obtain, for the
critically damped case:-
NOTE write it as a limiting condition.

Thus y'(0+) = 0 , so the system starts at rest, and y'(t) >0 when t>0,
so that y(t) has no maximum value, but increases monotonically to the
limit 1 as t tends to infty, i.e. there is no overshoot.

In the overdamped case, the
denominator will factor into two real linear terms in s, so that
the behavior in the t-domain will be an exponential approach
to u(t). The \sinh and \cosh laplace transforms provide the
neatest derivation. Let \omega"2=\omega2(1-\xi2). The behavior in
the t-domain is then

Now the \cosh(t)\rightarrow\infty, \sinh(t)\rightarrow\infty
as t\rightarrow\infty, so we have to satisfy ourselves that the system is
stable. Using the definition cosh(x) = (e^x+e^(-x))/2

\begin{Large}NOTE \end{Large}
If two functions have the same Laplace transform they must have the same
value ``almost everywhere" on the positive real axis.
This means that they may differ at finite number of points (or more
depending on the definition of integration employed,
eg. is it Lebegue Integration)

Steady-State

The control systems we have studied so far have had the property that,
provided they are stable, if they are given a step input, they approach
the new commanded value asymptotically, ie if x(t) = au(t) then
\lim_{t \rightarrow \infty} y(t) = 1 This behavior is clearly desirable
for a robot, but it is not always attainable. Suppose for example that
our link is loaded by gravity. Then in order for the link to be stationary,
the motor must exert a torque to resist the gravitational force, and therefore
there must be a non-zero current in the motor. But i = g(x-y)+k\dot{y}
so that a non-zero current can only exist with a stationary link if
there is a non-zero error.
Now if g is large, the error will be small, but we have remarked that
we may want to position a robot to one part in 10000, and that we cannot
have too large a gain or the system will be unstable. (Also if we are
digitising our sensor input, we will not be able to position our robot
exactly to one quantisation level, since to get any current we will have
to have at least one quantisation level of error). There is more than
one solution to this problem. For example we could :-
\begin{enumerate}
\item Calculate what the current required to overcome gravity is, and
add this into the command to the motors. In the case of gravity, this
quantity depends only on the joint angle y, so call it p(y)
Then we use the equation

i = g(x-y) +k\dot{y}+ p(y)$$ If we get the sums
right, then gravitational loads no longer produce any error. The approach can
be extended to handle dynamic loads arising from trying to track a changing
input (eg a ramp x(t)=a\,r(t)).
However it will not work for friction, since the value of frictional forces
cannot be predicted with any accuracy, and there is an inconvenience if the
robot is carrying a load, since the static (and dynamic) characteristics
of the load have to be known.
\item Consider the integral of the error with respect to time. If the error
is of a single sign, and non-zero, then its integral is a monotonic increasing
function. If x(t)-y(t)>\epsilon>0 where \epsilon is a constant, then
\int(x-y)dt will increase without limit. So if we employ a controller
for which i = gP(x-y) - gV\dot{y} + gI\int(x-y)dt,
we can be sure that if we give it
a step command, the error will not
remain greater than some constant
\epsilon. It may, of course, oscillate, and will do so if we choose too
large a value for gI.
\end{enumerate}

A controller which employs Position, its Derivative, and Integration is
known as a PID controller. The earlier controller which employed just
the position and derivative of y is known as a PD controller.

The limiting behavior of a linear controller can be derived using the
following Final Value Theorem.

If F(s) = {\cal L}(f(t)) then

\lim_{t\rightarrow\infty} f(t) = \lim_{ s\rightarrow 0} sF(s)$$

provided that F(s) is analytic in the right half plane including the
j axis. (For a rational F(s) this means that it must have no
poles in the right half plane). Thus, applying this theorem to our
PD controller with a unit step input, we get:-

If x(t) = u(t) then ,applying the final value theorem,
we see that y(t) tends to 1
as t\rightarrow\infty, provided the denominator has no roots in the right
half plane. Now if gI is large, we can see thet the denominator behaves as
s^3+gI which has roots s=\alpha e^{j(\pi/3+2k\pi/3)} where
\alpha^3=gI. k=1 gives the root s=-\alpha, while k=0 and k=2 give
two roots in the positive half plane. Thus for a large gI we may expect
the system to be unstable. If gI = 0 we have a root s=0 and the
two roots of the denominator of the system analysed earlier, which will be
complex in the underdamped case. As gI is increased from zero, the
real root will start to move along the negative real axis, and the
complex roots will initially remain in the left half plane, so we have a
stable system for small values of gI (to show that the real root moves this
way, hold all coefficients constant except gI and consider
\frac{ds}{dgI}, which is negative if s=0).

Practical Considerations

An analog velocity loop

The position of robot joints is usually controlled digitally since this
provides greater accuracy and flexibility. However, unless great computer
power is used to give short sample times,
it may be advantageous to provide velocity feedback through
an analog controller. Recall that a long sample time can mean instabilility.
If we restructure the equation

i = g(x-y) + k\dot{y}$$ as

v= g1(x-y),\;\;\;\; i = k(v-\dot{y})$$

Where g1 = g/k
We can regard the second equation as defining a velocity controller
which is fed a velocity command by the first. This arrangement is
more stable if the velocity control is performed by an analog controller
which by its nature imposes no significant delay,
and the position control is performed by a digital controller.
The precision of velocity control does not have to be
great for most robotic use, since precision of position is normally
what is required, so that the analog velocity loop is very adequate from
that point of view.

Saturation of the Controller

While digital controllers can readily represent a very large range of numbers,
the analog components of a system, the power amplifiers and the motors,
will only have a limited range in which their behavior is linear.
If the input to an amplifier goes outside the operational range the device
will saturate, that is to say it will give a constant output.
A motor (or at least vital parts of it) may well melt. In practice,
the gains of a controller are such, that if a step input is given to
move any significant fraction of the full range of movement of a joint,
the amplifiers will saturate. This may be acceptable in a single link device,
provided the amplifiers are properly designed to go into saturation, and
properly matched to the motors. However once a controller goes into
saturation, it is no longer controlling. And if we have multiple joints on a
robot, they will behave in a very uncoordinated manner if we give them
all large step inputs (eg. to move the robot from one place to another), so
that the path taken by the device will be unpredictable.

A practical solution to this problem, suitable for moderate performance
robots, is to interpolate
in joint space. If we have joints with angular or
linear displacements (y1,y2,...,yn) then we can choose command
functions (x1(t)...xn(t)) to go from a point (a1...an)
to (b1...bn) by a straight line.

\vec{x(t)} = (1-\lambda(t))\vec{a} + \lambda(t)\vec{b}$$

There remains a problem of choosing the function \lambda. One possibility
is to have \lambda(t)=lt for some constant l, which can be chosen by the
programmer of the robot. He will typically develop the program with a small
value of l and then increase it until the robot starts behaving improperly.
This is not perhaps a very scientific approach, but it is often done.
If we have a kinematic transform
from cartesian space to joint space,
the contoller can work out a straight line in cartesian space, and transform
that to a sequence of values for x(t). Note that there is no
guarantee
that y(t) will track x(t) exactly. Indeed, if \lambda(t)=lt
then it is impossible to obtain exact tracking, since a robot cannot
start moving instantaneously. However this form of interpolation does provide
a practical approach to controlling a robot which giving step inputs does not.

One of the most important aspects for the controllability of a robot is
what the lowest resonant frequency of vibration is. One can think of a robot
as a device that converts a control signal to a movement. Consider
for example a robot that is to pick up an object at point A and put it
down at point B, and repeat the action every second. If it has a resonant
frequency below 2hz it will be very difficult to control it to perform this
task.