Bugs

[From Bill Powers (971219.1014 MST)]

I've been sort of trailing along after Richard Kennaway on the bug project.
Recently he got a set of leg control systems working so that when the bug
model is run, the body settles down to a equilibrium position. This is just
a start, but a lot of problems had to be solved to get that far.

I have just got my own cruder model to that stage. My bug drops to the
surface from an initial height; the legs flex when the feet touch the
ground, bringing the body to a stop with a single overshoot (which is the
quickest possible settling time). One advantage of my version is that it's
programmed in Pascal and runs a heck of a lot faster than Richard's Java
version does (on my machine). I get about 30 frames per second from my
program, and between 1.5 and 2.5 frames per second with the Java program,
depending on whether I'm moving the mouse.

The next step for both of us will be to start adding control systems for
x,y, and z position, plus pitch, roll and yaw, all working at the same
time. The bug still won't walk but it will rise and fall, rock and roll,
and lean this way and that.

One of the problems I've had is what happens just as the feet touch the
ground. The transition from controlling the leg positions to controlling
the body position leads to problems like oscillating back and forth between
contact and no contact, and other stuff. While contemplating these
problems, I think I've learned some things about bugs. This is called
second-guessing evolution.

One thing is that the sticky feet of bugs are very handy for the modeler. I
can make the model so it requires a definite upward effort to unstick the
feet from the surface; this immediately solves the problem at the boundary
between contact and no contact.

That led me to think of the problem of walking with sticky feet. If the bug
has three feet on the ground at a given time, then at the end of a tripod
step, three feet have to come down to the ground, and press hard enough
while muscles lift the three feet already in contact and unstick them from
the ground. We can even estimate about what the stickiness must be.

First, we assume that the bug is designed to be able to walk upside down,
underneath a surface. This means that the stickiness of three feet has to
be enough to support the weight of the bug -- plus, presumably, some margin
to accommodate accelerations and wind forces and so on. If we say that
twice the bug's weight has to be supported by three feet, then it should
take a force of at least 2/3 of the weight to unstick one foot.

Unless the bug has some way of varying the stickiness of its feet (I assume
it does not), the same force is required to unstick a foot in the normal
foot-down position of walking. Thus the minimum downward force that a leg
must be able to exert is 1/3 the weight of the bug plus the force needed to
unstick a foot, or 2/3 of the weight, so one leg must be able to exert at
least enough force to support the whole weight -- plus however much more is
needed for vertical accelerations.

These thoughts reminded me of the problem of the fly jumping to escape from
a swatter. I think Bruce Abbott mentioned that the initial move is an
upward jump propelled by the two middle legs (after which flying begins).
So why the jump, and why the middle two legs? The answer may be the
stickiness of the feet as just worked out above.

If the bug just jumped upward by pushing downward with all the feet, then 6
feet would have to become unstuck, overcoming a force of at least 6 times
the weight of the body plus whatever is needed for upward velocity after
the unsticking. But if four feet were first unstuck by being retracted, and
then only 2 feet were used for the jump, the total force to be overcome
would be only twice the weight of the body (the legs would thrust down
hard, and immediately be pulled up with enough force to unstick them). One
would expect the middle two legs to be used simply because the object is to
achieve an upward translation, not a tumble. And the middle two legs would
probably have stronger muscles than the other legs, if normally used this
way to take flight.

The reason for the initial jump is also clear. Without the jump, and with
the bug standing on all six feet, the wings would have to provide enough
lift to unstick all six feet, and that would require a lift of about six
times the weight of the bug. This could easily be more force that the
aerodynamic design can produce while hovering (with no forward motion).
With the jump, the flight mechanism needs to exert a lift only equal to the
weight of the bug, plus a little. So I conclude that the jump is simply
part of the normal takeoff maneuver, not just for emergency situations. Is
this borne out by observation?

I presume that all bugs with sticky feet would have this problem in taking
flight. The worst case would be taking off from a normal foot-down position
with all six legs on the ground.

I know that this is not nearly as interesting a subject as dealing with
suicides or reconciling PCT with other theories, but it's a real PCT
problem and as they say, God (or Darwin) is in the details.

Best,

Bill P.

[From Bruce Abbott (971219.1610 EST)]

Bill Powers (971219.1014 MST) --

The reason for the initial jump is also clear. Without the jump, and with
the bug standing on all six feet, the wings would have to provide enough
lift to unstick all six feet, and that would require a lift of about six
times the weight of the bug. This could easily be more force that the
aerodynamic design can produce while hovering (with no forward motion).
With the jump, the flight mechanism needs to exert a lift only equal to the
weight of the bug, plus a little. So I conclude that the jump is simply
part of the normal takeoff maneuver, not just for emergency situations. Is
this borne out by observation?

Yes. The same muscles that extend the middle legs (thereby popping the bug
off whatever surface it is on) also pull on the "pan" to which the flight
muscles attach internally. This starts the flight-motor oscillator.
However, it is also the case that you can start and stop the flight-motor
simply by removing or replacing something the fly can "hold onto" via its
feet (such as a small ball of styrofoam). Whether pulling the surface away
from the feet triggers the same rapid contraction of the middle leg muscles
as when the fly jumps, I do not know.

I don't know by what means the bug's feet adhere to a surface, or whether
this "stickiness" is under the bug's control.

This is a nice project. Any chance I could see the code to date?

Regards,

Bruce

Hi, Bruce --

This is a nice project. Any chance I could see the code to date?

Sure. Want me to put you on the bug list?

Relevant source and units attached.

Billprogram bug7;
{Now using vectors unit}
{leg movements and body rotations OK as of bug6.
Leg forces on body now correct.}

uses dos, crt, graph, vectors, grutils, mouse;

const bw = 20.0; {body half-width, half-length,mm}
      bl = 80.0;
      leg1 = 50.0; {leg segment lengths,mm}
      leg2 = 50.0;
      lw1 = 4; {leg widths,mm}
      lw2 = 2;
      Jx = 0.1; {inertial constants}
      Jy = 0.2;
      Jz = 0.2;
      Mass = 0.2; {gm}
      dt = 0.02; {iteration time}
      pause = true;
      nopause = false;
      vzero: triple = (0.0,0.0,0.0);
      left = 203; right = 205; up = 200; down = 208; PgUp = 201; PgDn = 209;
      Ins = 210; Del = 211; EndKey = 207; Home = 199;
      f1 = 187; f2 = 188; f3 = 189; f4 = 190; f5 = 191;
      f6 = 192; f7 = 193; f8 = 194; f9 = 195; f10 = 196;

type
     tripleptr = ^triple;

{ polytype = array[1..40] of triple;

     poly2dtype = array[1..40] of record
                                    x2d,y2d: integer
                                   end;}

     legtype = array[1..6] of
           record
            upper,lower, {upper and lower leg}
            fr,fp,fy,fc:vector; {reaction and contact forces}
            M,F,ft: triple; {foot x y z}
            cont: triple; {original contact x y z}
            a1,sa1,ca1, {angle adjacent shoulder, functions}
            ftp,sfp,cfp, {angle to foot, funtions}
            sp,sy,kp, {shoulder pitch, yaw; knee pitch}
            rsp,rsy,rkp, {reference levels for joint angles}
            csp,ssp,csy,ssy,ckp,skp, {sines, cosines}
            t1,t2,t3,l1,l2: real;
            l3: vector; {lengths, leg segs}
            contact: boolean;
           end;

      bodytype = record
                  front,right,
                  rear,left: vector; {body outline}
                  bc: triple; {xyz location of body center}
                  ba: triple; {body angles: roll, pitch, yaw}
                  bx,by,bz: vector; {body basis directions}
                  cbr,sbr,cbp,sbp,cby,sby: real; {sines and cosines}
                  sdr,cdr,sdp,cdp,sdy,cdy: real; {sine, cosine of deltas}
                  aa,la: point; {angular, linear acceleration}
                  av,lv: point; {angular, linear velocity}
                  da,dl: point; {angular, linear deltas}
                 end;

var mleg,leg: legtype;
    mbody,body: bodytype;
    msx,msy: real;
    tfx,tfy,tfz,tmx,tmy,tmz: real; {total forces and moments}
    g,gv3: vector;
    i: integer;
    ch: char;
    first: boolean;
    iters: longint;
    numstr: string;
    l3sq: real;
    gain: real;
    image: array[1..34] of vector;
    xaxis,yaxis,zaxis: vector;

{Get character as byte; add 128 if function key etc.}
function getcod: byte;
var ch: byte;
begin
ch := ord(readkey);
if (ch = 0) and keypressed then
  ch := ord(readkey) or $80;
getcod := ch;
end;

{show a real number on the screen with label. For debugging}
procedure outreal(x,y: integer; s: string; r: real );
begin
str(r:7:3,numstr);
setfillstyle(0,0);
bar(x,y,x+8*(length(s+numstr)+2),y+10);
outtextxy(x,y,s + numstr);
end;

procedure makelegs; {load master leg record}
var i,j: integer;
begin
for i := 1 to 3 do
with mleg[i] do
begin
  with upper do
  begin
   mag := leg1;
   org[0] := 2*(bl/3.0)*(2 -i); {x}
   org[1] := 0.0; {y}
   org[2] := bw; {z}
   sp := 0.1;
   sy := 0.0;
   dc[0] := 0.0;
   dc[1] := sin(sp);
   dc[2] := cos(sp);
  end;
  with lower do
  begin
   mag := leg2;
   for j := 0 to 2 do
    org[j] := upper.org[j] + upper.mag*upper.dc[j];
   kp := -1.2;
   dc[0] := 0.0;
   dc[1] := sin(sp+kp);
   dc[2] := cos(sp+kp);
  end;
end;
for i := 4 to 6 do
with mleg[i] do
begin
  with upper do
  begin
   mag := leg1;
   org[0] := 2*(bl/3.0)*(i - 5); {x}
   org[1] := 0.0; {y}
   org[2] := -bw; {z}
   dc[0] := 0.0;
   dc[1] := sin(sp);
   dc[2] := -cos(sp);
  end;
  with lower do
  begin
   mag := leg2;
   for j := 0 to 2 do
    org[j] := upper.org[j] + upper.mag*upper.dc[j];
   kp := 1.2;
   dc[0] := 0.0;
   dc[1] := cos(sp+kp);
   dc[2] := -sin(sp+kp);
  end;
end;
for i := 1 to 6 do {locate foot}
with mleg[i] do
begin
   for j := 0 to 2 do
  ft[j] := lower.org[j] + lower.mag*lower.dc[j];
end;
end;
{load corners and center of master body arrays, set body angles}
procedure makebody;
var i: integer;
begin
with mbody do
begin
with right do
begin
  org[0] := bl; org[1] := 0.0; org[2] := bw;
  dc[0] := -1.0; dc[1] := 0.0; dc[2] := 0.0;
  mag := 2.0*bl;
end;
with rear do
begin
  org[0] := -bl; org[1] := 0.0; org[2] := bw;
  dc[0] := 0.0; dc[1] := 0.0; dc[2] := -1.0;
  mag := 2.0 * bw;
end;
with left do
begin
  org[0] := -bl; org[1] := 0.0; org[2] := -bw;
  dc[0] := 1.0; dc[1] := 0.0; dc[2] := 0.0;
  mag := 2.0 * bl;
end;
with front do
begin
  org[0] := bl; org[1] := 0.0; org[2] := -bw;
  dc[0] := 0.0; dc[1] := 0.0; dc[2] := 1.0;
  mag := 2.0 * bw;
end;
  bc[0] := 0.0; bc[1] := 150.0; bc[2] := 0.0; {body center}
  {body angles in order of roll, pitch, yaw}
  ba[0] := 0.0; cbr := 1.0; sbr := 0.0;
  ba[1] := 0.0; cbp := 1.0; sbp := 1.0;
  ba[2] := 0.0; cby := 1.0; sby := 0.0;

end;
end;

procedure makegrid;
const vd: real = 480;
      eh: real = 480;
var i,j: integer;
    x,y,z: real;
    xa,xb,ya,yb: integer;
begin
setcolor(darkgray);
for j := -50 to 50 do
for i := 0 to 1 do
     begin
      z := vd + 2000*i;
      y := 0.0;
      x := 30.0*j;
      xa := hcenter + round(vd*x/z);
      ya := vsize - round(((z - vd)*eh + vd*y)/z);
     if (j = -50) or (i = 0) then moveto(xa,ya) else lineto(xa,ya);
    end;
for i := 0 to 50 do
for j := -1 to 1 do
     begin
      x := 2000*j;
      y := 0.0;
      z := vd + 40.0*i;
      xa := hcenter + round(vd*x/z);
      ya := vsize - round(((z - vd)*eh + vd*y)/z);
     if (j = -1) or (i = 0) then moveto(xa,ya) else lineto(xa,ya);
    end;
setcolor(white);
end;

procedure showbug;
var i,j: integer;
begin
  for i := 1 to 16 do
   showvector(image[i]);
end;

procedure initialize;
var i: integer;
begin
initgraphics;
graphmode := getmaxmode;
setgraphmode(graphmode);
hsize := getmaxx + 1;
vsize := getmaxy + 1;
hcenter := hsize div 2;
vcenter := vsize div 2;
gain := 20000.0;
makebody;
makelegs;
first := true;
mousex := 0;
mousey := 0;
with g do
begin
  mag := 980; dc[0] := 0; dc[1] := -1.0; dc[2] := 0;
end;
for i := 1 to 6 do
with leg[i] do
  begin
   sp := 0.1;
   kp := -1.0;
   sy := 0.0;
  end;
xaxis.dc[0] := 1.0; xaxis.dc[1] := 0.0; xaxis.dc[2] := 0.0;
yaxis.dc[0] := 0.0; yaxis.dc[1] := 1.0; yaxis.dc[2] := 0.0;
zaxis.dc[0] := 0.0; zaxis.dc[1] := 0.0; zaxis.dc[2] := 1.0;
end;

function invtan(x,y: real): real;
begin
if (abs(x) = 0.0) and (abs(y) = 0.0) then invtan := 0.0
else
begin
  if (abs(x) = 0.0) then
   if y > 0.0 then invtan := pi/2.0 else invtan := -pi/2.0
  else
   if abs(y) = 0.0 then invtan := 0.0
   else
    if abs(y) < abs(x) then
     invtan := arctan(y/x)
    else
     invtan := pi/2.0 - arctan(x/y);
end;
end;

{derotate object to standard position}
procedure derotate3( var v: triple; {input point or dir cosines}
                    var o: triple; {origin of rotation}
                    r,p,yaw: real); {angles}

var i: integer;
    x,y,z,
    ct,st: real;

begin
for i := 0 to 2 do v[i] := v[i] - o[i]; {remove offset}

ct := cos(-yaw); st := sin(-yaw);
x := ct*v[0] + st*v[2]; {yaw: rotate about y}
z := -st*v[0] + ct*v[2];
v[0] := x; v[2] := z;

ct:= cos(-p); st := sin(-p);
x := ct*v[0] + st*v[1]; {pitch: rotate about z}
y := -st*v[0] + ct*v[1];
v[0] := x; v[1] := y;

ct := cos(-r); st := sin(-r);
y := ct*v[1] + st*v[2]; {roll: rotate about x}
z := -st*v[1] + ct*v[2];
v[1] := y; v[2] := z;

for i := 0 to 2 do v[i] := v[i] + o[i]; {restore offset}
end;

{rotate from standard position into real position}
procedure rotate3( var v: triple; {input point or dir cosines}
                   var o: triple; {origin of rotation}
                    r,p,yaw: real); {angles}

var i: integer;
    x,y,z,
    ct,st: real;

begin
for i := 0 to 2 do v[i] := v[i] - o[i]; {remove offset}

ct := cos(r); st := sin(r);
y := ct*v[1] + st*v[2]; {roll: rotate about x}
z := -st*v[1] + ct*v[2];
v[1] := y; v[2] := z;

ct:= cos(p); st := sin(p);
x := ct*v[0] + st*v[1]; {pitch: rotate about z}
y := -st*v[0] + ct*v[1];
v[0] := x; v[1] := y;

ct := cos(yaw); st := sin(yaw);
x := ct*v[0] + st*v[2]; {yaw: rotate about y}
z := -st*v[0] + ct*v[2];
v[0] := x; v[2] := z;

for i := 0 to 2 do v[i] := v[i] + o[i]; {restore offset}
end;

{rotate origin about center, rotate direction cosines}
procedure rotatevector( var vin: vector;
                        var vout: vector;
                        var o: triple; {origin of rotation}
                        r,p,yaw: real);
begin
vout := vin;
with vout do
  begin
   rotate3(org,o,r,p,yaw);
   rotate3(dc,vzero,r,p,yaw);
  end;
end;

{translate image of body and legs; show on screen}
procedure translate;
var i,j: integer;
begin
for i := 1 to 16 do
for j := 0 to 2 do
begin
  image[i].org[j] := image[i].org[j] + body.bc[j];
  if i in [1..4] then
  if image[i].org[1] < 0.0 then image[i].org[1] := 0.0;
end;
end;

{rotate body and legs about zero in standard position}
procedure bodypos;
var i,j: integer;
begin
with body do
begin
  rotatevector(front,image[1],vzero,ba[0],ba[2],ba[1]);
  rotatevector(right,image[2],vzero,ba[0],ba[2],ba[1]);
  rotatevector(rear,image[3],vzero,ba[0],ba[2],ba[1]);
  rotatevector(left,image[4],vzero,ba[0],ba[2],ba[1]);
end;
for i := 1 to 6 do
  with leg[i] do
  with body do
  begin {center} {r} {yaw} {p}
   rotatevector(upper, image[3+2*i], vzero, ba[0],ba[2],ba[1]);
   rotatevector(lower, image[4+2*i], vzero, ba[0],ba[2],ba[1]);
  end;
end;

{COMPUTING MOTIONS AND NEW BODY POSITION}

{Given forces and moments acting at center of body, calculate
linear and angular velocities, delta xyz and delta angles, and
new body and leg positions}

procedure movebug;
var i,j: integer;
    d1x,d2x,d1y,d2y,d1z,d2z: real; {temporary linear deltas}

begin
with body do
begin
{angular acceleration = total moment.mag about axis / moment.mag of inertia}
{Jx, Jy, and Jz are global constants}

with aa do
begin
  x := tmx/Jx; {angular acceleration}
  y := tmy/Jy;
  z := tmz/Jz;
  x := 0.0;
end;

with av do
begin
  x := x + aa.x*dt; {angular velocity}
  y := y + aa.y*dt;
  z := z + aa.z*dt;
end;

{Linear accelerations = total Force.mag along axis / Mass}
{Mass is a global constant}
with la do
with gv3 do
begin
  x := tfx/Mass + mag*dc[0]; {linear acceleration}
  y := tfy/Mass + mag*dc[1];
  z := tfz/Mass + mag*dc[2];
end;
if body.bc[1] < 2 then
i := i;
with lv do
begin
  x := x + la.x*dt; {linear velocity}
  y := y + la.y*dt;
  z := z + la.z*dt;
end;

with da do
begin
  x := (av.x + 0.5*aa.x*dt)*dt; {delta angle about axis}
  y := (av.y + 0.5*aa.y*dt)*dt;
  z := (av.z + 0.5*aa.z*dt)*dt;
end;

with dl do
begin
  x := (lv.x + 0.5*la.x*dt)*dt; {delta position}
  y := (lv.y + 0.5*la.y*dt)*dt;
  z := (lv.z + 0.5*la.z*dt)*dt;
end;
{rotate linear deltas to final position}
rotate3(triple(dl),vzero,ba[0],ba[2],ba[1]);
bc[0] := bc[0] + dl.x;
bc[1] := bc[1] + dl.y;
if bc[1] < 0.0 then bc[1] := 0.0; {hit bottom}
bc[2] := bc[2] + dl.z;
ba[0] := ba[0] + da.x;
ba[1] := ba[1] + da.y;
ba[2] := ba[2] + da.z;
end;
end; { of movebug}

{Find sp, kp, sy, a1 for leg with foot fixed at "cont"}

procedure computeleg;
var i,j: integer;
    u,v,w: real;
begin
for i := 1 to 6 do
  with leg[i] do
  if contact then
   begin
    ft := cont; {set foot to contact position}
    ft[1] := cont[1] - 2.0; {set it 2 pixels below ground level}
    makevector(upper.org,ft,l3);
    with l3 do
    begin
     sfp := (ft[1] - org[1])/mag; {sine foot pitch}
     u := sqrt(sqr(ft[0] - org[0]) + sqr(ft[2] - org[2]));
    end;
    cfp := u/l3.mag;
    ftp := invtan(cfp,sfp); {foot pitch angle, from shoulder}

    ca1 := (sqr(lower.mag) - sqr(upper.mag) - sqr(l3.mag))
           /(-2.0*upper.mag*l3.mag);
    sa1 := sqrt(1.0 - ca1*ca1);
    a1 := invtan(ca1,sa1); {interior angle at shoulder}
    sp := a1 + ftp; {shoulder pitch angle}

    ckp := (sqr(l3.mag) - sqr(upper.mag) - sqr(lower.mag))
            /(2.0*upper.mag*lower.mag); {exterior angle}
    skp := sqrt(1.0 - ckp*ckp);
    kp := -invtan(ckp,skp); {knee pitch angle}
   end;
end;

{Given pitch and yaw angles, construct leg in correct
position relative to body in standard position.}
procedure legpos;
var t,u,v,w: real;
    i,j: integer;
begin
for i := 1 to 6 do
with leg[i] do
begin
  t := sin(sy);
  if i < 4 then u := cos(sy) else u := -cos(sy);
  v := sin(sp);
  w := cos(sp);
  with upper do
  begin
   dc[0] := w*t; {x}
   dc[1] := v; {y}
   dc[2] := w*u; {z}
  end;
   v := sin(sp+kp);
   w := cos(sp+kp);
  with lower do
  begin
   dc[0] := w*t ;
   dc[1] := v;
   dc[2] := w*u;
  end;
  for j := 0 to 2 do
  begin
   lower.org[j] := upper.org[j] + upper.mag*upper.dc[j];
   ft[j] := lower.org[j] + lower.mag*lower.dc[j];
  end;
  if not contact then cont := ft;
end;
end;

procedure checkcontact;
var i,j: integer;
    foot: triple;
begin
for i := 1 to 6 do
with image[4+2*i] do
begin
  for j := 0 to 2 do
   foot[j] := org[j] + mag*dc[j];
  if (not leg[i].contact) and (foot[1] < 0.0) then
  with leg[i] do
  with body do
  begin
   cont := foot;
   {Move contact point to proper place in standard position}
   derotate3(cont, bc, ba[0],ba[2],ba[1]);
  end;
  with leg[i] do
  contact :=
     (foot[1] < 0.0) or ((f[1] > 0.5*mass*g.mag));
if leg[i].contact then leg[i].cont[1] := - body.bc[1];
  {derotate gravity}
  gv3 := g;
  with body do {rotate gravity vector for standard position}
  derotate3(gv3.dc,vzero,ba[0],ba[2],ba[1]);
end;
end;

procedure muscleforce; {compute forces with legs in standard positions}
const olde1: array[1..6] of real = (0,0,0,0,0,0);
      olde2: array[1..6] of real = (0,0,0,0,0,0);
      olde3: array[1..6] of real = (0,0,0,0,0,0);
      e1: array[1..6] of real = (0,0,0,0,0,0);
      e2: array[1..6] of real = (0,0,0,0,0,0);
      e3: array[1..6] of real = (0,0,0,0,0,0);
      damp: real = 10.0;
var i: integer;
begin
tfx := 0.0; tmx := 0.0;
tfy := 0.0; tmy := 0.0;
tfz := 0.0; tmz := 0.0;
for i := 1 to 6 do
  with leg[i] do
  if contact then
  begin
   olde1[i] := e1[i];
   e1[i] := rsp - sp;
   T1 := gain*(e1[i] + damp*(e1[1] - olde1[i]));
   olde2[i] := e2[i];
   e2[i] := rkp - kp;
   T2 := gain*(e2[i] + damp*(e2[1] - olde2[i]));
   olde3[i] := e3[i];
   e3[i] := rsy - sy;
   T3 := gain*(e3[i] + damp*(e3[1] - olde3[i]));
   makevector(upper.org,ft,l3);
   Fr := l3;
   Fr.mag := (T2 - T1)/l3.mag; {Radial force vector}
   Fr.org := ft;
   perp(lower,upper,Fy);
   Fy.mag := t3/l3.mag;
   Fy.org := ft;
   perp(Fr,Fy,Fp);
   Fp.mag := T1/l3.mag;
   Fp.org := ft;
   {compute reaction forces in x,y,z}
   F[0] := fp.mag*cosangle(Fp,xaxis) +
           fr.mag*cosangle(Fr,xaxis) +
           fy.mag*cosangle(Fy,xaxis);
   F[1] := fp.mag*cosangle(Fp,yaxis) +
           fr.mag*cosangle(Fr,yaxis) +
           fy.mag*cosangle(Fy,yaxis);
   F[2] := fp.mag*cosangle(Fp,zaxis) +
           fr.mag*cosangle(Fr,zaxis) +
           fy.mag*cosangle(Fy,zaxis);

   {Compute reaction moments rotated to match leg positions}
   M[0] := 0.1*(t3*sfp + t1)*csy; {dyne-cm}
   M[1] := -0.1*t3*cfp;
   M[2] := 0.1*(t3*sfp - t1)*ssy;

   tfx := tfx + F[0]; {total linear forces}
   tfy := tfy + F[1];
   tfz := tfz + F[2];

   with body do
   with upper do {total moments acting on body}
    begin
     tmx := tmx + M[0] + org[2]*f[1] - org[1]*F[2];
     tmy := tmy + M[1] + org[2]*F[0] - org[0]*F[2];
     tmz := tmz + M[2] + org[0]*F[1] - org[1]*F[0];
    end;
  end
else
  begin
   Fp.mag := 0.0; Fr.mag := 0.0; Fy.mag := 0.0;
   sp := rsp;
   kp := rkp;
   sy := rsy;
  end;
end;

begin
initialize;
i := 1;
setcolor(white);
setwritemode(XORput);
iters := 0;
body := mbody;
leg := mleg; {body and legs to master positions}
outtextxy(100,200,'PRESS KEY TO START');
ch := readkey;
clearviewport;
makegrid;
outtextxy(0,470,'Type ''q'' to quit');
mousey := 0;
mousex := 0;
for i := 1 to 6 do
with leg[i] do
begin
  rsp := 0.1;
  rkp := -1.2;
  rsy := 0.0;
   sp := rsp;
   kp := rkp;
   sy := rsy;
  contact := false;
end;
  body.ba[0] := 0.0;
  body.ba[2] := 0.0;
  body.ba[1] := 0.5;
  body.bc[1] := 100.0;
repeat
  readmouse;
  msy := mousey/10000;
  for i := 1 to 6 do
  with leg[i] do
  begin
   rsp := 0.1+ msy;
   rkp := -1.2 + 2.0*msy;
  end;
  computeleg; {if contact, compute sp,sy,kp,a1,l3,ft}
  legpos; {leg & foot position relative to body}
  muscleforce; {compute forces fp,fy,fr}
  showbug; {erase previous image}
  movebug;
  bodypos; {rotate body and leg images to real position}
  translate; {move according to body center coordinates}
  showbug;
  delay(20);
  checkcontact;
  if keypressed then ch := readkey;
  inc(iters);
  until ch = 'q';
  closegraph;
  writeln(iters:10,' iterations');
  ch := readkey;
end.unit vectors;

{vector operations. Right-hand coordinate system}

Interface

type triple = array[0..2] of real;
     point = record x,y,z: real end;
     vector = record
               org, {origin of vector}
               dc: triple; {direction cosines, x y z}
               mag: real; {magnitude}
              end;

procedure showvector(v:vector);
procedure makevector(var a,b: triple; var v: vector);
function cosangle(var a,b: vector): real;
procedure dircos(var a,b,c: triple);
procedure proj(var a,b,c: vector);
procedure perp(var a,b,c: vector);
procedure resultant(var a,b,c:vector);

implementation
uses dos,graph,grutils;

procedure showvector(v:vector);

var x1,y1,z1,x2,y2,z2,eh,vd: real;
    xa,ya,xb,yb,i: integer;
begin
eh := vsize;
vd := vsize;
with v do
begin
  x1 := org[0]; y1 := org[1]; z1 := org[2] + vd + 200;
  x2 := x1 + mag*dc[0];
  y2 := y1 + mag*dc[1];
  z2 := z1 + mag*dc[2];
end;
xa := hcenter + round(vd*x1/z1);
ya := vsize - round(((z1 - vd)*eh + vd*y1)/z1);
xb := hcenter + round(vd*x2/z2);
yb := vsize - round(((z2 - vd)*eh + vd*y2)/z2);
line(xa,ya,xb,yb);
end;

{given points A and B, compute origin, magnitude, and direction
cosines for vector V}
procedure makevector(var a,b: triple; var v: vector);
var i: integer;
begin
with v do
  begin
   mag := sqrt(sqr(b[0]-a[0]) + sqr(b[1]-a[1])+sqr(b[2]-a[2]));
   org := a;
   if mag > 0.0 then
   for i := 0 to 2 do
    dc[i] := (b[i] - a[i])/mag;
  end
end;

{Calculate direction cosines C given two points A and B. Sense is
in direction of second point. Cosines are in order of x,y,z.}
procedure dircos(var a,b,c: triple);
var i: integer;
    var u: real;
begin
u := 0.0;
for i := 0 to 2 do
  u := u + sqr(b[i] - a[i]);
u := sqrt(u);
for i := 0 to 2 do
  c[i] := (b[i] - a[i])/u;
end;

{Cosine of angle between vectors A and B, 1 to -1}
function cosangle(var a,b: vector): real;
var i: integer;
    c: real;
begin
   c := 0.0;
   for i := 0 to 2 do
   c := c + a.dc[i]*b.dc[i];
   cosangle := c;
end;

{project vector A onto direction of vector B; result, vector C.
Origin of C is origin of B}
procedure proj(var a,b,c: vector);
begin
c.mag := a.mag* cosangle(a,b);
c.dc := b.dc;
c.org := b.org;
end;

{construct direction cosines of vector C at right angles to
vectors A, B. If A rotates into B, C is in direction according to
right-hand rule. Origin of C is origin of A. Magnitude of C is 1.0}
procedure perp(var a,b,c: vector);
var i,j,k: integer;
begin
for i := 0 to 2 do
begin
  j := (i+1) mod 3;
  k := (i+2) mod 3;
  c.dc[i] := a.dc[j]*b.dc[k] - a.dc[k]*b.dc[j];
end;
  c.org := a.org;
  c.mag := 1.0;
end;

{ Add vectors. Origin of resultant is origin of first vector}
procedure resultant(var a,b,c:vector);
var i: integer;
begin
c.mag := sqrt(sqr(a.mag) + sqr(b.mag) + 2.0*a.mag*b.mag*cosangle(a,b));
for i := 0 to 2 do
  c.dc[i] := (a.mag * a.dc[i] + b.mag * b.dc[i])/c.mag;
c.org := a.org;
end;

begin
end.

unit GrUtils;
{ Graphics Utilities Unit }

interface

uses
  Graph, bgidriv;

var
  GraphDriver, GraphMode, Error: integer;
  hsize, hcenter, vsize, vcenter: integer;

procedure InitGraphics;
procedure Retrace;

implementation

const BGIDIR = '\tp\bgi';

procedure retrace;
begin
case Graphdriver of
  Ega,Vga,Ega64,EgaMono: begin
           while (port[$3da] and 8) = 8 do ;
           while (port[$3da] and 8) = 0 do ;
          end;
  HercMono: begin
           while (port[$3ba] and $80) = 0 do ;
           while (port[$3ba] and $80) = $80 do ;
           end;
  ELSE begin
          while (port[$3da] and 8) = 0 do ;
          while (port[$3da] and 8) = 8 do ;
       end;
end;
end;

procedure Abort(Msg : string);
begin
  Writeln(Msg, ': ', GraphErrorMsg(GraphResult));
  Halt(1);
end;

procedure InitGraphics; {ADAPTS TO HARDWARE}
begin
  { Register all the drivers }
  if RegisterBGIdriver(@CGADriverProc) < 0 then
    Abort('CGA');
  if RegisterBGIdriver(@EGAVGADriverProc) < 0 then
    Abort('EGA/VGA');

  GraphDriver := Detect; { autodetect the hardware }
  InitGraph(GraphDriver, GraphMode, ''); { activate graphics }
  if GraphResult <> grOk then { any errors? }
  begin
    Writeln('Graphics init error: ', GraphErrorMsg(GraphDriver));
    Halt(1);
  end;
  vsize := getmaxx; hsize := getmaxy;
  vcenter := (vsize + 1) div 2;
  hcenter := (hsize + 1) div 2;
end;

begin
end.

[From Bill Powers (971219.1955 MST)]

From: David Goldstein
Subject: Re: Bugs
Date: 12/19/97

Bill,

This is what real PCT is about? How do you think PCT will be advanced
by this?

If so: How many bugs can stand on the head of a pin, with our without
sticky feet?

Seriously, what is there about this problem which attracts you to it,
as a theoretician? I assume that your choice of this problem is a means
by which you are controlling for something, just can't figure out what?
Did the discussion on OCD, suicide, Bipolar Disorder, etc., really
drive you buggy? and hence we have the treatise on bugs?

There are several reasons for developing the bug model. One is that a
number of people have proposed bug models, but have used such sketchy and
simple physical models that the underlying theory is put to no real test.
When this model is done, I (we) will publish the basic physical part as a
test bed for all theories -- of course presenting our own while we're at it.

The second reason is that this was a second choice for me. What I really
want is a physical model of a human being (or a cat), with all the dynamics
in all the degrees of freedom properly expressed, and a good muscle model
to go with each DF. This would be another test bed for theories, one that
would challenge all open-loop, computed-output, and model-based theories --
as well as PCT, of course. Unfortunately, the mathematical modeling is
beyond me.

The bug model is much simpler, since we can neglect the mass of the legs
and such things as Coriolis forces, treating the body as single rigid
object, all without seriously misrepresenting the real behavior of the
physical system. I figure that if we can come up with a realistic bug model
that behaves the way real bugs do, this might interest people who do have
the mathematical skills required to help out with the model of the more
complex system. Richard Kennaway has started making noises about looking
for grant money to carry the project further.

To me, HPCT is to a complete model of behavior as an artist's pencil
sketches are to the finished oil painting. A truly enormous amount of
research has to be done to get beyond the sketch stage. By the time such
research gets to the levels that have some psychological pizzazz, I believe
that our very conception of the problems we're trying to solve will hae
changed fundamentally. What we guess now about how the higher levels work,
what people's problems really are, and what the solutions really are, will
look quaint and primitive in 100 years. So it behooves us, or at least some
of us, to pay some attention to laying firm foundations on which others can
build, rather than trying to jump into the most complex parts of the job
before anyone is prepared to handle them. We haven't carried psychology to
the level of 18th-century physics yet. Someone has to be working at the
foundations while others are trying to do something about life's problems
as they are understood right now.

No value judgement intended or implied. Life goes on; we have to deal with
it as best we can. But I don't. My sights are set on a target a long way off.

Best,

Bill P.