mechanizmos58

vpython_Tarapov_page_57_5_springs.py

Sep 18th, 2021
988
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1.  
  2. # Tarapov_page_57_5_springs.py
  3.  
  4. # 19_09_21
  5.  
  6. '''
  7. Vector and Tensor Analysis with Applications
  8. Borisenko, Tarapov, Silverman
  9. Exercise 21, page 57
  10.  
  11. "...springs that have relaxed lengths of zero":
  12. Introduction to Classical Mechanics
  13. by David Morin (2017), p. 124
  14.  
  15.  
  16. Classical Mechanics with Calculus of Variations...
  17. by Mark Levi (2014), pp. 126-127
  18. p. 126
  19. "linear zero length spring"
  20. p. 127
  21. "...the moment of inertia can be interpreted as the potential
  22.    energy of a set of linear zero length springs." #todo
  23.  
  24. Morin p. 690
  25.  
  26. F = ma vs. F = dp/dt
  27.  
  28. F("system") = dp/dt
  29.  
  30.            = m(dv/dt) + v(dm/dt)
  31.  
  32. cf. "The Momentum Principle"
  33. Matter & Interactions (2007)
  34. by Chabay, Sherwood, pp.45-
  35.  
  36.  
  37. '''
  38.  
  39. from vpython import *
  40. import sys
  41.  
  42.  
  43. def main():
  44.  
  45.     scene.width, scene.height = 1000, 800
  46.     scene.background = color.black
  47.     scene.align = 'left'
  48.     scene.autoscale = False
  49.     scene.autozoom = False
  50.     scene.userspin = True
  51.     scene.userzoom = True
  52.     scene.userpan = True
  53.     scene.resizable = False
  54.     scene.center = vec(0,0,0)
  55.     scene.range = 1.7
  56.  
  57.     RATE = 200
  58.  
  59.     PAUSED_AT_START = 0
  60.  
  61.     RANDOM_K = 0
  62.  
  63.     # Fdrag = -DRAG_CONSTANT * p
  64.     DRAG_CONSTANT = 0.6
  65.  
  66.     if 1:
  67.  
  68.         N = 5
  69.  
  70.         L = [vec(-0.3, 0.6, 0),
  71.              vec(0.8, 0.9, 0),
  72.              vec(1.2, 0, 0),
  73.              vec(0.7, -1, 0),
  74.              vec(-0.6, -1, 0)]
  75.  
  76.     else:
  77.  
  78.         N = 12
  79.  
  80.         L = []
  81.         v = vec(1,0,0)
  82.         L.append(v)
  83.         for _ in range(N-1):
  84.             v = v.rotate(radians(360./N), vec(0,0,1))
  85.             # ~ v.rotate_in_place(radians(360./N), vec(0,0,1))
  86.             L.append(v)
  87.  
  88.     # CENTER MASS, M
  89.     x = -1.0 + 2.0*random()
  90.     y = -1.0 + 2.0*random()
  91.     M = sphere(pos=vec(x,y,0), radius=0.05, color=vec(1,0,0))
  92.     M.p = vec(0,0,0)
  93.     M.mass = 1.0
  94.  
  95.     # MIRROR MASS, M2
  96.     dx = 1.7
  97.     M2 = simple_sphere(pos=M.pos - vec(dx,0,0), radius=0.007,
  98.         color=vec(1,0,0), make_trail = 1, trail_radius=0.004)
  99.  
  100.     # INITIALIZE MASS SPRING SYSTEM
  101.     points = []
  102.     for i in range(N):
  103.         point = simple_sphere(pos=L[i], radius=0.02, color=vec(1,1,1))
  104.         point.m = 1.0
  105.         points.append(point)
  106.  
  107.     springs = []
  108.     for i in range(N):
  109.         spring = helix(pos=M.pos, axis=points[i].pos - M.pos,
  110.                   coils=18, radius=0.02, thickness=0.005, opacity=0.7)
  111.         if RANDOM_K:
  112.             spring.k = random() * 2.0
  113.         else:
  114.             spring.k = 2.0
  115.         spring.L0 = 0.0 # "zero-length spring"
  116.         springs.append(spring)
  117.  
  118.  
  119.     '''
  120.    CALCULATE EQUILIBRIUM POSITION OF M
  121.    ===================================
  122.    This formula is accurate only for springs
  123.        that have "relaxed lengths of zero" (Morin)
  124.    Formula from:
  125.    Vector and Tensor Analysis with Applications
  126.    Borisenko, Tarapov, Silverman
  127.    Exercise 21, page 57
  128.  
  129.    This formula is analogous to the Center of Mass
  130.    formula (substituting m for k)
  131.    '''
  132.     R = vec(0.,0.,0.)
  133.     K = 0.
  134.     for i in range(N):
  135.         k = springs[i].k
  136.         R += k * points[i].pos
  137.         K += k
  138.     R /= K
  139.     sphere(pos=R, radius=0.02, color=color.green)
  140.  
  141.  
  142.     def pause_run():
  143.         g.paused = not g.paused
  144.         pause.text='  RUN   ' if g.paused else 'PAUSE'
  145.     scene.append_to_caption('\n   ')
  146.     pause = button(bind=pause_run)
  147.     g.paused = not PAUSED_AT_START
  148.     pause_run()
  149.  
  150.  
  151.     def key(e):
  152.         s = e.key
  153.         if s in [' ', '\n']: pause_run()
  154.     scene.bind('keydown', key)
  155.  
  156.  
  157.     # IMPART SOME INITIAL MOMENTUM
  158.     M.p = vec(0.5, 0.5,0)
  159.  
  160.     t = 0
  161.     dt = 0.002
  162.     while t < 30:
  163.  
  164.         rate(RATE)
  165.  
  166.         if g.paused: continue
  167.  
  168.         for i in range(N):
  169.  
  170.             r = M.pos - points[i].pos
  171.  
  172.             Fs = -springs[i].k * (r.mag - springs[i].L0) * r.hat
  173.             Fdrag = -DRAG_CONSTANT * M.p
  174.  
  175.             M.p += (Fs + Fdrag) * dt
  176.  
  177.             M.pos += M.p / M.mass * dt
  178.  
  179.             springs[i].pos = M.pos
  180.             springs[i].axis = points[i].pos - M.pos
  181.  
  182.             M2.pos = M.pos - vec(dx, 0, 0)
  183.  
  184.             t += dt
  185.  
  186.  
  187.     print()
  188.     s1 = '   Predicted Equilibrium point: '
  189.     s2 = '   Actual Equilibrium point:    '
  190.     print(s1, R)
  191.     print(s2, M.pos)
  192.  
  193.     scene.append_to_caption('\n\n')
  194.     scene.append_to_caption(s1 + str(R) + '\n')
  195.     scene.append_to_caption('\n')
  196.     scene.append_to_caption(s2 + '   ' + str(M.pos) + '\n')
  197.  
  198.  
  199. class g: ...
  200.  
  201. if __name__ == '__main__':
  202.     sys.exit(main())
  203.  
  204.  
RAW Paste Data