test_frame.py 25 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660
  1. from sympy.core.numbers import pi
  2. from sympy.core.symbol import symbols
  3. from sympy.functions.elementary.trigonometric import (cos, sin)
  4. from sympy.matrices.dense import (eye, zeros)
  5. from sympy.matrices.immutable import ImmutableDenseMatrix as Matrix
  6. from sympy.simplify.simplify import simplify
  7. from sympy.physics.vector import (ReferenceFrame, Vector, CoordinateSym,
  8. dynamicsymbols, time_derivative, express,
  9. dot)
  10. from sympy.physics.vector.frame import _check_frame
  11. from sympy.physics.vector.vector import VectorTypeError
  12. from sympy.testing.pytest import raises
  13. import warnings
  14. Vector.simp = True
  15. def test_dict_list():
  16. A = ReferenceFrame('A')
  17. B = ReferenceFrame('B')
  18. C = ReferenceFrame('C')
  19. D = ReferenceFrame('D')
  20. E = ReferenceFrame('E')
  21. F = ReferenceFrame('F')
  22. B.orient_axis(A, A.x, 1.0)
  23. C.orient_axis(B, B.x, 1.0)
  24. D.orient_axis(C, C.x, 1.0)
  25. assert D._dict_list(A, 0) == [D, C, B, A]
  26. E.orient_axis(D, D.x, 1.0)
  27. assert C._dict_list(A, 0) == [C, B, A]
  28. assert C._dict_list(E, 0) == [C, D, E]
  29. # only 0, 1, 2 permitted for second argument
  30. raises(ValueError, lambda: C._dict_list(E, 5))
  31. # no connecting path
  32. raises(ValueError, lambda: F._dict_list(A, 0))
  33. def test_coordinate_vars():
  34. """Tests the coordinate variables functionality"""
  35. A = ReferenceFrame('A')
  36. assert CoordinateSym('Ax', A, 0) == A[0]
  37. assert CoordinateSym('Ax', A, 1) == A[1]
  38. assert CoordinateSym('Ax', A, 2) == A[2]
  39. raises(ValueError, lambda: CoordinateSym('Ax', A, 3))
  40. q = dynamicsymbols('q')
  41. qd = dynamicsymbols('q', 1)
  42. assert isinstance(A[0], CoordinateSym) and \
  43. isinstance(A[0], CoordinateSym) and \
  44. isinstance(A[0], CoordinateSym)
  45. assert A.variable_map(A) == {A[0]:A[0], A[1]:A[1], A[2]:A[2]}
  46. assert A[0].frame == A
  47. B = A.orientnew('B', 'Axis', [q, A.z])
  48. assert B.variable_map(A) == {B[2]: A[2], B[1]: -A[0]*sin(q) + A[1]*cos(q),
  49. B[0]: A[0]*cos(q) + A[1]*sin(q)}
  50. assert A.variable_map(B) == {A[0]: B[0]*cos(q) - B[1]*sin(q),
  51. A[1]: B[0]*sin(q) + B[1]*cos(q), A[2]: B[2]}
  52. assert time_derivative(B[0], A) == -A[0]*sin(q)*qd + A[1]*cos(q)*qd
  53. assert time_derivative(B[1], A) == -A[0]*cos(q)*qd - A[1]*sin(q)*qd
  54. assert time_derivative(B[2], A) == 0
  55. assert express(B[0], A, variables=True) == A[0]*cos(q) + A[1]*sin(q)
  56. assert express(B[1], A, variables=True) == -A[0]*sin(q) + A[1]*cos(q)
  57. assert express(B[2], A, variables=True) == A[2]
  58. assert time_derivative(A[0]*A.x + A[1]*A.y + A[2]*A.z, B) == A[1]*qd*A.x - A[0]*qd*A.y
  59. assert time_derivative(B[0]*B.x + B[1]*B.y + B[2]*B.z, A) == - B[1]*qd*B.x + B[0]*qd*B.y
  60. assert express(B[0]*B[1]*B[2], A, variables=True) == \
  61. A[2]*(-A[0]*sin(q) + A[1]*cos(q))*(A[0]*cos(q) + A[1]*sin(q))
  62. assert (time_derivative(B[0]*B[1]*B[2], A) -
  63. (A[2]*(-A[0]**2*cos(2*q) -
  64. 2*A[0]*A[1]*sin(2*q) +
  65. A[1]**2*cos(2*q))*qd)).trigsimp() == 0
  66. assert express(B[0]*B.x + B[1]*B.y + B[2]*B.z, A) == \
  67. (B[0]*cos(q) - B[1]*sin(q))*A.x + (B[0]*sin(q) + \
  68. B[1]*cos(q))*A.y + B[2]*A.z
  69. assert express(B[0]*B.x + B[1]*B.y + B[2]*B.z, A, variables=True) == \
  70. A[0]*A.x + A[1]*A.y + A[2]*A.z
  71. assert express(A[0]*A.x + A[1]*A.y + A[2]*A.z, B) == \
  72. (A[0]*cos(q) + A[1]*sin(q))*B.x + \
  73. (-A[0]*sin(q) + A[1]*cos(q))*B.y + A[2]*B.z
  74. assert express(A[0]*A.x + A[1]*A.y + A[2]*A.z, B, variables=True) == \
  75. B[0]*B.x + B[1]*B.y + B[2]*B.z
  76. N = B.orientnew('N', 'Axis', [-q, B.z])
  77. assert N.variable_map(A) == {N[0]: A[0], N[2]: A[2], N[1]: A[1]}
  78. C = A.orientnew('C', 'Axis', [q, A.x + A.y + A.z])
  79. mapping = A.variable_map(C)
  80. assert mapping[A[0]] == 2*C[0]*cos(q)/3 + C[0]/3 - 2*C[1]*sin(q + pi/6)/3 +\
  81. C[1]/3 - 2*C[2]*cos(q + pi/3)/3 + C[2]/3
  82. assert mapping[A[1]] == -2*C[0]*cos(q + pi/3)/3 + \
  83. C[0]/3 + 2*C[1]*cos(q)/3 + C[1]/3 - 2*C[2]*sin(q + pi/6)/3 + C[2]/3
  84. assert mapping[A[2]] == -2*C[0]*sin(q + pi/6)/3 + C[0]/3 - \
  85. 2*C[1]*cos(q + pi/3)/3 + C[1]/3 + 2*C[2]*cos(q)/3 + C[2]/3
  86. def test_ang_vel():
  87. q1, q2, q3, q4 = dynamicsymbols('q1 q2 q3 q4')
  88. q1d, q2d, q3d, q4d = dynamicsymbols('q1 q2 q3 q4', 1)
  89. N = ReferenceFrame('N')
  90. A = N.orientnew('A', 'Axis', [q1, N.z])
  91. B = A.orientnew('B', 'Axis', [q2, A.x])
  92. C = B.orientnew('C', 'Axis', [q3, B.y])
  93. D = N.orientnew('D', 'Axis', [q4, N.y])
  94. u1, u2, u3 = dynamicsymbols('u1 u2 u3')
  95. assert A.ang_vel_in(N) == (q1d)*A.z
  96. assert B.ang_vel_in(N) == (q2d)*B.x + (q1d)*A.z
  97. assert C.ang_vel_in(N) == (q3d)*C.y + (q2d)*B.x + (q1d)*A.z
  98. A2 = N.orientnew('A2', 'Axis', [q4, N.y])
  99. assert N.ang_vel_in(N) == 0
  100. assert N.ang_vel_in(A) == -q1d*N.z
  101. assert N.ang_vel_in(B) == -q1d*A.z - q2d*B.x
  102. assert N.ang_vel_in(C) == -q1d*A.z - q2d*B.x - q3d*B.y
  103. assert N.ang_vel_in(A2) == -q4d*N.y
  104. assert A.ang_vel_in(N) == q1d*N.z
  105. assert A.ang_vel_in(A) == 0
  106. assert A.ang_vel_in(B) == - q2d*B.x
  107. assert A.ang_vel_in(C) == - q2d*B.x - q3d*B.y
  108. assert A.ang_vel_in(A2) == q1d*N.z - q4d*N.y
  109. assert B.ang_vel_in(N) == q1d*A.z + q2d*A.x
  110. assert B.ang_vel_in(A) == q2d*A.x
  111. assert B.ang_vel_in(B) == 0
  112. assert B.ang_vel_in(C) == -q3d*B.y
  113. assert B.ang_vel_in(A2) == q1d*A.z + q2d*A.x - q4d*N.y
  114. assert C.ang_vel_in(N) == q1d*A.z + q2d*A.x + q3d*B.y
  115. assert C.ang_vel_in(A) == q2d*A.x + q3d*C.y
  116. assert C.ang_vel_in(B) == q3d*B.y
  117. assert C.ang_vel_in(C) == 0
  118. assert C.ang_vel_in(A2) == q1d*A.z + q2d*A.x + q3d*B.y - q4d*N.y
  119. assert A2.ang_vel_in(N) == q4d*A2.y
  120. assert A2.ang_vel_in(A) == q4d*A2.y - q1d*N.z
  121. assert A2.ang_vel_in(B) == q4d*N.y - q1d*A.z - q2d*A.x
  122. assert A2.ang_vel_in(C) == q4d*N.y - q1d*A.z - q2d*A.x - q3d*B.y
  123. assert A2.ang_vel_in(A2) == 0
  124. C.set_ang_vel(N, u1*C.x + u2*C.y + u3*C.z)
  125. assert C.ang_vel_in(N) == (u1)*C.x + (u2)*C.y + (u3)*C.z
  126. assert N.ang_vel_in(C) == (-u1)*C.x + (-u2)*C.y + (-u3)*C.z
  127. assert C.ang_vel_in(D) == (u1)*C.x + (u2)*C.y + (u3)*C.z + (-q4d)*D.y
  128. assert D.ang_vel_in(C) == (-u1)*C.x + (-u2)*C.y + (-u3)*C.z + (q4d)*D.y
  129. q0 = dynamicsymbols('q0')
  130. q0d = dynamicsymbols('q0', 1)
  131. E = N.orientnew('E', 'Quaternion', (q0, q1, q2, q3))
  132. assert E.ang_vel_in(N) == (
  133. 2 * (q1d * q0 + q2d * q3 - q3d * q2 - q0d * q1) * E.x +
  134. 2 * (q2d * q0 + q3d * q1 - q1d * q3 - q0d * q2) * E.y +
  135. 2 * (q3d * q0 + q1d * q2 - q2d * q1 - q0d * q3) * E.z)
  136. F = N.orientnew('F', 'Body', (q1, q2, q3), 313)
  137. assert F.ang_vel_in(N) == ((sin(q2)*sin(q3)*q1d + cos(q3)*q2d)*F.x +
  138. (sin(q2)*cos(q3)*q1d - sin(q3)*q2d)*F.y + (cos(q2)*q1d + q3d)*F.z)
  139. G = N.orientnew('G', 'Axis', (q1, N.x + N.y))
  140. assert G.ang_vel_in(N) == q1d * (N.x + N.y).normalize()
  141. assert N.ang_vel_in(G) == -q1d * (N.x + N.y).normalize()
  142. def test_dcm():
  143. q1, q2, q3, q4 = dynamicsymbols('q1 q2 q3 q4')
  144. N = ReferenceFrame('N')
  145. A = N.orientnew('A', 'Axis', [q1, N.z])
  146. B = A.orientnew('B', 'Axis', [q2, A.x])
  147. C = B.orientnew('C', 'Axis', [q3, B.y])
  148. D = N.orientnew('D', 'Axis', [q4, N.y])
  149. E = N.orientnew('E', 'Space', [q1, q2, q3], '123')
  150. assert N.dcm(C) == Matrix([
  151. [- sin(q1) * sin(q2) * sin(q3) + cos(q1) * cos(q3), - sin(q1) *
  152. cos(q2), sin(q1) * sin(q2) * cos(q3) + sin(q3) * cos(q1)], [sin(q1) *
  153. cos(q3) + sin(q2) * sin(q3) * cos(q1), cos(q1) * cos(q2), sin(q1) *
  154. sin(q3) - sin(q2) * cos(q1) * cos(q3)], [- sin(q3) * cos(q2), sin(q2),
  155. cos(q2) * cos(q3)]])
  156. # This is a little touchy. Is it ok to use simplify in assert?
  157. test_mat = D.dcm(C) - Matrix(
  158. [[cos(q1) * cos(q3) * cos(q4) - sin(q3) * (- sin(q4) * cos(q2) +
  159. sin(q1) * sin(q2) * cos(q4)), - sin(q2) * sin(q4) - sin(q1) *
  160. cos(q2) * cos(q4), sin(q3) * cos(q1) * cos(q4) + cos(q3) * (- sin(q4) *
  161. cos(q2) + sin(q1) * sin(q2) * cos(q4))], [sin(q1) * cos(q3) +
  162. sin(q2) * sin(q3) * cos(q1), cos(q1) * cos(q2), sin(q1) * sin(q3) -
  163. sin(q2) * cos(q1) * cos(q3)], [sin(q4) * cos(q1) * cos(q3) -
  164. sin(q3) * (cos(q2) * cos(q4) + sin(q1) * sin(q2) * sin(q4)), sin(q2) *
  165. cos(q4) - sin(q1) * sin(q4) * cos(q2), sin(q3) * sin(q4) * cos(q1) +
  166. cos(q3) * (cos(q2) * cos(q4) + sin(q1) * sin(q2) * sin(q4))]])
  167. assert test_mat.expand() == zeros(3, 3)
  168. assert E.dcm(N) == Matrix(
  169. [[cos(q2)*cos(q3), sin(q3)*cos(q2), -sin(q2)],
  170. [sin(q1)*sin(q2)*cos(q3) - sin(q3)*cos(q1), sin(q1)*sin(q2)*sin(q3) +
  171. cos(q1)*cos(q3), sin(q1)*cos(q2)], [sin(q1)*sin(q3) +
  172. sin(q2)*cos(q1)*cos(q3), - sin(q1)*cos(q3) + sin(q2)*sin(q3)*cos(q1),
  173. cos(q1)*cos(q2)]])
  174. def test_w_diff_dcm1():
  175. # Ref:
  176. # Dynamics Theory and Applications, Kane 1985
  177. # Sec. 2.1 ANGULAR VELOCITY
  178. A = ReferenceFrame('A')
  179. B = ReferenceFrame('B')
  180. c11, c12, c13 = dynamicsymbols('C11 C12 C13')
  181. c21, c22, c23 = dynamicsymbols('C21 C22 C23')
  182. c31, c32, c33 = dynamicsymbols('C31 C32 C33')
  183. c11d, c12d, c13d = dynamicsymbols('C11 C12 C13', level=1)
  184. c21d, c22d, c23d = dynamicsymbols('C21 C22 C23', level=1)
  185. c31d, c32d, c33d = dynamicsymbols('C31 C32 C33', level=1)
  186. DCM = Matrix([
  187. [c11, c12, c13],
  188. [c21, c22, c23],
  189. [c31, c32, c33]
  190. ])
  191. B.orient(A, 'DCM', DCM)
  192. b1a = (B.x).express(A)
  193. b2a = (B.y).express(A)
  194. b3a = (B.z).express(A)
  195. # Equation (2.1.1)
  196. B.set_ang_vel(A, B.x*(dot((b3a).dt(A), B.y))
  197. + B.y*(dot((b1a).dt(A), B.z))
  198. + B.z*(dot((b2a).dt(A), B.x)))
  199. # Equation (2.1.21)
  200. expr = ( (c12*c13d + c22*c23d + c32*c33d)*B.x
  201. + (c13*c11d + c23*c21d + c33*c31d)*B.y
  202. + (c11*c12d + c21*c22d + c31*c32d)*B.z)
  203. assert B.ang_vel_in(A) - expr == 0
  204. def test_w_diff_dcm2():
  205. q1, q2, q3 = dynamicsymbols('q1:4')
  206. N = ReferenceFrame('N')
  207. A = N.orientnew('A', 'axis', [q1, N.x])
  208. B = A.orientnew('B', 'axis', [q2, A.y])
  209. C = B.orientnew('C', 'axis', [q3, B.z])
  210. DCM = C.dcm(N).T
  211. D = N.orientnew('D', 'DCM', DCM)
  212. # Frames D and C are the same ReferenceFrame,
  213. # since they have equal DCM respect to frame N.
  214. # Therefore, D and C should have same angle velocity in N.
  215. assert D.dcm(N) == C.dcm(N) == Matrix([
  216. [cos(q2)*cos(q3), sin(q1)*sin(q2)*cos(q3) +
  217. sin(q3)*cos(q1), sin(q1)*sin(q3) -
  218. sin(q2)*cos(q1)*cos(q3)], [-sin(q3)*cos(q2),
  219. -sin(q1)*sin(q2)*sin(q3) + cos(q1)*cos(q3),
  220. sin(q1)*cos(q3) + sin(q2)*sin(q3)*cos(q1)],
  221. [sin(q2), -sin(q1)*cos(q2), cos(q1)*cos(q2)]])
  222. assert (D.ang_vel_in(N) - C.ang_vel_in(N)).express(N).simplify() == 0
  223. def test_orientnew_respects_parent_class():
  224. class MyReferenceFrame(ReferenceFrame):
  225. pass
  226. B = MyReferenceFrame('B')
  227. C = B.orientnew('C', 'Axis', [0, B.x])
  228. assert isinstance(C, MyReferenceFrame)
  229. def test_orientnew_respects_input_indices():
  230. N = ReferenceFrame('N')
  231. q1 = dynamicsymbols('q1')
  232. A = N.orientnew('a', 'Axis', [q1, N.z])
  233. #modify default indices:
  234. minds = [x+'1' for x in N.indices]
  235. B = N.orientnew('b', 'Axis', [q1, N.z], indices=minds)
  236. assert N.indices == A.indices
  237. assert B.indices == minds
  238. def test_orientnew_respects_input_latexs():
  239. N = ReferenceFrame('N')
  240. q1 = dynamicsymbols('q1')
  241. A = N.orientnew('a', 'Axis', [q1, N.z])
  242. #build default and alternate latex_vecs:
  243. def_latex_vecs = [(r"\mathbf{\hat{%s}_%s}" % (A.name.lower(),
  244. A.indices[0])), (r"\mathbf{\hat{%s}_%s}" %
  245. (A.name.lower(), A.indices[1])),
  246. (r"\mathbf{\hat{%s}_%s}" % (A.name.lower(),
  247. A.indices[2]))]
  248. name = 'b'
  249. indices = [x+'1' for x in N.indices]
  250. new_latex_vecs = [(r"\mathbf{\hat{%s}_{%s}}" % (name.lower(),
  251. indices[0])), (r"\mathbf{\hat{%s}_{%s}}" %
  252. (name.lower(), indices[1])),
  253. (r"\mathbf{\hat{%s}_{%s}}" % (name.lower(),
  254. indices[2]))]
  255. B = N.orientnew(name, 'Axis', [q1, N.z], latexs=new_latex_vecs)
  256. assert A.latex_vecs == def_latex_vecs
  257. assert B.latex_vecs == new_latex_vecs
  258. assert B.indices != indices
  259. def test_orientnew_respects_input_variables():
  260. N = ReferenceFrame('N')
  261. q1 = dynamicsymbols('q1')
  262. A = N.orientnew('a', 'Axis', [q1, N.z])
  263. #build non-standard variable names
  264. name = 'b'
  265. new_variables = ['notb_'+x+'1' for x in N.indices]
  266. B = N.orientnew(name, 'Axis', [q1, N.z], variables=new_variables)
  267. for j,var in enumerate(A.varlist):
  268. assert var.name == A.name + '_' + A.indices[j]
  269. for j,var in enumerate(B.varlist):
  270. assert var.name == new_variables[j]
  271. def test_issue_10348():
  272. u = dynamicsymbols('u:3')
  273. I = ReferenceFrame('I')
  274. I.orientnew('A', 'space', u, 'XYZ')
  275. def test_issue_11503():
  276. A = ReferenceFrame("A")
  277. A.orientnew("B", "Axis", [35, A.y])
  278. C = ReferenceFrame("C")
  279. A.orient(C, "Axis", [70, C.z])
  280. def test_partial_velocity():
  281. N = ReferenceFrame('N')
  282. A = ReferenceFrame('A')
  283. u1, u2 = dynamicsymbols('u1, u2')
  284. A.set_ang_vel(N, u1 * A.x + u2 * N.y)
  285. assert N.partial_velocity(A, u1) == -A.x
  286. assert N.partial_velocity(A, u1, u2) == (-A.x, -N.y)
  287. assert A.partial_velocity(N, u1) == A.x
  288. assert A.partial_velocity(N, u1, u2) == (A.x, N.y)
  289. assert N.partial_velocity(N, u1) == 0
  290. assert A.partial_velocity(A, u1) == 0
  291. def test_issue_11498():
  292. A = ReferenceFrame('A')
  293. B = ReferenceFrame('B')
  294. # Identity transformation
  295. A.orient(B, 'DCM', eye(3))
  296. assert A.dcm(B) == Matrix([[1, 0, 0], [0, 1, 0], [0, 0, 1]])
  297. assert B.dcm(A) == Matrix([[1, 0, 0], [0, 1, 0], [0, 0, 1]])
  298. # x -> y
  299. # y -> -z
  300. # z -> -x
  301. A.orient(B, 'DCM', Matrix([[0, 1, 0], [0, 0, -1], [-1, 0, 0]]))
  302. assert B.dcm(A) == Matrix([[0, 1, 0], [0, 0, -1], [-1, 0, 0]])
  303. assert A.dcm(B) == Matrix([[0, 0, -1], [1, 0, 0], [0, -1, 0]])
  304. assert B.dcm(A).T == A.dcm(B)
  305. def test_reference_frame():
  306. raises(TypeError, lambda: ReferenceFrame(0))
  307. raises(TypeError, lambda: ReferenceFrame('N', 0))
  308. raises(ValueError, lambda: ReferenceFrame('N', [0, 1]))
  309. raises(TypeError, lambda: ReferenceFrame('N', [0, 1, 2]))
  310. raises(TypeError, lambda: ReferenceFrame('N', ['a', 'b', 'c'], 0))
  311. raises(ValueError, lambda: ReferenceFrame('N', ['a', 'b', 'c'], [0, 1]))
  312. raises(TypeError, lambda: ReferenceFrame('N', ['a', 'b', 'c'], [0, 1, 2]))
  313. raises(TypeError, lambda: ReferenceFrame('N', ['a', 'b', 'c'],
  314. ['a', 'b', 'c'], 0))
  315. raises(ValueError, lambda: ReferenceFrame('N', ['a', 'b', 'c'],
  316. ['a', 'b', 'c'], [0, 1]))
  317. raises(TypeError, lambda: ReferenceFrame('N', ['a', 'b', 'c'],
  318. ['a', 'b', 'c'], [0, 1, 2]))
  319. N = ReferenceFrame('N')
  320. assert N[0] == CoordinateSym('N_x', N, 0)
  321. assert N[1] == CoordinateSym('N_y', N, 1)
  322. assert N[2] == CoordinateSym('N_z', N, 2)
  323. raises(ValueError, lambda: N[3])
  324. N = ReferenceFrame('N', ['a', 'b', 'c'])
  325. assert N['a'] == N.x
  326. assert N['b'] == N.y
  327. assert N['c'] == N.z
  328. raises(ValueError, lambda: N['d'])
  329. assert str(N) == 'N'
  330. A = ReferenceFrame('A')
  331. B = ReferenceFrame('B')
  332. q0, q1, q2, q3 = symbols('q0 q1 q2 q3')
  333. raises(TypeError, lambda: A.orient(B, 'DCM', 0))
  334. raises(TypeError, lambda: B.orient(N, 'Space', [q1, q2, q3], '222'))
  335. raises(TypeError, lambda: B.orient(N, 'Axis', [q1, N.x + 2 * N.y], '222'))
  336. raises(TypeError, lambda: B.orient(N, 'Axis', q1))
  337. raises(IndexError, lambda: B.orient(N, 'Axis', [q1]))
  338. raises(TypeError, lambda: B.orient(N, 'Quaternion', [q0, q1, q2, q3], '222'))
  339. raises(TypeError, lambda: B.orient(N, 'Quaternion', q0))
  340. raises(TypeError, lambda: B.orient(N, 'Quaternion', [q0, q1, q2]))
  341. raises(NotImplementedError, lambda: B.orient(N, 'Foo', [q0, q1, q2]))
  342. raises(TypeError, lambda: B.orient(N, 'Body', [q1, q2], '232'))
  343. raises(TypeError, lambda: B.orient(N, 'Space', [q1, q2], '232'))
  344. N.set_ang_acc(B, 0)
  345. assert N.ang_acc_in(B) == Vector(0)
  346. N.set_ang_vel(B, 0)
  347. assert N.ang_vel_in(B) == Vector(0)
  348. def test_check_frame():
  349. raises(VectorTypeError, lambda: _check_frame(0))
  350. def test_dcm_diff_16824():
  351. # NOTE : This is a regression test for the bug introduced in PR 14758,
  352. # identified in 16824, and solved by PR 16828.
  353. # This is the solution to Problem 2.2 on page 264 in Kane & Lenvinson's
  354. # 1985 book.
  355. q1, q2, q3 = dynamicsymbols('q1:4')
  356. s1 = sin(q1)
  357. c1 = cos(q1)
  358. s2 = sin(q2)
  359. c2 = cos(q2)
  360. s3 = sin(q3)
  361. c3 = cos(q3)
  362. dcm = Matrix([[c2*c3, s1*s2*c3 - s3*c1, c1*s2*c3 + s3*s1],
  363. [c2*s3, s1*s2*s3 + c3*c1, c1*s2*s3 - c3*s1],
  364. [-s2, s1*c2, c1*c2]])
  365. A = ReferenceFrame('A')
  366. B = ReferenceFrame('B')
  367. B.orient(A, 'DCM', dcm)
  368. AwB = B.ang_vel_in(A)
  369. alpha2 = s3*c2*q1.diff() + c3*q2.diff()
  370. beta2 = s1*c2*q3.diff() + c1*q2.diff()
  371. assert simplify(AwB.dot(A.y) - alpha2) == 0
  372. assert simplify(AwB.dot(B.y) - beta2) == 0
  373. def test_orient_explicit():
  374. A = ReferenceFrame('A')
  375. B = ReferenceFrame('B')
  376. A.orient_explicit(B, eye(3))
  377. assert A.dcm(B) == Matrix([[1, 0, 0], [0, 1, 0], [0, 0, 1]])
  378. def test_orient_axis():
  379. A = ReferenceFrame('A')
  380. B = ReferenceFrame('B')
  381. A.orient_axis(B,-B.x, 1)
  382. A1 = A.dcm(B)
  383. A.orient_axis(B, B.x, -1)
  384. A2 = A.dcm(B)
  385. A.orient_axis(B, 1, -B.x)
  386. A3 = A.dcm(B)
  387. assert A1 == A2
  388. assert A2 == A3
  389. raises(TypeError, lambda: A.orient_axis(B, 1, 1))
  390. def test_orient_body():
  391. A = ReferenceFrame('A')
  392. B = ReferenceFrame('B')
  393. B.orient_body_fixed(A, (1,1,0), 'XYX')
  394. assert B.dcm(A) == Matrix([[cos(1), sin(1)**2, -sin(1)*cos(1)], [0, cos(1), sin(1)], [sin(1), -sin(1)*cos(1), cos(1)**2]])
  395. def test_orient_body_advanced():
  396. q1, q2, q3 = dynamicsymbols('q1:4')
  397. c1, c2, c3 = symbols('c1:4')
  398. u1, u2, u3 = dynamicsymbols('q1:4', 1)
  399. # Test with everything as dynamicsymbols
  400. A, B = ReferenceFrame('A'), ReferenceFrame('B')
  401. B.orient_body_fixed(A, (q1, q2, q3), 'zxy')
  402. assert A.dcm(B) == Matrix([
  403. [-sin(q1) * sin(q2) * sin(q3) + cos(q1) * cos(q3), -sin(q1) * cos(q2),
  404. sin(q1) * sin(q2) * cos(q3) + sin(q3) * cos(q1)],
  405. [sin(q1) * cos(q3) + sin(q2) * sin(q3) * cos(q1), cos(q1) * cos(q2),
  406. sin(q1) * sin(q3) - sin(q2) * cos(q1) * cos(q3)],
  407. [-sin(q3) * cos(q2), sin(q2), cos(q2) * cos(q3)]])
  408. assert B.ang_vel_in(A).to_matrix(B) == Matrix([
  409. [-sin(q3) * cos(q2) * u1 + cos(q3) * u2],
  410. [sin(q2) * u1 + u3],
  411. [sin(q3) * u2 + cos(q2) * cos(q3) * u1]])
  412. # Test with constant symbol
  413. A, B = ReferenceFrame('A'), ReferenceFrame('B')
  414. B.orient_body_fixed(A, (q1, c2, q3), 131)
  415. assert A.dcm(B) == Matrix([
  416. [cos(c2), -sin(c2) * cos(q3), sin(c2) * sin(q3)],
  417. [sin(c2) * cos(q1), -sin(q1) * sin(q3) + cos(c2) * cos(q1) * cos(q3),
  418. -sin(q1) * cos(q3) - sin(q3) * cos(c2) * cos(q1)],
  419. [sin(c2) * sin(q1), sin(q1) * cos(c2) * cos(q3) + sin(q3) * cos(q1),
  420. -sin(q1) * sin(q3) * cos(c2) + cos(q1) * cos(q3)]])
  421. assert B.ang_vel_in(A).to_matrix(B) == Matrix([
  422. [cos(c2) * u1 + u3],
  423. [-sin(c2) * cos(q3) * u1],
  424. [sin(c2) * sin(q3) * u1]])
  425. # Test all symbols not time dependent
  426. A, B = ReferenceFrame('A'), ReferenceFrame('B')
  427. B.orient_body_fixed(A, (c1, c2, c3), 123)
  428. assert B.ang_vel_in(A) == Vector(0)
  429. def test_orient_space_advanced():
  430. # space fixed is in the end like body fixed only in opposite order
  431. q1, q2, q3 = dynamicsymbols('q1:4')
  432. c1, c2, c3 = symbols('c1:4')
  433. u1, u2, u3 = dynamicsymbols('q1:4', 1)
  434. # Test with everything as dynamicsymbols
  435. A, B = ReferenceFrame('A'), ReferenceFrame('B')
  436. B.orient_space_fixed(A, (q3, q2, q1), 'yxz')
  437. assert A.dcm(B) == Matrix([
  438. [-sin(q1) * sin(q2) * sin(q3) + cos(q1) * cos(q3), -sin(q1) * cos(q2),
  439. sin(q1) * sin(q2) * cos(q3) + sin(q3) * cos(q1)],
  440. [sin(q1) * cos(q3) + sin(q2) * sin(q3) * cos(q1), cos(q1) * cos(q2),
  441. sin(q1) * sin(q3) - sin(q2) * cos(q1) * cos(q3)],
  442. [-sin(q3) * cos(q2), sin(q2), cos(q2) * cos(q3)]])
  443. assert B.ang_vel_in(A).to_matrix(B) == Matrix([
  444. [-sin(q3) * cos(q2) * u1 + cos(q3) * u2],
  445. [sin(q2) * u1 + u3],
  446. [sin(q3) * u2 + cos(q2) * cos(q3) * u1]])
  447. # Test with constant symbol
  448. A, B = ReferenceFrame('A'), ReferenceFrame('B')
  449. B.orient_space_fixed(A, (q3, c2, q1), 131)
  450. assert A.dcm(B) == Matrix([
  451. [cos(c2), -sin(c2) * cos(q3), sin(c2) * sin(q3)],
  452. [sin(c2) * cos(q1), -sin(q1) * sin(q3) + cos(c2) * cos(q1) * cos(q3),
  453. -sin(q1) * cos(q3) - sin(q3) * cos(c2) * cos(q1)],
  454. [sin(c2) * sin(q1), sin(q1) * cos(c2) * cos(q3) + sin(q3) * cos(q1),
  455. -sin(q1) * sin(q3) * cos(c2) + cos(q1) * cos(q3)]])
  456. assert B.ang_vel_in(A).to_matrix(B) == Matrix([
  457. [cos(c2) * u1 + u3],
  458. [-sin(c2) * cos(q3) * u1],
  459. [sin(c2) * sin(q3) * u1]])
  460. # Test all symbols not time dependent
  461. A, B = ReferenceFrame('A'), ReferenceFrame('B')
  462. B.orient_space_fixed(A, (c1, c2, c3), 123)
  463. assert B.ang_vel_in(A) == Vector(0)
  464. def test_orient_body_simple_ang_vel():
  465. """This test ensures that the simplest form of that linear system solution
  466. is returned, thus the == for the expression comparison."""
  467. psi, theta, phi = dynamicsymbols('psi, theta, varphi')
  468. t = dynamicsymbols._t
  469. A = ReferenceFrame('A')
  470. B = ReferenceFrame('B')
  471. B.orient_body_fixed(A, (psi, theta, phi), 'ZXZ')
  472. A_w_B = B.ang_vel_in(A)
  473. assert A_w_B.args[0][1] == B
  474. assert A_w_B.args[0][0][0] == (sin(theta)*sin(phi)*psi.diff(t) +
  475. cos(phi)*theta.diff(t))
  476. assert A_w_B.args[0][0][1] == (sin(theta)*cos(phi)*psi.diff(t) -
  477. sin(phi)*theta.diff(t))
  478. assert A_w_B.args[0][0][2] == cos(theta)*psi.diff(t) + phi.diff(t)
  479. def test_orient_space():
  480. A = ReferenceFrame('A')
  481. B = ReferenceFrame('B')
  482. B.orient_space_fixed(A, (0,0,0), '123')
  483. assert B.dcm(A) == Matrix([[1, 0, 0], [0, 1, 0], [0, 0, 1]])
  484. def test_orient_quaternion():
  485. A = ReferenceFrame('A')
  486. B = ReferenceFrame('B')
  487. B.orient_quaternion(A, (0,0,0,0))
  488. assert B.dcm(A) == Matrix([[0, 0, 0], [0, 0, 0], [0, 0, 0]])
  489. def test_looped_frame_warning():
  490. A = ReferenceFrame('A')
  491. B = ReferenceFrame('B')
  492. C = ReferenceFrame('C')
  493. a, b, c = symbols('a b c')
  494. B.orient_axis(A, A.x, a)
  495. C.orient_axis(B, B.x, b)
  496. with warnings.catch_warnings(record = True) as w:
  497. warnings.simplefilter("always")
  498. A.orient_axis(C, C.x, c)
  499. assert issubclass(w[-1].category, UserWarning)
  500. assert 'Loops are defined among the orientation of frames. ' + \
  501. 'This is likely not desired and may cause errors in your calculations.' in str(w[-1].message)
  502. def test_frame_dict():
  503. A = ReferenceFrame('A')
  504. B = ReferenceFrame('B')
  505. C = ReferenceFrame('C')
  506. a, b, c = symbols('a b c')
  507. B.orient_axis(A, A.x, a)
  508. assert A._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(a), -sin(a)],[0, sin(a), cos(a)]])}
  509. assert B._dcm_dict == {A: Matrix([[1, 0, 0],[0, cos(a), sin(a)],[0, -sin(a), cos(a)]])}
  510. assert C._dcm_dict == {}
  511. B.orient_axis(C, C.x, b)
  512. # Previous relation is not wiped
  513. assert A._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(a), -sin(a)],[0, sin(a), cos(a)]])}
  514. assert B._dcm_dict == {A: Matrix([[1, 0, 0],[0, cos(a), sin(a)],[0, -sin(a), cos(a)]]), \
  515. C: Matrix([[1, 0, 0],[0, cos(b), sin(b)],[0, -sin(b), cos(b)]])}
  516. assert C._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(b), -sin(b)],[0, sin(b), cos(b)]])}
  517. A.orient_axis(B, B.x, c)
  518. # Previous relation is updated
  519. assert B._dcm_dict == {C: Matrix([[1, 0, 0],[0, cos(b), sin(b)],[0, -sin(b), cos(b)]]),\
  520. A: Matrix([[1, 0, 0],[0, cos(c), -sin(c)],[0, sin(c), cos(c)]])}
  521. assert A._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(c), sin(c)],[0, -sin(c), cos(c)]])}
  522. assert C._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(b), -sin(b)],[0, sin(b), cos(b)]])}
  523. def test_dcm_cache_dict():
  524. A = ReferenceFrame('A')
  525. B = ReferenceFrame('B')
  526. C = ReferenceFrame('C')
  527. D = ReferenceFrame('D')
  528. a, b, c = symbols('a b c')
  529. B.orient_axis(A, A.x, a)
  530. C.orient_axis(B, B.x, b)
  531. D.orient_axis(C, C.x, c)
  532. assert D._dcm_dict == {C: Matrix([[1, 0, 0],[0, cos(c), sin(c)],[0, -sin(c), cos(c)]])}
  533. assert C._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(b), sin(b)],[0, -sin(b), cos(b)]]), \
  534. D: Matrix([[1, 0, 0],[0, cos(c), -sin(c)],[0, sin(c), cos(c)]])}
  535. assert B._dcm_dict == {A: Matrix([[1, 0, 0],[0, cos(a), sin(a)],[0, -sin(a), cos(a)]]), \
  536. C: Matrix([[1, 0, 0],[0, cos(b), -sin(b)],[0, sin(b), cos(b)]])}
  537. assert A._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(a), -sin(a)],[0, sin(a), cos(a)]])}
  538. assert D._dcm_dict == D._dcm_cache
  539. D.dcm(A) # Check calculated dcm relation is stored in _dcm_cache and not in _dcm_dict
  540. assert list(A._dcm_cache.keys()) == [A, B, D]
  541. assert list(D._dcm_cache.keys()) == [C, A]
  542. assert list(A._dcm_dict.keys()) == [B]
  543. assert list(D._dcm_dict.keys()) == [C]
  544. assert A._dcm_dict != A._dcm_cache
  545. A.orient_axis(B, B.x, b) # _dcm_cache of A is wiped out and new relation is stored.
  546. assert A._dcm_dict == {B: Matrix([[1, 0, 0],[0, cos(b), sin(b)],[0, -sin(b), cos(b)]])}
  547. assert A._dcm_dict == A._dcm_cache
  548. assert B._dcm_dict == {C: Matrix([[1, 0, 0],[0, cos(b), -sin(b)],[0, sin(b), cos(b)]]), \
  549. A: Matrix([[1, 0, 0],[0, cos(b), -sin(b)],[0, sin(b), cos(b)]])}