[From Bill Powers (970324.0100 MST)}
OK, after a few hours' sleep let's see if we can approach the problem
systematicallyl.
[Hans Blom, 970324]
If you want to be obtuse and go on arguing about the Ding an Sich, then I'm
going to call a halt to this discussion on the grounds that you're not
serious about it. You, the engineer, can analyze the physical theodolite,
and you, the engineer, can analyze the controller. If you can't see any
difference between them, we have no basis for talking to each other.
Assuming that you're interested, I will go through your program.
program theodolite_model_1;
const
J = 1.0; {theodolite characteristic}
var
x, {perceived position}
v, {perceived velocity}
u, {controller's output}
r, {reference level}
d, {disturbance estimate}
p, {prediction of next observation}
dt, {time increment; var so it can be changed}
t: real; {time; for display/loop control purposes}
This list of variables is unacceptable: the controller's internal model of
the physical theodolite is not shown separately, as it must be. What we need
are the following:
x {actual position}
xmod {modeled position}
v {actual velocity}
vmod {modeled velocity}
d {actual disturbance}
dmod {modeled disturbance}
I am substituting "dmod" for what you call "d", and "d" for what you call
"disturbance," for uniformity of notation.
I have agreed that the controller's internal model is to be accurate, but
this doesn't mean it behaves the same as the real theodolite. The FORM of
the model is accurate, and the PARAMETERS are accurate, but the internal
model is driven by u plus the MODELED disturbance, while the theodolite is
driven by u plus the ACTUAL disturbance. If the modeled disturbance differs
from the actual disturbance, there may be a difference (however momentary)
between x and xmod, and between v and vmod. That difference will propagate
to the next iteration, with results to be discovered.
When you do calculations inside the controller, the values of position and
velocity you use MUST be xmod and vmod, not x and v. By using x and v, you
are substituting the actual values of these variables for the modeled
values. This prevents any differences from developing between the model and
the actual theodolite. This strategy allows using the assumed perfection of
the model as a premise in proving that the model is perfect.
The following procedure is acceptable, using x and v for the actual values
of position and velocity:
procedure observe;
{this procedure implements the theodolite motion equations; it
generates the NEXT observation(s) x(t+dt) and v(t+dt), given u(t)
and d(t); thus it also shifts the time by dt; note that aa, vv and
xx are LOCAL variables, not visible outside this function; declaring
them as "const" saves their value from run to run; aa, vv and xx are
initially known to be zero}
function disturbance: real;
{this function implements the disturbance, which is known only
to procedure observe, the theodolite's motion simulator}
begin
if (t < 4.0) or (t > 5.0) then
disturbance := 0.0
else
disturbance := 10.0
end;
const
aa: real = 0.0;
vv: real = 0.0;
xx: real = 0.0;
begin
aa := (u + disturbance) / J;
xx := xx + vv*dt + 0.5*aa*dt*dt;
vv := vv + aa*dt;
t := t + dt;
x := xx; v := vv {the controller may use these}
end;
The following procedure is not acceptable: it uses the values of x and v
where it should use xmod and vmod. Also the value of the reference signal it
taken from the future, at time t+dt, which has not occurred yet:
procedure control;
{this function implements the controller; it knows the motion equa-
tions and J, is allowed to observe x(t) and v(t), and delivers u(t)}
begin
u := 2.0 * J * (reference (t+dt) - x - v*dt)/(dt*dt) - d;
end;
The body of this procedure should read, in my notation,
u := 2.0 * J * (reference (t) - xmod - vmod*dt)/(dt*dt) - dmod;
Now for the main loop:
repeat
writeln (t:12:3, r:12:6, x:12:6, v:12:6, d:12:6, u:12:6);
The following line is superfluous, since you use the "reference" function
inside the next line, the "control" function. Its main purpose here is to
shift the next value of the reference signal down one line in the printout,
to make it line up with the values of x: This is what creates the false
appearance of x matching the reference signal perfectly.
r := reference (t+dt); {specify desired x}
The controller, next, must use the corrected values given above:
control; {controller computes u}
Now comes the prediction step: the controller's internal model is used to
predict the next value, not of x but of xmod:
{predict next observation}
p := x + v*dt + 0.5 * (u + d) * (dt*dt) / J;
This line should read
xmod := xmod + vmod*dt + 0.5*(u+dmod)*(dt*dt)/J
The next statement computes the actual position and velocity based on the
same u as in the "prediction" statement, but using the actual disturbance.
Note that time does not shift by dt where your comment is, but AFTER the
computations in "observe" have been done (check the "observe" function):
{time shifts by dt here}
observe; {observe x (and v?)}
The next statement is now incorrect, since it fails to take into account
possible differences between x and xmod, and v and vmod. And in my notation,
the "d" should be "dmod."
{update disturbance estimate}
d := d + 2.0 * J * (x-p) / (dt*dt)
To derive the correct form, you must use the equations for the predicted and
actual values of x(t+dt):
xmod(t+dt) = xmod(t) + vmod(t)*dt + 0.5*(u(t)+dmod(t))*dt*dt/J, and
x(t+dt) = x(t) + v(t)*dt + 0.5*(u(t)+ d(t))*dt*dt/J
Subtracting the first equation from the second and solving for dmod(t), we
get the difference between the predicted and actual values of x(t+1), which
is the error in d. Adding this difference to d gives us a corrected value of
d based on the _current_ iteration (the best we can do). This leads to the
expression I used in my program:
dmod :=
dmod + 2*J/dt/dt*(x-xmod-(oldx-oldxmod)-(oldv-oldvmod)*dt);
In your computation of d, a subtle trick occurs:
d := d + 2.0 * J * (x-p) / (dt*dt)
What is "x-p" as you have computed it? x comes from the "observe" function,
where (substituting for the temporary variables) we have
x := x + v*dt + 0.5*(u + disturbance)*dt*dt/J,
where "disturbance" means the true disturbance.
The expression for p as you wrote it is
p := x + v*dt + 0.5 * (u + d) * (dt*dt) / J;
where "d" is the modeled disturbance.
Note that x and v in this equation are the _actual_ values, and so are
identical to the same values in the other equation. If we now compute x - p,
the first two terms disappear, leaving
(x-p) = 0.5*dt*dt/J*(u + disturbance - u - d)
Substituting this into your
d := d + 2.0 * J * (x-p) / (dt*dt)
we get
d := d + disturbance - d, so
d := disturbance.
If course if the controller's xmod and vmod actually did equal x and v
respectively, then this would be true: the correction of the modeled
disturbance would be perfect. So your development uses the following logic:
1. Assert that xmod = x and vmod = v.
2. Compute the value of d that will make xmod = x (which you have already
caused to be true by substituting x for xmod).
If xmod is NOT equal to x, and vmod is NOT equal to v, then the correction
of d done above will be wrong. And of course they are not equal; if they
were, d would already be correct!
To sum up:
In your development, you assume from the start that the controller's model
and the theodolite will behave identically, simply because the same equation
is used. But this assumption overlooks the fact that in the controller, the
modeled disturbance can differ from the real one, for at least one iteration
after the real one changes. This assumption effectively forces position and
velocity in the model to match the real position and velocity when they do
not and cannot match. Thus the effects of this momentary mismatch are not
caught by your model; they do not propagate around the loop as they should.
The true effect of the momentary mismatch is that there is always a small
unbalance between the disturbance and the output that is supposed to cancel
it. This shows up in the velocity term, which instead of going to zero, goes
to a small nonzero value that causes a continual departure of position from
the reference position. This value gets smaller as dt gets smaller, but it
never becomes zero, and the disturbance always causes a continuing and
unlimited departure of x from the reference value.
To this we must add the fact that your controller (in either your form or
the corrected form) produces an output that alternates on every iteration
between extreme values, with an amplitude that grows exactly as dt
decreases. As dt decreases toward zero, the behavior of your model does not
approach some asymptote, but diverges toward infinity. This is clearly not a
well-behaved controller, if one can call it a controller at all.
···
---------------------------------
If there is one thing I am sure of in addition to death and taxes, it is
that you will reject my analysis and assert that yours is correct. I assure
you that if this is what happens, I will be finished with this discussion.
So I look forward to getting more sleep in the future.
Best,
Bill P.