test_utils.py 8.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202
  1. from sympy.core.numbers import comp, Rational
  2. from sympy.physics.optics.utils import (refraction_angle, fresnel_coefficients,
  3. deviation, brewster_angle, critical_angle, lens_makers_formula,
  4. mirror_formula, lens_formula, hyperfocal_distance,
  5. transverse_magnification)
  6. from sympy.physics.optics.medium import Medium
  7. from sympy.physics.units import e0
  8. from sympy.core.numbers import oo
  9. from sympy.core.symbol import symbols
  10. from sympy.functions.elementary.miscellaneous import sqrt
  11. from sympy.matrices.dense import Matrix
  12. from sympy.geometry.point import Point3D
  13. from sympy.geometry.line import Ray3D
  14. from sympy.geometry.plane import Plane
  15. from sympy.testing.pytest import raises
  16. ae = lambda a, b, n: comp(a, b, 10**-n)
  17. def test_refraction_angle():
  18. n1, n2 = symbols('n1, n2')
  19. m1 = Medium('m1')
  20. m2 = Medium('m2')
  21. r1 = Ray3D(Point3D(-1, -1, 1), Point3D(0, 0, 0))
  22. i = Matrix([1, 1, 1])
  23. n = Matrix([0, 0, 1])
  24. normal_ray = Ray3D(Point3D(0, 0, 0), Point3D(0, 0, 1))
  25. P = Plane(Point3D(0, 0, 0), normal_vector=[0, 0, 1])
  26. assert refraction_angle(r1, 1, 1, n) == Matrix([
  27. [ 1],
  28. [ 1],
  29. [-1]])
  30. assert refraction_angle([1, 1, 1], 1, 1, n) == Matrix([
  31. [ 1],
  32. [ 1],
  33. [-1]])
  34. assert refraction_angle((1, 1, 1), 1, 1, n) == Matrix([
  35. [ 1],
  36. [ 1],
  37. [-1]])
  38. assert refraction_angle(i, 1, 1, [0, 0, 1]) == Matrix([
  39. [ 1],
  40. [ 1],
  41. [-1]])
  42. assert refraction_angle(i, 1, 1, (0, 0, 1)) == Matrix([
  43. [ 1],
  44. [ 1],
  45. [-1]])
  46. assert refraction_angle(i, 1, 1, normal_ray) == Matrix([
  47. [ 1],
  48. [ 1],
  49. [-1]])
  50. assert refraction_angle(i, 1, 1, plane=P) == Matrix([
  51. [ 1],
  52. [ 1],
  53. [-1]])
  54. assert refraction_angle(r1, 1, 1, plane=P) == \
  55. Ray3D(Point3D(0, 0, 0), Point3D(1, 1, -1))
  56. assert refraction_angle(r1, m1, 1.33, plane=P) == \
  57. Ray3D(Point3D(0, 0, 0), Point3D(Rational(100, 133), Rational(100, 133), -789378201649271*sqrt(3)/1000000000000000))
  58. assert refraction_angle(r1, 1, m2, plane=P) == \
  59. Ray3D(Point3D(0, 0, 0), Point3D(1, 1, -1))
  60. assert refraction_angle(r1, n1, n2, plane=P) == \
  61. Ray3D(Point3D(0, 0, 0), Point3D(n1/n2, n1/n2, -sqrt(3)*sqrt(-2*n1**2/(3*n2**2) + 1)))
  62. assert refraction_angle(r1, 1.33, 1, plane=P) == 0 # TIR
  63. assert refraction_angle(r1, 1, 1, normal_ray) == \
  64. Ray3D(Point3D(0, 0, 0), direction_ratio=[1, 1, -1])
  65. assert ae(refraction_angle(0.5, 1, 2), 0.24207, 5)
  66. assert ae(refraction_angle(0.5, 2, 1), 1.28293, 5)
  67. raises(ValueError, lambda: refraction_angle(r1, m1, m2, normal_ray, P))
  68. raises(TypeError, lambda: refraction_angle(m1, m1, m2)) # can add other values for arg[0]
  69. raises(TypeError, lambda: refraction_angle(r1, m1, m2, None, i))
  70. raises(TypeError, lambda: refraction_angle(r1, m1, m2, m2))
  71. def test_fresnel_coefficients():
  72. assert all(ae(i, j, 5) for i, j in zip(
  73. fresnel_coefficients(0.5, 1, 1.33),
  74. [0.11163, -0.17138, 0.83581, 0.82862]))
  75. assert all(ae(i, j, 5) for i, j in zip(
  76. fresnel_coefficients(0.5, 1.33, 1),
  77. [-0.07726, 0.20482, 1.22724, 1.20482]))
  78. m1 = Medium('m1')
  79. m2 = Medium('m2', n=2)
  80. assert all(ae(i, j, 5) for i, j in zip(
  81. fresnel_coefficients(0.3, m1, m2),
  82. [0.31784, -0.34865, 0.65892, 0.65135]))
  83. ans = [[-0.23563, -0.97184], [0.81648, -0.57738]]
  84. got = fresnel_coefficients(0.6, m2, m1)
  85. for i, j in zip(got, ans):
  86. for a, b in zip(i.as_real_imag(), j):
  87. assert ae(a, b, 5)
  88. def test_deviation():
  89. n1, n2 = symbols('n1, n2')
  90. r1 = Ray3D(Point3D(-1, -1, 1), Point3D(0, 0, 0))
  91. n = Matrix([0, 0, 1])
  92. i = Matrix([-1, -1, -1])
  93. normal_ray = Ray3D(Point3D(0, 0, 0), Point3D(0, 0, 1))
  94. P = Plane(Point3D(0, 0, 0), normal_vector=[0, 0, 1])
  95. assert deviation(r1, 1, 1, normal=n) == 0
  96. assert deviation(r1, 1, 1, plane=P) == 0
  97. assert deviation(r1, 1, 1.1, plane=P).evalf(3) + 0.119 < 1e-3
  98. assert deviation(i, 1, 1.1, normal=normal_ray).evalf(3) + 0.119 < 1e-3
  99. assert deviation(r1, 1.33, 1, plane=P) is None # TIR
  100. assert deviation(r1, 1, 1, normal=[0, 0, 1]) == 0
  101. assert deviation([-1, -1, -1], 1, 1, normal=[0, 0, 1]) == 0
  102. assert ae(deviation(0.5, 1, 2), -0.25793, 5)
  103. assert ae(deviation(0.5, 2, 1), 0.78293, 5)
  104. def test_brewster_angle():
  105. m1 = Medium('m1', n=1)
  106. m2 = Medium('m2', n=1.33)
  107. assert ae(brewster_angle(m1, m2), 0.93, 2)
  108. m1 = Medium('m1', permittivity=e0, n=1)
  109. m2 = Medium('m2', permittivity=e0, n=1.33)
  110. assert ae(brewster_angle(m1, m2), 0.93, 2)
  111. assert ae(brewster_angle(1, 1.33), 0.93, 2)
  112. def test_critical_angle():
  113. m1 = Medium('m1', n=1)
  114. m2 = Medium('m2', n=1.33)
  115. assert ae(critical_angle(m2, m1), 0.85, 2)
  116. def test_lens_makers_formula():
  117. n1, n2 = symbols('n1, n2')
  118. m1 = Medium('m1', permittivity=e0, n=1)
  119. m2 = Medium('m2', permittivity=e0, n=1.33)
  120. assert lens_makers_formula(n1, n2, 10, -10) == 5.0*n2/(n1 - n2)
  121. assert ae(lens_makers_formula(m1, m2, 10, -10), -20.15, 2)
  122. assert ae(lens_makers_formula(1.33, 1, 10, -10), 15.15, 2)
  123. def test_mirror_formula():
  124. u, v, f = symbols('u, v, f')
  125. assert mirror_formula(focal_length=f, u=u) == f*u/(-f + u)
  126. assert mirror_formula(focal_length=f, v=v) == f*v/(-f + v)
  127. assert mirror_formula(u=u, v=v) == u*v/(u + v)
  128. assert mirror_formula(u=oo, v=v) == v
  129. assert mirror_formula(u=oo, v=oo) is oo
  130. assert mirror_formula(focal_length=oo, u=u) == -u
  131. assert mirror_formula(u=u, v=oo) == u
  132. assert mirror_formula(focal_length=oo, v=oo) is oo
  133. assert mirror_formula(focal_length=f, v=oo) == f
  134. assert mirror_formula(focal_length=oo, v=v) == -v
  135. assert mirror_formula(focal_length=oo, u=oo) is oo
  136. assert mirror_formula(focal_length=f, u=oo) == f
  137. assert mirror_formula(focal_length=oo, u=u) == -u
  138. raises(ValueError, lambda: mirror_formula(focal_length=f, u=u, v=v))
  139. def test_lens_formula():
  140. u, v, f = symbols('u, v, f')
  141. assert lens_formula(focal_length=f, u=u) == f*u/(f + u)
  142. assert lens_formula(focal_length=f, v=v) == f*v/(f - v)
  143. assert lens_formula(u=u, v=v) == u*v/(u - v)
  144. assert lens_formula(u=oo, v=v) == v
  145. assert lens_formula(u=oo, v=oo) is oo
  146. assert lens_formula(focal_length=oo, u=u) == u
  147. assert lens_formula(u=u, v=oo) == -u
  148. assert lens_formula(focal_length=oo, v=oo) is -oo
  149. assert lens_formula(focal_length=oo, v=v) == v
  150. assert lens_formula(focal_length=f, v=oo) == -f
  151. assert lens_formula(focal_length=oo, u=oo) is oo
  152. assert lens_formula(focal_length=oo, u=u) == u
  153. assert lens_formula(focal_length=f, u=oo) == f
  154. raises(ValueError, lambda: lens_formula(focal_length=f, u=u, v=v))
  155. def test_hyperfocal_distance():
  156. f, N, c = symbols('f, N, c')
  157. assert hyperfocal_distance(f=f, N=N, c=c) == f**2/(N*c)
  158. assert ae(hyperfocal_distance(f=0.5, N=8, c=0.0033), 9.47, 2)
  159. def test_transverse_magnification():
  160. si, so = symbols('si, so')
  161. assert transverse_magnification(si, so) == -si/so
  162. assert transverse_magnification(30, 15) == -2
  163. def test_lens_makers_formula_thick_lens():
  164. n1, n2 = symbols('n1, n2')
  165. m1 = Medium('m1', permittivity=e0, n=1)
  166. m2 = Medium('m2', permittivity=e0, n=1.33)
  167. assert ae(lens_makers_formula(m1, m2, 10, -10, d=1), -19.82, 2)
  168. assert lens_makers_formula(n1, n2, 1, -1, d=0.1) == n2/((2.0 - (0.1*n1 - 0.1*n2)/n1)*(n1 - n2))
  169. def test_lens_makers_formula_plano_lens():
  170. n1, n2 = symbols('n1, n2')
  171. m1 = Medium('m1', permittivity=e0, n=1)
  172. m2 = Medium('m2', permittivity=e0, n=1.33)
  173. assert ae(lens_makers_formula(m1, m2, 10, oo), -40.30, 2)
  174. assert lens_makers_formula(n1, n2, 10, oo) == 10.0*n2/(n1 - n2)