# Category Archives: engineering

## Matrix inversion, row ops, program added

The row operation matrix inversion method is so neat and ingenious, and it has the same operations for all dimensions of matrix.

Here is a step by step approach, where firstly the dimension is chosen, then the first of the buttons is selected (Start). After which the buttons are selected in order. The states of the left and right matrices are displayed at each stage, and finally the identity matrix appears on the left, and the inverse of the original matrix appears on the right.

http://mathcomesalive.com/mathsite/matrix%20prog%203.html

The first display shows the original matrix on the left. Nothing has been done yet.

The second one is the 2 x 2 matrix inverted.

The third is the 4 x 4 matrix inverted.

The matrices can be altered with the file name for the application, “matrix prog 3.html”. using right click and “view page source”. You are then on your own, with javascript !!!!!!!!

## Part 5: James Watt, the steam engine, and the first automatic control system

I was going to do my own version of the James Watt flyball governor, but I found a very readable, illustrated account of the life and inventions of James Watt. Here it is:

http://www.mainlesson.com/display.php?author=bachman&book=inventors&story=watt

In the picture below you can see the flyball governor, driven by the output shaft of the engine. The link from the governor to the steam valve however is not clear, but it is there.

Here is a detailed view of the governor.

## Part 4: Tuning the feedback controller

Our first order process is described by the equation yn+1 = ayn + kxn , where yn is the process output now (at time n),  xn is the process input between time n and time n+1, and yn+1 is the process output at time n+1.
a is the coefficient determining how quickly the process settles after a change in the input, and k is related to the process steady state gain (ratio of settled output to constant input).

If we set the input to be a constant x and the output settled value to be a constant y, then

y = ay + kx, and solving for y/x we get y/x = k/(1-a), the actual steady state gain.

In what follows the k will represent the actual steady state gain, the old k divided by (1-a)

The fist two plots show the process alone, with input set to 6 at time zero. In the first the steady state gain is set to 1, and in the second it is set to 2

Now we look a t the process under “direct” control, where the input is determined only by the chosen setpoint value. The two equations are

yn+1 = ayn + kxn  and  xn = Adn  (direct control: h is zero)

To obtain a controlled system with overall steady state gain equal to 1 (settled output  equal to desired output) it is easy to see that  has to be equal to (1-a)/k

It is not so obvious how the choice of h affects the performance of the controlled system. To do this we observe that the complete system is described entirely by the process equation and the controller equation together, and we can eliminate the xn from the two equations to get yn+1 = ayn + k(Adn + h(yn – dn))

which rearranged is yn+1 = (a + kh)yn + k(A – h)dn

Substituting   (1-a)/k for A, as found above, gives  yn+1 = (a + kh)yn + (1 – a – kh)dn

which has the required steady state gain of 1.

This final equation has the SAME structure as the process equation,
with a + kh in place of a

So now we will see how the value of the “a” coefficient affects the dynamic response of the system.

If h > 0 the controlled system will respond slower than with h = 0, and if h < 0 it will respond faster:

Setpoint changes were made at time 20 and at time 40

Congratulations if you got this far. This introduction to computer controlled processes has been kept as simple as possible, while using just the minimum amount of really basic math. The difficulties are in the interpretation and meaning of the various equations, and this something which is studiously avoided in school math. Such a shame.

Now you can run the program yourself, and play with the coefficients. It is a webpage with javascript:  http://mathcomesalive.com/mathsite/firstordersiml.html

Aspects and theoretical stuff which follow this (not here !) include the  backward shift operator z and its use in forming the transfer function of the system, behaviour of systems with wave form inputs to assess frequency response, representation of systems in matrix form (state space), non-linear systems and limit cycles, optimal control, adaptive control, and more…..

1 Comment

## Part 3: Computer control of real dynamic systems

Let us return to the industrial steam engine driving a range of equipment, each piece in the charge of an operator. In order to keep the factory working properly it is necessary to keep the speed of the steam engine reasonably constant, at what is called the desired output, or the set point. A human controller can achieve this fairly well, by direct observation of the speed and experience of the effect of changing the steam flow rate, but automatically? Well, James Watt solved the problem with a mechanical device (details later). We will now see how to “do it with a computer”.

Diagram of the physical system. There will be a valve on the steam line (not shown).

The information flow diagram. The load is not measured, and may vary.

Now we have a controller in the picture, fed with two pieces of information, the current speed and the desired speed. The controller can be human, mechanical or computer based. It has a set of rules to figure out the appropriate steam flow rate.

Now we have computer control. The output of the system is measured (sampled) at regular intervals, the computer calculates the required flow rate , sends the corresponding value to the valve actuating mechanism, which holds this value until the next sampling time, when the process is repeated. Using n=1, 2, 3, ….. for the times at which the output is sampled and the controller does its bit, we have at time n the speed measure yn, the desired speed dn and the difference between them (or error in the speed), and they are used to calculate the system input xn with the formula shown in the controller box.

The equation in the steam engine box is our fairly simple first order linear system model as described in the previous post. It is a good idea to have a model of the process dynamics for many reasons (!), one of which is that we can do experiments on the whole controlled system by simulation rather than on the real thing.

Looking at the controller function (formula), xn = Adn + h(yn – dn) , it is “obvious” that if the coefficient h is zero we have direct control (no feedback), and so the value of the coefficient A must be chosen accordingly (see next post).

It is not so obvious how the choice of h affects the performance of the controlled system. To do this we observe that the complete system is described entirely by the process equation and the controller equation together, and we can eliminate the xn from the two equations to get yn+1 = ayn + k(Adn + h(yn – dn))

which rearranged is yn+1 = (a + h)yn + k(A – h)dn

and has the SAME structure as the process equation, with a + h in place of a

In the next post, on tuning our controller, we will see how the value of the “a” coefficient affects the dynamic response of the system.

1 Comment

## Time Series and Discrete Processes, continued

A simple function has an input, call it x, and an output, f(x). Sometimes the input is t, for time, and the output is the value f(t), observed at time t. With such a function, or process, we can look at what happened in the past, but we cannot know for certain what will happen in the future. Of course, we may have a (mathematical) model of the process, which allows us to make predictions about the future. These processes are “simple”, insofar as the future behaviour depends only on the value of the time variable.

An example is the height of a projectile, s(t) = ut – 0.5gt2

Borrowed from dy/dan “Will it hit the hoop?”

Far more interesting are what I shall call “dynamic” processes, in which the future behaviour depends on how the process output got to its current value (at time now!). This implies that in the simplest case the rate of change of the process output is involved.

This is an early industrial static steam engine, designed by James Watt. The input is steam flow, the output is the rotational speed of the flywheel.

An example is a DC (direct current) electric motor driving a load. The way its speed changes over time depends not only on the voltage applied, but also on the speed at the moment (now).

This setup can be modeled by a simple first order differential equation,  but since computers came along it became clear that a simpler model using difference equations would do the job as well.

Let time be seen as equally spaced points, labeled …n-1, n, n+1, … and the speed of the motor at time n be vn.

Then the model of the process  is vn+1 = avn + kun, where un is the  voltage input at time n. The values of a and k are to be determined from theoretical or observed behaviour of the motor.This model is very good if we choose to keep the input voltage constant for the duration of the time interval.

Now for the amusing bit. I had taught business students about exponentially weighted averages numerous times previous to my sabbatical year studying control systems, and when I finally became happy with the discrete models as described above I realized in one of those moments that don’t happen very often that the equations were THE SAME in both cases. A bit of rewriting needed !

ewma(n) = (1-k)*ewma(n-1) + k*x(n)
vn+1 = avn + kun

Next post will look at feedback and computer control of dynamic systems.

## Time Series and Discrete Processes, or Calculus Takes a Back Seat

In the business and manufacturing world data keeps on coming in. What is it telling me? What is it telling me about the past? Can it tell me something believable about the future. These streams of data form the food for “time series analysis”, procedures for condensing
the data into manageable amounts.

Here is a plot of quarterly Apple I-phone sales. (ignore the obvious error, this is common)

I see several things, one is a pronounced and fairly regular variation according to the quarter, another is an apparent tendency for the sales to be increasing over time, and the third is what I will call the jiggling up and down, more formally called the unpredictable or random behaviour.
The quarterly variation can be removed from the data by “seasonal analysis”, a standard method which I am not going to deal wirh (look it up. Here’s a link:

https://onlinecourses.science.psu.edu/stat510/node/69
The gradual increase, if real, can be found by fitting a straight line to the data using linear regression, and then subtracting the line value from the observed value for each data point. This gives the “residuals”, which are the unpredictable bits.

Some time series have no pattern, and can be viewed as entirely residuals.
We cannot presume that the residuals are just “random”, so we need a process of smoothing the data to a) reveal any short term drifting, and b) some help with prediction.
The simple way is that of “moving average”, the average of the last k residuals. Let’s take k = 5, why not !
Then symbolically, with the most recent five data items now as x(n), x(n-1),..,x(n-4), the moving average is

ma = (x(n)+x(n-1)+x(n-2)+x(n-3)+x(n-4))/5

If we write it as ma = x(n)/5 + x(n-1)/5 + x(n-2)/5 + x(n-3)/5 + x(n-4)/5 we can see that this is a weighted sum of the last five items, each with a weight of 1/5 (the weights add up to 1).
This is fine for smoothing the data, but not very good for finding out anything about future short term behaviour.
It was figured out around 60 years ago that if the weights got smaller as the data points became further back in time then things might be better.

Consider taking the first weight, that applied to x(n) as 1/5, or 0.2 so the next one back is 0.2 squared and so on. These won’t add up to 1 which is needed for an average, so we fix it.
0.2 + 0.2^2 + 0.2^3 + … is equal to 0.2/(1-0.2) or 0.2/0.8 so the initial weights all need to be divided by the 0.8, and we get the “exponentially weighted moving average” ewma.

ewma(n) =
x(n)*0.8+x(n-1)*0.8*0.2+x(n-2)*0.8*0.04+x(n-3)*0.8*0.008+x(n-)*0.8*0.0016 + …… where n is the time variable, in the form 1, 2, 3, 4, …

A quick look at this shows that ewma(n-1) is buried in the right hand side, and we can see that ewma(n) = 0.2*ewma(n-1) + 0.8*x(n),
which makes the progressive calculations very simple.
In words, the next value of the ewma is a weighted average of the previous value of the ewma and the new data value. (0.2 + 0.8 = 1)

The weighting we have just worked out is not very useful, as it pays too much attention to the new data, and the ewma will therefore be almost as jumpy as the data. Better values are in the range 0.2 to 0.05 as you can see in the following pictures, in which the k value is the weighting of the data value x, and the weighting of the ewma value is then (1-k):

So the general form is ewma(n) = (1-k)*ewma(n-1) + k*x(n)

With k=0.5 we do not get a very good idea of the general level of the sequence, as the ewma value is halfway between the data and the previous ewma value, so we try k=0.1

Much better, in spite of my having increased the random range from 4 to 10.

Going back to random=4, but giving the data an upward drift the ewma picks up the direction of the drift, but fails to provide a good estimate of current data value. This can be fixed by modifying the ewma formula. A similar situation arises when the data values rise at first and then start to fall (second picture)

Common sense says that to reduce the tracking error for a drifting sequence it is enough to increase the value of k. But that does not get rid of the error. We need a measure of error which gets bigger and bigger as long as the error persists. Well, one thing certainly gets bigger is the total error, so let us use it and see what happens.

Writing the error at time n as err(n), and the sum of errors to date as totalerr(n) we have

err(n) = x(n) – ewma(n), and totalerr(n) = totalerr(n-1) + err(n)

Then we can incorporate the accumulated error into the ewma formula by adding a small multiple of totalerr(n) to get

ewma(n) =  ewma(n) = (1-k)*ewma(n-1) + k*x(n) + h*totalerr(n)

In the first example below h = 0.05, and in the second h = 0.03, as things were too lively with h=0.05.

A good article on moving averages is:
http://www.mcoscillator.com/learning_center/kb/market_history_and_background/who_first_came_up_with_moving_averages/

In my next post I will show how the exponentially weighted moving average can be described in the language of discrete time feedback systems models, as used in modern control systems, and with luck I will get as far as the z-transfer function idea.

## Real problems with conic sections (ellipse, parabola) part two

So suppose we have a parabolic curve and we want to find out stuff about it.

Its equation … Oh, we have no axes.

Its focus … That would be nice, but it is a bit out of reach.

Its axis, in fact its axis of symmetry … Fold it in half? But how?

Try the method of part one, with the ellipse. (previous post)

This looks promising. I even get another axis, for my coordinate system, if I really want the equation.

Now, analysis of the standard equation for a parabola (see later) says that a line at 45 deg to the axis, as shown, cuts the parabola at a point four focal lengths from the axis. In the picture, marked on the “vertical”axis, this is the length DH

So I need a point one quarter of the way from D to H. Easy !

and then the circle center D, with radius DH/4 cuts the axis of the parabola at the focal point (the focus).

Even better, we get the directrix as well …

and now for the mathy bit (well, you do the algebra, I did the picture)

Yes, I know that this one points up and the previous one pointed to the right !

All diagrams were created with my geometrical construction program, GEOSTRUCT

You will find it here:

www.mathcomesalive.com/geostruct/geostructforbrowser1.html

Filed under conic sections, construction, engineering, geometry, teaching

## Real problems with conic sections (ellipse, parabola)

So there is an oval hole in a metal casting. It’s supposed to be an elliptical hole. Is it ????? How can we find out ?????
A good start would be to find the line which would be the major axis if it was elliptical. This turns out to be an engineering problem, not a mathematical one (I cannot see a way!). If the oval curve has an axis of symmetry then the method below will find it:

Firstly, get a computer picture of the oval.
Take two circles, of different radii, and push them along until each one touches the oval in two places.

The line joining the two centers will be the axis of symmetry if there is one (this can be shown mathematically).

The equation of an ellipse uses the lengths of the major and minor axes. Do it !

The closeness to elliptic can be assessed in various ways. Think of one.

next…..finding the focus of a parabolic shape

Filed under conic sections, conics, education, engineering, geometry, math, teaching

## An almost real life geometry problem

I needed to move a point around a circle, in a computer graphics application, using the mouse pointer. It is clearly not sensible to have mouse pointer on the point all the time, so the problem was

“For a point anywhere, where is the point both on the circle and on the radial line?”

It may help to see the situation without the coordinate grid on show:

This is a problem with many ways to a solution, some of them incredibly messy !

Filed under engineering, geometry, math, teaching

## SOLVING PROBLEMS THAT MATTER, NOT GIVING ANSWERS

This is not mine, so read the original at