[From Bill Powers (970226.0945 MST)]

Can't leave well enough alone. At the end is the source code for a new

version of the PCT theodolite program; I'll also substitute the .exe file on

my ftp page, at www.frontier.net/users/powers_w. Just replace the old copies

with the new ones.

This version allows you to change dt from the keyboard. The upper limit is

set so you will see some computational oscillations when dt gets too large

for the gains of the control loops. There is no lower limit, so if you've

got an hour or two you can set dt to 1e-10 if you want. The + and - keys

double and halve the value of dt. Striking any key cuts the run short, so

you can raise and lower dt quickly. Also, the legends now include the

instructions (!) -- q to quit, space for another run, any key during the run

cuts it short.

Best,

Bill P.

## ···

====================================================================

program theodPCT;

{corrected for long integer in "for" loop to allow shorter dt}

{instructions added to screen}

{allows changing dt from keyboard}

uses dos, crt, graph, grutils;

var j,u,x,rx,v,rv,d,t,dt,dv,dx,gv,gx,r1,r2,d1,d2: real;

i: longint;

k,xscale: integer;

ch: char;

numstr: string;

begin

initgraphics;

dt := 0.001; {time resolution}

r1 := 3.0; {start of ref signal change}

r2 := 6.0; {end of ref signal change}

gv := 500.0; {gain of position control loop}

gx := 40.0; {gain of velocity control loop}

randomize;

repeat

rv := 0.0; {initialize velocity ref level & other variables}

x := 0.0;

v := 0.0;

t := 0.0;

d1 := 1.5 + 6.0*random; {start of disturbance}

d2 := d1 + 2.0; {end of disturbance}

xscale := round(9.0/dt/450.0); {x scaling for plots}

clearviewport;

setcolor(white);

outtextxy(0,0,'Pointing angle');

setcolor(lightcyan);

outtextxy(0,15,'Output torque');

setcolor(lightred);

outtextxy(0,30,'Disturbance');

setcolor(white);

str(dt:8:6,numstr);

outtextxy(0,45,'dt = ' + numstr + ' sec');

outtextxy(0,240 - 85,'1 radian');

outtextxy(0,240 - 104,'---------->');

setcolor(lightcyan);

outtextxy(0,240 - 124,'400 newton-meters');

setcolor(white);

outtextxy(0,420,'q to quit, space to repeat, any key to cut short');

outtextxy(0,435,'+,- to double or halve dt');

{ITERATION LOOP FOR RUN}

for i := 1 to round(9.0/dt) do

begin

if keypressed then i := round(9.0/dt);

{CHECK FOR CHANGES IN REF SIGNAL, DISTURANCE}

if (t > r1) and (t < r2)

then rx := 1 else rx := 0.0; {POSITION REF SIGNAL, RADIANS}

if (t > d1) and (t < d2) then

d := 100 else d := 0; {DISTURBANCE, NEWTON-METERS}

{CONTROL SYSTEM EQUATIONS}

rv := gx* (rx - x); {velocity ref level = output of position control}

u := gv*(rv - v); {output force = output of velocity control}

{LIMITS ON OUTPUT}

if u > 750.0 then u := 750.0 else {limit output to +/- 500 n-m}

if u < -750.0 then u := -750.0;

{ENVIRONMENT EQUATIONS}

v := v + (u + d)*dt; {integrate torque to get angular velocity}

x := x + v*dt; {integrate angular velocity to get angle}

{PLOT VARIABLES}

putpixel(i div xscale, 240,white);

putpixel(i div xscale, 240 - round(x*100),white);

putpixel(i div xscale, 240 - round(d/4),lightred);

putpixel(i div xscale, 240 - round(u/4),lightcyan);

{ADVANCE TIME BY DT}

t := t + dt;

end; {END OF RUN LOOP}

ch := readkey; {pause to view}

if ch in ['+','='] then dt := 2.0*dt;

if ch in ['-','_'] then dt := 0.5*dt;

if dt > 0.004 then dt := 0.004;

until ch = 'q'; {repeat on any key but 'q'}

closegraph;

end.