[python] Re: Bicycling physics
- From: Dirk Bonne <dirk_bonne@xxxxxxxxxx>
- To: python@xxxxxxxxxxxxx
- Date: Sat, 12 Mar 2005 14:19:42 +0100
|
Hi Ray, when I ran your program I got a little different output than the one you copied pasted in your mail. I had also to change the for loop in the gtrail (although it prints out trails going from 10 -> -25cm, the program uses -40 -> -48cm). With 10->-25 the values printed out get quite wierd (no self centering effect): steerAngle 20.0
pivot ground trail
10 5 0 -5 -10 -15 -20 -25
40 -4.077 -3.532 -2.98 -2.421 -1.855 -1.281 -0.699 -0.111
41 -4.109 -3.562 -3.007 -2.444 -1.873 -1.295 -0.708 -0.114
42 -4.139 -3.589 -3.031 -2.466 -1.892 -1.31 -0.719 -0.12
43 -4.165 -3.614 -3.055 -2.487 -1.91 -1.325 -0.731 -0.129
44 -4.189 -3.637 -3.077 -2.507 -1.929 -1.342 -0.745 -0.14
45 -4.209 -3.658 -3.097 -2.527 -1.947 -1.359 -0.761 -0.153
46 -4.226 -3.675 -3.115 -2.545 -1.966 -1.376 -0.777 -0.168
47 -4.24 -3.69 -3.131 -2.562 -1.983 -1.394 -0.795 -0.185
48 -4.249 -3.702 -3.145 -2.578 -2.0 -1.412 -0.813 -0.204
49 -4.255 -3.711 -3.156 -2.591 -2.016 -1.43 -0.833 -0.224
50 -4.257 -3.716 -3.165 -2.604 -2.031 -1.447 -0.852 -0.246
51 -4.254 -3.718 -3.171 -2.614 -2.045 -1.465 -0.873 -0.269
52 -4.248 -3.717 -3.175 -2.622 -2.057 -1.481 -0.893 -0.293
53 -4.236 -3.711 -3.175 -2.627 -2.068 -1.497 -0.914 -0.319
54 -4.221 -3.702 -3.172 -2.631 -2.077 -1.512 -0.935 -0.344
55 -4.2 -3.689 -3.166 -2.631 -2.085 -1.526 -0.955 -0.371
56 -4.175 -3.671 -3.156 -2.629 -2.09 -1.539 -0.975 -0.398
57 -4.145 -3.649 -3.142 -2.624 -2.093 -1.55 -0.994 -0.424
58 -4.109 -3.623 -3.125 -2.616 -2.094 -1.559 -1.012 -0.451
59 -4.069 -3.592 -3.104 -2.604 -2.092 -1.567 -1.029 -0.478
60 -4.023 -3.557 -3.079 -2.589 -2.087 -1.572 -1.045 -0.504
61 -3.972 -3.517 -3.05 -2.571 -2.08 -1.576 -1.059 -0.529
62 -3.916 -3.471 -3.016 -2.548 -2.069 -1.577 -1.072 -0.554
63 -3.854 -3.421 -2.978 -2.522 -2.055 -1.575 -1.083 -0.577
64 -3.786 -3.366 -2.935 -2.493 -2.038 -1.571 -1.092 -0.599
65 -3.713 -3.306 -2.888 -2.458 -2.017 -1.564 -1.098 -0.619
66 -3.634 -3.241 -2.836 -2.42 -1.993 -1.554 -1.102 -0.638
67 -3.55 -3.17 -2.779 -2.378 -1.965 -1.54 -1.103 -0.654
68 -3.46 -3.094 -2.717 -2.33 -1.932 -1.523 -1.102 -0.668
69 -3.364 -3.012 -2.651 -2.279 -1.896 -1.503 -1.097 -0.68
70 -3.262 -2.925 -2.579 -2.222 -1.856 -1.478 -1.089 -0.689
71 -3.154 -2.832 -2.502 -2.161 -1.811 -1.45 -1.078 -0.695
72 -3.04 -2.734 -2.419 -2.095 -1.762 -1.418 -1.063 -0.698
73 -2.921 -2.631 -2.332 -2.024 -1.708 -1.381 -1.045 -0.698
74 -2.795 -2.521 -2.239 -1.948 -1.649 -1.34 -1.022 -0.694
75 -2.664 -2.406 -2.141 -1.867 -1.586 -1.295 -0.995 -0.686
76 -2.527 -2.286 -2.037 -1.781 -1.517 -1.245 -0.964 -0.675
77 -2.384 -2.16 -1.928 -1.69 -1.444 -1.191 -0.929 -0.659
78 -2.235 -2.028 -1.814 -1.593 -1.366 -1.131 -0.889 -0.639
79 -2.081 -1.89 -1.694 -1.491 -1.282 -1.067 -0.845 -0.615
To be able to compare, I adapted my program (octave). Here is the
output of my program:octave> flevo2 aa 10.0 5.0 0.0 -5.0 -10.0 -15.0 -20.0 -25.0 40.0 -0.26 -0.12 0.01 0.15 0.29 0.42 0.56 0.69 41.0 -0.26 -0.12 0.01 0.15 0.29 0.42 0.56 0.70 42.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.56 0.70 43.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.57 0.70 44.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.57 0.71 45.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.57 0.71 46.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.57 0.71 47.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.57 0.71 48.0 -0.26 -0.12 0.01 0.15 0.29 0.43 0.56 0.70 49.0 -0.26 -0.12 0.01 0.15 0.29 0.42 0.56 0.70 50.0 -0.26 -0.12 0.01 0.15 0.29 0.42 0.56 0.70 51.0 -0.25 -0.12 0.01 0.15 0.28 0.42 0.56 0.69 52.0 -0.25 -0.12 0.01 0.15 0.28 0.42 0.55 0.69 53.0 -0.25 -0.12 0.01 0.15 0.28 0.41 0.55 0.68 54.0 -0.25 -0.12 0.01 0.14 0.28 0.41 0.54 0.67 55.0 -0.25 -0.12 0.01 0.14 0.27 0.40 0.53 0.67 56.0 -0.24 -0.12 0.01 0.14 0.27 0.40 0.53 0.66 57.0 -0.24 -0.11 0.01 0.14 0.26 0.39 0.52 0.65 58.0 -0.24 -0.11 0.01 0.14 0.26 0.39 0.51 0.64 59.0 -0.23 -0.11 0.01 0.13 0.26 0.38 0.50 0.63 60.0 -0.23 -0.11 0.01 0.13 0.25 0.37 0.49 0.62 61.0 -0.22 -0.11 0.01 0.13 0.25 0.36 0.48 0.60 62.0 -0.22 -0.10 0.01 0.12 0.24 0.36 0.47 0.59 63.0 -0.21 -0.10 0.01 0.12 0.23 0.35 0.46 0.58 64.0 -0.21 -0.10 0.01 0.12 0.23 0.34 0.45 0.56 65.0 -0.20 -0.10 0.01 0.11 0.22 0.33 0.44 0.54 66.0 -0.20 -0.09 0.01 0.11 0.21 0.32 0.42 0.53 67.0 -0.19 -0.09 0.01 0.11 0.21 0.31 0.41 0.51 68.0 -0.18 -0.09 0.01 0.10 0.20 0.30 0.40 0.49 69.0 -0.18 -0.09 0.01 0.10 0.19 0.29 0.38 0.48 70.0 -0.17 -0.08 0.01 0.09 0.18 0.27 0.37 0.46 71.0 -0.16 -0.08 0.01 0.09 0.18 0.26 0.35 0.44 72.0 -0.16 -0.08 0.00 0.09 0.17 0.25 0.33 0.42 73.0 -0.15 -0.07 0.00 0.08 0.16 0.24 0.32 0.40 74.0 -0.14 -0.07 0.00 0.08 0.15 0.23 0.30 0.38 75.0 -0.13 -0.07 0.00 0.07 0.14 0.21 0.28 0.36 76.0 -0.13 -0.06 0.00 0.07 0.13 0.20 0.27 0.33 77.0 -0.12 -0.06 0.00 0.06 0.13 0.19 0.25 0.31 78.0 -0.11 -0.05 0.00 0.06 0.12 0.17 0.23 0.29 79.0 -0.10 -0.05 0.00 0.05 0.11 0.16 0.21 0.27 I think the calculation of the frontAxle* values is too simplified. msssHT.py uses a linear formula for that, but that makes a rather large error. Instead of: k=steerAngle*2/pi I would propose k=1 - cos(steerAngle) This has considarable effect. For example at 20degrees k is only 0.0603 instead of 0.222 The second problem is in the steerDropAngle. The correct formula is rather complicated. I switch to a matrix calculation here. if you define the unrotated front wheel as [xt,yt,zt] = [R*cos(t), 0, R*sin(t)] where t in [0, 2*pi[ and R is the front wheel radius then the formula for the front wheel rotated in the new steering axle point would be (choke!): [cos(HA)*cos(HA)*R*cos(t) + cos(SA)*(1-cos(HA)*cos(HA))*R*cos(t) -sin(HA)*cos(HA)*R*sin(t) + cos(SA)*sin(HA)*cos(HA)*R*sin(t), -sin(SA)*sin(HA)*R*cos(t)- sin(SA)*cos(HA)*R*sin(t), -cos(HA)*sin(HA)*R*cos(t) + cos(SA)*cos(HA)*sin(HA)*R*cos(t) + sin(HA)*sin(HA)*R*sin(t) + cos(SA)*(1-sin(HA)*sin(HA))*R*sin(t)] where HA= headangle, SA is steerAngle This formula follows from substituting the formulas at: http://www.euclideanspace.com/maths/algebra/matrix/orthogonal/rotation/ The minimum z would be at (follows out of the derivation of the z component in parameter t tg(t) = ((sin(HA)*sin(HA) + cos(SA)*(1-sin(HA)*sin(HA))) / (cos(SA)*cos(HA)*sin(HA) - cos(HA)*sin(HA))) This t you need to put back into the largish formula above to find the Z/X value of the lowest point of the wheel. This values can then be used to calculate the steerDropAngle. As a test, I adapted your program (molested it actaully), but there are still a difference between the results. The new results of your program are here:
For lower pivot angles the results are now more or less the same as in
my program, but at larger pivot angles there is a diviation. But at
least the values that are output are sane now (negative trail give a
rise at steering angle).I haven't checked the reason for the deviation (and like i said it is not impossible that I am wrong. Calculating all this I find quite tricky). Attached is the adapted program. Tip: you might have a look at quaternions. They simplify working in three dimensional space a great deal. AFAIK, there is a quaternion package in the python programming language. cheers Dirk mtb@xxxxxxx wrote: Hi Dirk, |
#!/usr/bin/python
#Boa:PyApp:main
import math
modules ={}
def main():
wheelBase = 72.5
steerAngle = 20. * math.pi / 180.0
rearMassZ = 28.
rearMassX = 67.
frontWheelRadius = 33.
rearWheelRadius = 33.
## calc angle from raxle to rmass, from horizontal
rmass2RaxleAngle = math.asin((rearMassZ-rearWheelRadius)/rearMassX)
## calc dist from raxle to rmass
rmass2Raxle = ((rearWheelRadius-rearMassZ)**2+rearMassX**2)**.5
## angle front to rear axle
## a wheelie is a pos+ angle
angleBtwnAxles = math.atan((frontWheelRadius-rearWheelRadius)/wheelBase)
print 'steerAngle', steerAngle/math.pi*180.0
print 'pivot ground trail '
print ' ', '\t10', '\t5', '\t0', '\t-5', '\t-10', '\t-15', '\t-20',
'\t-25'
for angle in range(40, 71, 1):
headAngle = angle * math.pi / 180.0 # deg from horiz
print '\n',int(round(headAngle*180/math.pi)), ' ',
## iterate over the ground trail
for gtrail in range(10, -30, -5):
# print "gtrail=", gtrail, "\n";
mechTrail = frontWheelRadius*math.cos(headAngle) -
gtrail*math.sin(headAngle)
# print "mechTrail=", mechTrail, "\n";
## rear wheel does not move or roll, frame rotates around the axle
## turnAnglePnt is imaginary point the axle rotates about with
steer,
## is fixed with respect to rear frame, and is the center of a torus
## swept by the wheel
turnAnglePntZ =
frontWheelRadius+(math.cos(math.pi-headAngle)*mechTrail)#
turnAnglePntX = wheelBase-(math.sin(math.pi-headAngle)*mechTrail)#
turnAnglePnt2RearAxle =
((frontWheelRadius-turnAnglePntZ)**2+turnAnglePntX**2)**.5#
## front axle moves along the x-z plane with steer, sweeping
through the torus
## steerAngle==pi at 180
frontAxleZ = -math.cos(headAngle)*mechTrail*(1 -
math.cos(steerAngle))
frontAxleX = wheelBase - math.sin(headAngle)*mechTrail*(1 -
math.cos(steerAngle))
# print "frontAxleZ=", frontAxleZ, "\n"
# print "frontAxleX=", frontAxleX, "\n"
## angle between initial place, 0 steer, and current steer, frame
not moved yet...
## relative to rear axle height; add any existing angle
## now, rotate the steer, calc angle due to steer without trail -
holding frame fixed
## max fWheel drop is prop. to steer Angle==90
SA=steerAngle
HA=headAngle
R=frontWheelRadius
cosHA=math.cos(HA)
sinHA=math.sin(HA)
cosSA=math.cos(SA)
sinSA=math.sin(HA)
t1=(sinHA*sinHA + cosSA*(1-sinHA*sinHA))
t2=(cosSA*cosHA*sinHA - cosHA*sinHA)
if(t2 == 0):
t=(-math.pi/2)
else:
t= math.atan(t1 / t2)
Z=-cosHA*sinHA*R*math.cos(t) + cosSA*cosHA*sinHA*R*math.cos(t) +
sinHA*sinHA*R*math.sin(t) + cosSA*(1-sinHA*sinHA)*R*math.sin(t)
# print "Z=", Z, "\n";
X=frontAxleX-cosHA*cosHA*R*math.cos(t) +
cosSA*(1-cosHA*cosHA)*R*math.cos(t) -sinHA*cosHA*R*math.sin(t) +
cosSA*sinHA*cosHA*R*math.sin(t)
# print "X=", X, "\n";
combinedDropAngle = math.atan((R+frontAxleZ+Z)/X)
#steerDropAngle = math.asin(
(math.sin(steerAngle)*math.cos(headAngle)*frontWheelRadius) /
(wheelBase+math.sin(steerAngle)*math.cos(headAngle)*frontWheelRadius) )
## subtract the amount the wheel drops due to flop, from the amount
## the frame rises due to the whole wheel sweeping about the center
of the torus
## frontAxleAngle-steerDropAngle
## calc rearMassHeight
## so, hypotenous == rmass2Raxle
newrmass2RaxleAngle = rmass2RaxleAngle-combinedDropAngle
newrearMassZ =
rearWheelRadius+rmass2Raxle*math.sin(newrmass2RaxleAngle)#
#print '\t', round(newrearMassZ-rearMassZ, 2),
#print '\t', round(rearPitchAngle/math.pi * 180.0, 2),
#print ' ', round(frontAxle2Ground, 1),
#print ' ', round(turnAnglePntZ, 1),
print '\t', round(newrearMassZ-rearMassZ, 3),
if __name__ == '__main__':
main()
- Follow-Ups:
- [python] Re: Bicycling physics
- From: Olaf Johansson
- [python] Re: Bicycling physics
- From: Dirk Bonne
- [python] Re: Bicycling physics
- References:
- [python] Re: Bicycling physics
- From: mtb
- [python] Re: Bicycling physics
