test_functions.py 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509
  1. from sympy.core.numbers import pi
  2. from sympy.core.singleton import S
  3. from sympy.core.symbol import symbols
  4. from sympy.functions.elementary.miscellaneous import sqrt
  5. from sympy.functions.elementary.trigonometric import (cos, sin)
  6. from sympy.integrals.integrals import Integral
  7. from sympy.physics.vector import Dyadic, Point, ReferenceFrame, Vector
  8. from sympy.physics.vector.functions import (cross, dot, express,
  9. time_derivative,
  10. kinematic_equations, outer,
  11. partial_velocity,
  12. get_motion_params, dynamicsymbols)
  13. from sympy.simplify import trigsimp
  14. from sympy.testing.pytest import raises
  15. q1, q2, q3, q4, q5 = symbols('q1 q2 q3 q4 q5')
  16. N = ReferenceFrame('N')
  17. A = N.orientnew('A', 'Axis', [q1, N.z])
  18. B = A.orientnew('B', 'Axis', [q2, A.x])
  19. C = B.orientnew('C', 'Axis', [q3, B.y])
  20. def test_dot():
  21. assert dot(A.x, A.x) == 1
  22. assert dot(A.x, A.y) == 0
  23. assert dot(A.x, A.z) == 0
  24. assert dot(A.y, A.x) == 0
  25. assert dot(A.y, A.y) == 1
  26. assert dot(A.y, A.z) == 0
  27. assert dot(A.z, A.x) == 0
  28. assert dot(A.z, A.y) == 0
  29. assert dot(A.z, A.z) == 1
  30. def test_dot_different_frames():
  31. assert dot(N.x, A.x) == cos(q1)
  32. assert dot(N.x, A.y) == -sin(q1)
  33. assert dot(N.x, A.z) == 0
  34. assert dot(N.y, A.x) == sin(q1)
  35. assert dot(N.y, A.y) == cos(q1)
  36. assert dot(N.y, A.z) == 0
  37. assert dot(N.z, A.x) == 0
  38. assert dot(N.z, A.y) == 0
  39. assert dot(N.z, A.z) == 1
  40. assert trigsimp(dot(N.x, A.x + A.y)) == sqrt(2)*cos(q1 + pi/4)
  41. assert trigsimp(dot(N.x, A.x + A.y)) == trigsimp(dot(A.x + A.y, N.x))
  42. assert dot(A.x, C.x) == cos(q3)
  43. assert dot(A.x, C.y) == 0
  44. assert dot(A.x, C.z) == sin(q3)
  45. assert dot(A.y, C.x) == sin(q2)*sin(q3)
  46. assert dot(A.y, C.y) == cos(q2)
  47. assert dot(A.y, C.z) == -sin(q2)*cos(q3)
  48. assert dot(A.z, C.x) == -cos(q2)*sin(q3)
  49. assert dot(A.z, C.y) == sin(q2)
  50. assert dot(A.z, C.z) == cos(q2)*cos(q3)
  51. def test_cross():
  52. assert cross(A.x, A.x) == 0
  53. assert cross(A.x, A.y) == A.z
  54. assert cross(A.x, A.z) == -A.y
  55. assert cross(A.y, A.x) == -A.z
  56. assert cross(A.y, A.y) == 0
  57. assert cross(A.y, A.z) == A.x
  58. assert cross(A.z, A.x) == A.y
  59. assert cross(A.z, A.y) == -A.x
  60. assert cross(A.z, A.z) == 0
  61. def test_cross_different_frames():
  62. assert cross(N.x, A.x) == sin(q1)*A.z
  63. assert cross(N.x, A.y) == cos(q1)*A.z
  64. assert cross(N.x, A.z) == -sin(q1)*A.x - cos(q1)*A.y
  65. assert cross(N.y, A.x) == -cos(q1)*A.z
  66. assert cross(N.y, A.y) == sin(q1)*A.z
  67. assert cross(N.y, A.z) == cos(q1)*A.x - sin(q1)*A.y
  68. assert cross(N.z, A.x) == A.y
  69. assert cross(N.z, A.y) == -A.x
  70. assert cross(N.z, A.z) == 0
  71. assert cross(N.x, A.x) == sin(q1)*A.z
  72. assert cross(N.x, A.y) == cos(q1)*A.z
  73. assert cross(N.x, A.x + A.y) == sin(q1)*A.z + cos(q1)*A.z
  74. assert cross(A.x + A.y, N.x) == -sin(q1)*A.z - cos(q1)*A.z
  75. assert cross(A.x, C.x) == sin(q3)*C.y
  76. assert cross(A.x, C.y) == -sin(q3)*C.x + cos(q3)*C.z
  77. assert cross(A.x, C.z) == -cos(q3)*C.y
  78. assert cross(C.x, A.x) == -sin(q3)*C.y
  79. assert cross(C.y, A.x).express(C).simplify() == sin(q3)*C.x - cos(q3)*C.z
  80. assert cross(C.z, A.x) == cos(q3)*C.y
  81. def test_operator_match():
  82. """Test that the output of dot, cross, outer functions match
  83. operator behavior.
  84. """
  85. A = ReferenceFrame('A')
  86. v = A.x + A.y
  87. d = v | v
  88. zerov = Vector(0)
  89. zerod = Dyadic(0)
  90. # dot products
  91. assert d & d == dot(d, d)
  92. assert d & zerod == dot(d, zerod)
  93. assert zerod & d == dot(zerod, d)
  94. assert d & v == dot(d, v)
  95. assert v & d == dot(v, d)
  96. assert d & zerov == dot(d, zerov)
  97. assert zerov & d == dot(zerov, d)
  98. raises(TypeError, lambda: dot(d, S.Zero))
  99. raises(TypeError, lambda: dot(S.Zero, d))
  100. raises(TypeError, lambda: dot(d, 0))
  101. raises(TypeError, lambda: dot(0, d))
  102. assert v & v == dot(v, v)
  103. assert v & zerov == dot(v, zerov)
  104. assert zerov & v == dot(zerov, v)
  105. raises(TypeError, lambda: dot(v, S.Zero))
  106. raises(TypeError, lambda: dot(S.Zero, v))
  107. raises(TypeError, lambda: dot(v, 0))
  108. raises(TypeError, lambda: dot(0, v))
  109. # cross products
  110. raises(TypeError, lambda: cross(d, d))
  111. raises(TypeError, lambda: cross(d, zerod))
  112. raises(TypeError, lambda: cross(zerod, d))
  113. assert d ^ v == cross(d, v)
  114. assert v ^ d == cross(v, d)
  115. assert d ^ zerov == cross(d, zerov)
  116. assert zerov ^ d == cross(zerov, d)
  117. assert zerov ^ d == cross(zerov, d)
  118. raises(TypeError, lambda: cross(d, S.Zero))
  119. raises(TypeError, lambda: cross(S.Zero, d))
  120. raises(TypeError, lambda: cross(d, 0))
  121. raises(TypeError, lambda: cross(0, d))
  122. assert v ^ v == cross(v, v)
  123. assert v ^ zerov == cross(v, zerov)
  124. assert zerov ^ v == cross(zerov, v)
  125. raises(TypeError, lambda: cross(v, S.Zero))
  126. raises(TypeError, lambda: cross(S.Zero, v))
  127. raises(TypeError, lambda: cross(v, 0))
  128. raises(TypeError, lambda: cross(0, v))
  129. # outer products
  130. raises(TypeError, lambda: outer(d, d))
  131. raises(TypeError, lambda: outer(d, zerod))
  132. raises(TypeError, lambda: outer(zerod, d))
  133. raises(TypeError, lambda: outer(d, v))
  134. raises(TypeError, lambda: outer(v, d))
  135. raises(TypeError, lambda: outer(d, zerov))
  136. raises(TypeError, lambda: outer(zerov, d))
  137. raises(TypeError, lambda: outer(zerov, d))
  138. raises(TypeError, lambda: outer(d, S.Zero))
  139. raises(TypeError, lambda: outer(S.Zero, d))
  140. raises(TypeError, lambda: outer(d, 0))
  141. raises(TypeError, lambda: outer(0, d))
  142. assert v | v == outer(v, v)
  143. assert v | zerov == outer(v, zerov)
  144. assert zerov | v == outer(zerov, v)
  145. raises(TypeError, lambda: outer(v, S.Zero))
  146. raises(TypeError, lambda: outer(S.Zero, v))
  147. raises(TypeError, lambda: outer(v, 0))
  148. raises(TypeError, lambda: outer(0, v))
  149. def test_express():
  150. assert express(Vector(0), N) == Vector(0)
  151. assert express(S.Zero, N) is S.Zero
  152. assert express(A.x, C) == cos(q3)*C.x + sin(q3)*C.z
  153. assert express(A.y, C) == sin(q2)*sin(q3)*C.x + cos(q2)*C.y - \
  154. sin(q2)*cos(q3)*C.z
  155. assert express(A.z, C) == -sin(q3)*cos(q2)*C.x + sin(q2)*C.y + \
  156. cos(q2)*cos(q3)*C.z
  157. assert express(A.x, N) == cos(q1)*N.x + sin(q1)*N.y
  158. assert express(A.y, N) == -sin(q1)*N.x + cos(q1)*N.y
  159. assert express(A.z, N) == N.z
  160. assert express(A.x, A) == A.x
  161. assert express(A.y, A) == A.y
  162. assert express(A.z, A) == A.z
  163. assert express(A.x, B) == B.x
  164. assert express(A.y, B) == cos(q2)*B.y - sin(q2)*B.z
  165. assert express(A.z, B) == sin(q2)*B.y + cos(q2)*B.z
  166. assert express(A.x, C) == cos(q3)*C.x + sin(q3)*C.z
  167. assert express(A.y, C) == sin(q2)*sin(q3)*C.x + cos(q2)*C.y - \
  168. sin(q2)*cos(q3)*C.z
  169. assert express(A.z, C) == -sin(q3)*cos(q2)*C.x + sin(q2)*C.y + \
  170. cos(q2)*cos(q3)*C.z
  171. # Check to make sure UnitVectors get converted properly
  172. assert express(N.x, N) == N.x
  173. assert express(N.y, N) == N.y
  174. assert express(N.z, N) == N.z
  175. assert express(N.x, A) == (cos(q1)*A.x - sin(q1)*A.y)
  176. assert express(N.y, A) == (sin(q1)*A.x + cos(q1)*A.y)
  177. assert express(N.z, A) == A.z
  178. assert express(N.x, B) == (cos(q1)*B.x - sin(q1)*cos(q2)*B.y +
  179. sin(q1)*sin(q2)*B.z)
  180. assert express(N.y, B) == (sin(q1)*B.x + cos(q1)*cos(q2)*B.y -
  181. sin(q2)*cos(q1)*B.z)
  182. assert express(N.z, B) == (sin(q2)*B.y + cos(q2)*B.z)
  183. assert express(N.x, C) == (
  184. (cos(q1)*cos(q3) - sin(q1)*sin(q2)*sin(q3))*C.x -
  185. sin(q1)*cos(q2)*C.y +
  186. (sin(q3)*cos(q1) + sin(q1)*sin(q2)*cos(q3))*C.z)
  187. assert express(N.y, C) == (
  188. (sin(q1)*cos(q3) + sin(q2)*sin(q3)*cos(q1))*C.x +
  189. cos(q1)*cos(q2)*C.y +
  190. (sin(q1)*sin(q3) - sin(q2)*cos(q1)*cos(q3))*C.z)
  191. assert express(N.z, C) == (-sin(q3)*cos(q2)*C.x + sin(q2)*C.y +
  192. cos(q2)*cos(q3)*C.z)
  193. assert express(A.x, N) == (cos(q1)*N.x + sin(q1)*N.y)
  194. assert express(A.y, N) == (-sin(q1)*N.x + cos(q1)*N.y)
  195. assert express(A.z, N) == N.z
  196. assert express(A.x, A) == A.x
  197. assert express(A.y, A) == A.y
  198. assert express(A.z, A) == A.z
  199. assert express(A.x, B) == B.x
  200. assert express(A.y, B) == (cos(q2)*B.y - sin(q2)*B.z)
  201. assert express(A.z, B) == (sin(q2)*B.y + cos(q2)*B.z)
  202. assert express(A.x, C) == (cos(q3)*C.x + sin(q3)*C.z)
  203. assert express(A.y, C) == (sin(q2)*sin(q3)*C.x + cos(q2)*C.y -
  204. sin(q2)*cos(q3)*C.z)
  205. assert express(A.z, C) == (-sin(q3)*cos(q2)*C.x + sin(q2)*C.y +
  206. cos(q2)*cos(q3)*C.z)
  207. assert express(B.x, N) == (cos(q1)*N.x + sin(q1)*N.y)
  208. assert express(B.y, N) == (-sin(q1)*cos(q2)*N.x +
  209. cos(q1)*cos(q2)*N.y + sin(q2)*N.z)
  210. assert express(B.z, N) == (sin(q1)*sin(q2)*N.x -
  211. sin(q2)*cos(q1)*N.y + cos(q2)*N.z)
  212. assert express(B.x, A) == A.x
  213. assert express(B.y, A) == (cos(q2)*A.y + sin(q2)*A.z)
  214. assert express(B.z, A) == (-sin(q2)*A.y + cos(q2)*A.z)
  215. assert express(B.x, B) == B.x
  216. assert express(B.y, B) == B.y
  217. assert express(B.z, B) == B.z
  218. assert express(B.x, C) == (cos(q3)*C.x + sin(q3)*C.z)
  219. assert express(B.y, C) == C.y
  220. assert express(B.z, C) == (-sin(q3)*C.x + cos(q3)*C.z)
  221. assert express(C.x, N) == (
  222. (cos(q1)*cos(q3) - sin(q1)*sin(q2)*sin(q3))*N.x +
  223. (sin(q1)*cos(q3) + sin(q2)*sin(q3)*cos(q1))*N.y -
  224. sin(q3)*cos(q2)*N.z)
  225. assert express(C.y, N) == (
  226. -sin(q1)*cos(q2)*N.x + cos(q1)*cos(q2)*N.y + sin(q2)*N.z)
  227. assert express(C.z, N) == (
  228. (sin(q3)*cos(q1) + sin(q1)*sin(q2)*cos(q3))*N.x +
  229. (sin(q1)*sin(q3) - sin(q2)*cos(q1)*cos(q3))*N.y +
  230. cos(q2)*cos(q3)*N.z)
  231. assert express(C.x, A) == (cos(q3)*A.x + sin(q2)*sin(q3)*A.y -
  232. sin(q3)*cos(q2)*A.z)
  233. assert express(C.y, A) == (cos(q2)*A.y + sin(q2)*A.z)
  234. assert express(C.z, A) == (sin(q3)*A.x - sin(q2)*cos(q3)*A.y +
  235. cos(q2)*cos(q3)*A.z)
  236. assert express(C.x, B) == (cos(q3)*B.x - sin(q3)*B.z)
  237. assert express(C.y, B) == B.y
  238. assert express(C.z, B) == (sin(q3)*B.x + cos(q3)*B.z)
  239. assert express(C.x, C) == C.x
  240. assert express(C.y, C) == C.y
  241. assert express(C.z, C) == C.z == (C.z)
  242. # Check to make sure Vectors get converted back to UnitVectors
  243. assert N.x == express((cos(q1)*A.x - sin(q1)*A.y), N).simplify()
  244. assert N.y == express((sin(q1)*A.x + cos(q1)*A.y), N).simplify()
  245. assert N.x == express((cos(q1)*B.x - sin(q1)*cos(q2)*B.y +
  246. sin(q1)*sin(q2)*B.z), N).simplify()
  247. assert N.y == express((sin(q1)*B.x + cos(q1)*cos(q2)*B.y -
  248. sin(q2)*cos(q1)*B.z), N).simplify()
  249. assert N.z == express((sin(q2)*B.y + cos(q2)*B.z), N).simplify()
  250. """
  251. These don't really test our code, they instead test the auto simplification
  252. (or lack thereof) of SymPy.
  253. assert N.x == express((
  254. (cos(q1)*cos(q3)-sin(q1)*sin(q2)*sin(q3))*C.x -
  255. sin(q1)*cos(q2)*C.y +
  256. (sin(q3)*cos(q1)+sin(q1)*sin(q2)*cos(q3))*C.z), N)
  257. assert N.y == express((
  258. (sin(q1)*cos(q3) + sin(q2)*sin(q3)*cos(q1))*C.x +
  259. cos(q1)*cos(q2)*C.y +
  260. (sin(q1)*sin(q3) - sin(q2)*cos(q1)*cos(q3))*C.z), N)
  261. assert N.z == express((-sin(q3)*cos(q2)*C.x + sin(q2)*C.y +
  262. cos(q2)*cos(q3)*C.z), N)
  263. """
  264. assert A.x == express((cos(q1)*N.x + sin(q1)*N.y), A).simplify()
  265. assert A.y == express((-sin(q1)*N.x + cos(q1)*N.y), A).simplify()
  266. assert A.y == express((cos(q2)*B.y - sin(q2)*B.z), A).simplify()
  267. assert A.z == express((sin(q2)*B.y + cos(q2)*B.z), A).simplify()
  268. assert A.x == express((cos(q3)*C.x + sin(q3)*C.z), A).simplify()
  269. # Tripsimp messes up here too.
  270. #print express((sin(q2)*sin(q3)*C.x + cos(q2)*C.y -
  271. # sin(q2)*cos(q3)*C.z), A)
  272. assert A.y == express((sin(q2)*sin(q3)*C.x + cos(q2)*C.y -
  273. sin(q2)*cos(q3)*C.z), A).simplify()
  274. assert A.z == express((-sin(q3)*cos(q2)*C.x + sin(q2)*C.y +
  275. cos(q2)*cos(q3)*C.z), A).simplify()
  276. assert B.x == express((cos(q1)*N.x + sin(q1)*N.y), B).simplify()
  277. assert B.y == express((-sin(q1)*cos(q2)*N.x +
  278. cos(q1)*cos(q2)*N.y + sin(q2)*N.z), B).simplify()
  279. assert B.z == express((sin(q1)*sin(q2)*N.x -
  280. sin(q2)*cos(q1)*N.y + cos(q2)*N.z), B).simplify()
  281. assert B.y == express((cos(q2)*A.y + sin(q2)*A.z), B).simplify()
  282. assert B.z == express((-sin(q2)*A.y + cos(q2)*A.z), B).simplify()
  283. assert B.x == express((cos(q3)*C.x + sin(q3)*C.z), B).simplify()
  284. assert B.z == express((-sin(q3)*C.x + cos(q3)*C.z), B).simplify()
  285. """
  286. assert C.x == express((
  287. (cos(q1)*cos(q3)-sin(q1)*sin(q2)*sin(q3))*N.x +
  288. (sin(q1)*cos(q3)+sin(q2)*sin(q3)*cos(q1))*N.y -
  289. sin(q3)*cos(q2)*N.z), C)
  290. assert C.y == express((
  291. -sin(q1)*cos(q2)*N.x + cos(q1)*cos(q2)*N.y + sin(q2)*N.z), C)
  292. assert C.z == express((
  293. (sin(q3)*cos(q1)+sin(q1)*sin(q2)*cos(q3))*N.x +
  294. (sin(q1)*sin(q3)-sin(q2)*cos(q1)*cos(q3))*N.y +
  295. cos(q2)*cos(q3)*N.z), C)
  296. """
  297. assert C.x == express((cos(q3)*A.x + sin(q2)*sin(q3)*A.y -
  298. sin(q3)*cos(q2)*A.z), C).simplify()
  299. assert C.y == express((cos(q2)*A.y + sin(q2)*A.z), C).simplify()
  300. assert C.z == express((sin(q3)*A.x - sin(q2)*cos(q3)*A.y +
  301. cos(q2)*cos(q3)*A.z), C).simplify()
  302. assert C.x == express((cos(q3)*B.x - sin(q3)*B.z), C).simplify()
  303. assert C.z == express((sin(q3)*B.x + cos(q3)*B.z), C).simplify()
  304. def test_time_derivative():
  305. #The use of time_derivative for calculations pertaining to scalar
  306. #fields has been tested in test_coordinate_vars in test_essential.py
  307. A = ReferenceFrame('A')
  308. q = dynamicsymbols('q')
  309. qd = dynamicsymbols('q', 1)
  310. B = A.orientnew('B', 'Axis', [q, A.z])
  311. d = A.x | A.x
  312. assert time_derivative(d, B) == (-qd) * (A.y | A.x) + \
  313. (-qd) * (A.x | A.y)
  314. d1 = A.x | B.y
  315. assert time_derivative(d1, A) == - qd*(A.x|B.x)
  316. assert time_derivative(d1, B) == - qd*(A.y|B.y)
  317. d2 = A.x | B.x
  318. assert time_derivative(d2, A) == qd*(A.x|B.y)
  319. assert time_derivative(d2, B) == - qd*(A.y|B.x)
  320. d3 = A.x | B.z
  321. assert time_derivative(d3, A) == 0
  322. assert time_derivative(d3, B) == - qd*(A.y|B.z)
  323. q1, q2, q3, q4 = dynamicsymbols('q1 q2 q3 q4')
  324. q1d, q2d, q3d, q4d = dynamicsymbols('q1 q2 q3 q4', 1)
  325. q1dd, q2dd, q3dd, q4dd = dynamicsymbols('q1 q2 q3 q4', 2)
  326. C = B.orientnew('C', 'Axis', [q4, B.x])
  327. v1 = q1 * A.z
  328. v2 = q2*A.x + q3*B.y
  329. v3 = q1*A.x + q2*A.y + q3*A.z
  330. assert time_derivative(B.x, C) == 0
  331. assert time_derivative(B.y, C) == - q4d*B.z
  332. assert time_derivative(B.z, C) == q4d*B.y
  333. assert time_derivative(v1, B) == q1d*A.z
  334. assert time_derivative(v1, C) == - q1*sin(q)*q4d*A.x + \
  335. q1*cos(q)*q4d*A.y + q1d*A.z
  336. assert time_derivative(v2, A) == q2d*A.x - q3*qd*B.x + q3d*B.y
  337. assert time_derivative(v2, C) == q2d*A.x - q2*qd*A.y + \
  338. q2*sin(q)*q4d*A.z + q3d*B.y - q3*q4d*B.z
  339. assert time_derivative(v3, B) == (q2*qd + q1d)*A.x + \
  340. (-q1*qd + q2d)*A.y + q3d*A.z
  341. assert time_derivative(d, C) == - qd*(A.y|A.x) + \
  342. sin(q)*q4d*(A.z|A.x) - qd*(A.x|A.y) + sin(q)*q4d*(A.x|A.z)
  343. raises(ValueError, lambda: time_derivative(B.x, C, order=0.5))
  344. raises(ValueError, lambda: time_derivative(B.x, C, order=-1))
  345. def test_get_motion_methods():
  346. #Initialization
  347. t = dynamicsymbols._t
  348. s1, s2, s3 = symbols('s1 s2 s3')
  349. S1, S2, S3 = symbols('S1 S2 S3')
  350. S4, S5, S6 = symbols('S4 S5 S6')
  351. t1, t2 = symbols('t1 t2')
  352. a, b, c = dynamicsymbols('a b c')
  353. ad, bd, cd = dynamicsymbols('a b c', 1)
  354. a2d, b2d, c2d = dynamicsymbols('a b c', 2)
  355. v0 = S1*N.x + S2*N.y + S3*N.z
  356. v01 = S4*N.x + S5*N.y + S6*N.z
  357. v1 = s1*N.x + s2*N.y + s3*N.z
  358. v2 = a*N.x + b*N.y + c*N.z
  359. v2d = ad*N.x + bd*N.y + cd*N.z
  360. v2dd = a2d*N.x + b2d*N.y + c2d*N.z
  361. #Test position parameter
  362. assert get_motion_params(frame = N) == (0, 0, 0)
  363. assert get_motion_params(N, position=v1) == (0, 0, v1)
  364. assert get_motion_params(N, position=v2) == (v2dd, v2d, v2)
  365. #Test velocity parameter
  366. assert get_motion_params(N, velocity=v1) == (0, v1, v1 * t)
  367. assert get_motion_params(N, velocity=v1, position=v0, timevalue1=t1) == \
  368. (0, v1, v0 + v1*(t - t1))
  369. answer = get_motion_params(N, velocity=v1, position=v2, timevalue1=t1)
  370. answer_expected = (0, v1, v1*t - v1*t1 + v2.subs(t, t1))
  371. assert answer == answer_expected
  372. answer = get_motion_params(N, velocity=v2, position=v0, timevalue1=t1)
  373. integral_vector = Integral(a, (t, t1, t))*N.x + Integral(b, (t, t1, t))*N.y \
  374. + Integral(c, (t, t1, t))*N.z
  375. answer_expected = (v2d, v2, v0 + integral_vector)
  376. assert answer == answer_expected
  377. #Test acceleration parameter
  378. assert get_motion_params(N, acceleration=v1) == \
  379. (v1, v1 * t, v1 * t**2/2)
  380. assert get_motion_params(N, acceleration=v1, velocity=v0,
  381. position=v2, timevalue1=t1, timevalue2=t2) == \
  382. (v1, (v0 + v1*t - v1*t2),
  383. -v0*t1 + v1*t**2/2 + v1*t2*t1 - \
  384. v1*t1**2/2 + t*(v0 - v1*t2) + \
  385. v2.subs(t, t1))
  386. assert get_motion_params(N, acceleration=v1, velocity=v0,
  387. position=v01, timevalue1=t1, timevalue2=t2) == \
  388. (v1, v0 + v1*t - v1*t2,
  389. -v0*t1 + v01 + v1*t**2/2 + \
  390. v1*t2*t1 - v1*t1**2/2 + \
  391. t*(v0 - v1*t2))
  392. answer = get_motion_params(N, acceleration=a*N.x, velocity=S1*N.x,
  393. position=S2*N.x, timevalue1=t1, timevalue2=t2)
  394. i1 = Integral(a, (t, t2, t))
  395. answer_expected = (a*N.x, (S1 + i1)*N.x, \
  396. (S2 + Integral(S1 + i1, (t, t1, t)))*N.x)
  397. assert answer == answer_expected
  398. def test_kin_eqs():
  399. q0, q1, q2, q3 = dynamicsymbols('q0 q1 q2 q3')
  400. q0d, q1d, q2d, q3d = dynamicsymbols('q0 q1 q2 q3', 1)
  401. u1, u2, u3 = dynamicsymbols('u1 u2 u3')
  402. ke = kinematic_equations([u1,u2,u3], [q1,q2,q3], 'body', 313)
  403. assert ke == kinematic_equations([u1,u2,u3], [q1,q2,q3], 'body', '313')
  404. kds = kinematic_equations([u1, u2, u3], [q0, q1, q2, q3], 'quaternion')
  405. assert kds == [-0.5 * q0 * u1 - 0.5 * q2 * u3 + 0.5 * q3 * u2 + q1d,
  406. -0.5 * q0 * u2 + 0.5 * q1 * u3 - 0.5 * q3 * u1 + q2d,
  407. -0.5 * q0 * u3 - 0.5 * q1 * u2 + 0.5 * q2 * u1 + q3d,
  408. 0.5 * q1 * u1 + 0.5 * q2 * u2 + 0.5 * q3 * u3 + q0d]
  409. raises(ValueError, lambda: kinematic_equations([u1, u2, u3], [q0, q1, q2], 'quaternion'))
  410. raises(ValueError, lambda: kinematic_equations([u1, u2, u3], [q0, q1, q2, q3], 'quaternion', '123'))
  411. raises(ValueError, lambda: kinematic_equations([u1, u2, u3], [q0, q1, q2, q3], 'foo'))
  412. raises(TypeError, lambda: kinematic_equations(u1, [q0, q1, q2, q3], 'quaternion'))
  413. raises(TypeError, lambda: kinematic_equations([u1], [q0, q1, q2, q3], 'quaternion'))
  414. raises(TypeError, lambda: kinematic_equations([u1, u2, u3], q0, 'quaternion'))
  415. raises(ValueError, lambda: kinematic_equations([u1, u2, u3], [q0, q1, q2, q3], 'body'))
  416. raises(ValueError, lambda: kinematic_equations([u1, u2, u3], [q0, q1, q2, q3], 'space'))
  417. raises(ValueError, lambda: kinematic_equations([u1, u2, u3], [q0, q1, q2], 'body', '222'))
  418. assert kinematic_equations([0, 0, 0], [q0, q1, q2], 'space') == [S.Zero, S.Zero, S.Zero]
  419. def test_partial_velocity():
  420. q1, q2, q3, u1, u2, u3 = dynamicsymbols('q1 q2 q3 u1 u2 u3')
  421. u4, u5 = dynamicsymbols('u4, u5')
  422. r = symbols('r')
  423. N = ReferenceFrame('N')
  424. Y = N.orientnew('Y', 'Axis', [q1, N.z])
  425. L = Y.orientnew('L', 'Axis', [q2, Y.x])
  426. R = L.orientnew('R', 'Axis', [q3, L.y])
  427. R.set_ang_vel(N, u1 * L.x + u2 * L.y + u3 * L.z)
  428. C = Point('C')
  429. C.set_vel(N, u4 * L.x + u5 * (Y.z ^ L.x))
  430. Dmc = C.locatenew('Dmc', r * L.z)
  431. Dmc.v2pt_theory(C, N, R)
  432. vel_list = [Dmc.vel(N), C.vel(N), R.ang_vel_in(N)]
  433. u_list = [u1, u2, u3, u4, u5]
  434. assert (partial_velocity(vel_list, u_list, N) ==
  435. [[- r*L.y, r*L.x, 0, L.x, cos(q2)*L.y - sin(q2)*L.z],
  436. [0, 0, 0, L.x, cos(q2)*L.y - sin(q2)*L.z],
  437. [L.x, L.y, L.z, 0, 0]])
  438. # Make sure that partial velocities can be computed regardless if the
  439. # orientation between frames is defined or not.
  440. A = ReferenceFrame('A')
  441. B = ReferenceFrame('B')
  442. v = u4 * A.x + u5 * B.y
  443. assert partial_velocity((v, ), (u4, u5), A) == [[A.x, B.y]]
  444. raises(TypeError, lambda: partial_velocity(Dmc.vel(N), u_list, N))
  445. raises(TypeError, lambda: partial_velocity(vel_list, u1, N))
  446. def test_dynamicsymbols():
  447. #Tests to check the assumptions applied to dynamicsymbols
  448. f1 = dynamicsymbols('f1')
  449. f2 = dynamicsymbols('f2', real=True)
  450. f3 = dynamicsymbols('f3', positive=True)
  451. f4, f5 = dynamicsymbols('f4,f5', commutative=False)
  452. f6 = dynamicsymbols('f6', integer=True)
  453. assert f1.is_real is None
  454. assert f2.is_real
  455. assert f3.is_positive
  456. assert f4*f5 != f5*f4
  457. assert f6.is_integer