图片解析应用
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

335 lines
11 KiB

  1. from sympy.core import cacheit, Dummy, Ne, Integer, Rational, S, Wild
  2. from sympy.functions import binomial, sin, cos, Piecewise
  3. from .integrals import integrate
  4. # TODO sin(a*x)*cos(b*x) -> sin((a+b)x) + sin((a-b)x) ?
  5. # creating, each time, Wild's and sin/cos/Mul is expensive. Also, our match &
  6. # subs are very slow when not cached, and if we create Wild each time, we
  7. # effectively block caching.
  8. #
  9. # so we cache the pattern
  10. # need to use a function instead of lamda since hash of lambda changes on
  11. # each call to _pat_sincos
  12. def _integer_instance(n):
  13. return isinstance(n, Integer)
  14. @cacheit
  15. def _pat_sincos(x):
  16. a = Wild('a', exclude=[x])
  17. n, m = [Wild(s, exclude=[x], properties=[_integer_instance])
  18. for s in 'nm']
  19. pat = sin(a*x)**n * cos(a*x)**m
  20. return pat, a, n, m
  21. _u = Dummy('u')
  22. def trigintegrate(f, x, conds='piecewise'):
  23. """
  24. Integrate f = Mul(trig) over x.
  25. Examples
  26. ========
  27. >>> from sympy import sin, cos, tan, sec
  28. >>> from sympy.integrals.trigonometry import trigintegrate
  29. >>> from sympy.abc import x
  30. >>> trigintegrate(sin(x)*cos(x), x)
  31. sin(x)**2/2
  32. >>> trigintegrate(sin(x)**2, x)
  33. x/2 - sin(x)*cos(x)/2
  34. >>> trigintegrate(tan(x)*sec(x), x)
  35. 1/cos(x)
  36. >>> trigintegrate(sin(x)*tan(x), x)
  37. -log(sin(x) - 1)/2 + log(sin(x) + 1)/2 - sin(x)
  38. References
  39. ==========
  40. .. [1] http://en.wikibooks.org/wiki/Calculus/Integration_techniques
  41. See Also
  42. ========
  43. sympy.integrals.integrals.Integral.doit
  44. sympy.integrals.integrals.Integral
  45. """
  46. pat, a, n, m = _pat_sincos(x)
  47. f = f.rewrite('sincos')
  48. M = f.match(pat)
  49. if M is None:
  50. return
  51. n, m = M[n], M[m]
  52. if n.is_zero and m.is_zero:
  53. return x
  54. zz = x if n.is_zero else S.Zero
  55. a = M[a]
  56. if n.is_odd or m.is_odd:
  57. u = _u
  58. n_, m_ = n.is_odd, m.is_odd
  59. # take smallest n or m -- to choose simplest substitution
  60. if n_ and m_:
  61. # Make sure to choose the positive one
  62. # otherwise an incorrect integral can occur.
  63. if n < 0 and m > 0:
  64. m_ = True
  65. n_ = False
  66. elif m < 0 and n > 0:
  67. n_ = True
  68. m_ = False
  69. # Both are negative so choose the smallest n or m
  70. # in absolute value for simplest substitution.
  71. elif (n < 0 and m < 0):
  72. n_ = n > m
  73. m_ = not (n > m)
  74. # Both n and m are odd and positive
  75. else:
  76. n_ = (n < m) # NB: careful here, one of the
  77. m_ = not (n < m) # conditions *must* be true
  78. # n m u=C (n-1)/2 m
  79. # S(x) * C(x) dx --> -(1-u^2) * u du
  80. if n_:
  81. ff = -(1 - u**2)**((n - 1)/2) * u**m
  82. uu = cos(a*x)
  83. # n m u=S n (m-1)/2
  84. # S(x) * C(x) dx --> u * (1-u^2) du
  85. elif m_:
  86. ff = u**n * (1 - u**2)**((m - 1)/2)
  87. uu = sin(a*x)
  88. fi = integrate(ff, u) # XXX cyclic deps
  89. fx = fi.subs(u, uu)
  90. if conds == 'piecewise':
  91. return Piecewise((fx / a, Ne(a, 0)), (zz, True))
  92. return fx / a
  93. # n & m are both even
  94. #
  95. # 2k 2m 2l 2l
  96. # we transform S (x) * C (x) into terms with only S (x) or C (x)
  97. #
  98. # example:
  99. # 100 4 100 2 2 100 4 2
  100. # S (x) * C (x) = S (x) * (1-S (x)) = S (x) * (1 + S (x) - 2*S (x))
  101. #
  102. # 104 102 100
  103. # = S (x) - 2*S (x) + S (x)
  104. # 2k
  105. # then S is integrated with recursive formula
  106. # take largest n or m -- to choose simplest substitution
  107. n_ = (abs(n) > abs(m))
  108. m_ = (abs(m) > abs(n))
  109. res = S.Zero
  110. if n_:
  111. # 2k 2 k i 2i
  112. # C = (1 - S ) = sum(i, (-) * B(k, i) * S )
  113. if m > 0:
  114. for i in range(0, m//2 + 1):
  115. res += (S.NegativeOne**i * binomial(m//2, i) *
  116. _sin_pow_integrate(n + 2*i, x))
  117. elif m == 0:
  118. res = _sin_pow_integrate(n, x)
  119. else:
  120. # m < 0 , |n| > |m|
  121. # /
  122. # |
  123. # | m n
  124. # | cos (x) sin (x) dx =
  125. # |
  126. # |
  127. #/
  128. # /
  129. # |
  130. # -1 m+1 n-1 n - 1 | m+2 n-2
  131. # ________ cos (x) sin (x) + _______ | cos (x) sin (x) dx
  132. # |
  133. # m + 1 m + 1 |
  134. # /
  135. res = (Rational(-1, m + 1) * cos(x)**(m + 1) * sin(x)**(n - 1) +
  136. Rational(n - 1, m + 1) *
  137. trigintegrate(cos(x)**(m + 2)*sin(x)**(n - 2), x))
  138. elif m_:
  139. # 2k 2 k i 2i
  140. # S = (1 - C ) = sum(i, (-) * B(k, i) * C )
  141. if n > 0:
  142. # / /
  143. # | |
  144. # | m n | -m n
  145. # | cos (x)*sin (x) dx or | cos (x) * sin (x) dx
  146. # | |
  147. # / /
  148. #
  149. # |m| > |n| ; m, n >0 ; m, n belong to Z - {0}
  150. # n 2
  151. # sin (x) term is expanded here in terms of cos (x),
  152. # and then integrated.
  153. #
  154. for i in range(0, n//2 + 1):
  155. res += (S.NegativeOne**i * binomial(n//2, i) *
  156. _cos_pow_integrate(m + 2*i, x))
  157. elif n == 0:
  158. # /
  159. # |
  160. # | 1
  161. # | _ _ _
  162. # | m
  163. # | cos (x)
  164. # /
  165. #
  166. res = _cos_pow_integrate(m, x)
  167. else:
  168. # n < 0 , |m| > |n|
  169. # /
  170. # |
  171. # | m n
  172. # | cos (x) sin (x) dx =
  173. # |
  174. # |
  175. #/
  176. # /
  177. # |
  178. # 1 m-1 n+1 m - 1 | m-2 n+2
  179. # _______ cos (x) sin (x) + _______ | cos (x) sin (x) dx
  180. # |
  181. # n + 1 n + 1 |
  182. # /
  183. res = (Rational(1, n + 1) * cos(x)**(m - 1)*sin(x)**(n + 1) +
  184. Rational(m - 1, n + 1) *
  185. trigintegrate(cos(x)**(m - 2)*sin(x)**(n + 2), x))
  186. else:
  187. if m == n:
  188. ##Substitute sin(2x)/2 for sin(x)cos(x) and then Integrate.
  189. res = integrate((sin(2*x)*S.Half)**m, x)
  190. elif (m == -n):
  191. if n < 0:
  192. # Same as the scheme described above.
  193. # the function argument to integrate in the end will
  194. # be 1, this cannot be integrated by trigintegrate.
  195. # Hence use sympy.integrals.integrate.
  196. res = (Rational(1, n + 1) * cos(x)**(m - 1) * sin(x)**(n + 1) +
  197. Rational(m - 1, n + 1) *
  198. integrate(cos(x)**(m - 2) * sin(x)**(n + 2), x))
  199. else:
  200. res = (Rational(-1, m + 1) * cos(x)**(m + 1) * sin(x)**(n - 1) +
  201. Rational(n - 1, m + 1) *
  202. integrate(cos(x)**(m + 2)*sin(x)**(n - 2), x))
  203. if conds == 'piecewise':
  204. return Piecewise((res.subs(x, a*x) / a, Ne(a, 0)), (zz, True))
  205. return res.subs(x, a*x) / a
  206. def _sin_pow_integrate(n, x):
  207. if n > 0:
  208. if n == 1:
  209. #Recursion break
  210. return -cos(x)
  211. # n > 0
  212. # / /
  213. # | |
  214. # | n -1 n-1 n - 1 | n-2
  215. # | sin (x) dx = ______ cos (x) sin (x) + _______ | sin (x) dx
  216. # | |
  217. # | n n |
  218. #/ /
  219. #
  220. #
  221. return (Rational(-1, n) * cos(x) * sin(x)**(n - 1) +
  222. Rational(n - 1, n) * _sin_pow_integrate(n - 2, x))
  223. if n < 0:
  224. if n == -1:
  225. ##Make sure this does not come back here again.
  226. ##Recursion breaks here or at n==0.
  227. return trigintegrate(1/sin(x), x)
  228. # n < 0
  229. # / /
  230. # | |
  231. # | n 1 n+1 n + 2 | n+2
  232. # | sin (x) dx = _______ cos (x) sin (x) + _______ | sin (x) dx
  233. # | |
  234. # | n + 1 n + 1 |
  235. #/ /
  236. #
  237. return (Rational(1, n + 1) * cos(x) * sin(x)**(n + 1) +
  238. Rational(n + 2, n + 1) * _sin_pow_integrate(n + 2, x))
  239. else:
  240. #n == 0
  241. #Recursion break.
  242. return x
  243. def _cos_pow_integrate(n, x):
  244. if n > 0:
  245. if n == 1:
  246. #Recursion break.
  247. return sin(x)
  248. # n > 0
  249. # / /
  250. # | |
  251. # | n 1 n-1 n - 1 | n-2
  252. # | sin (x) dx = ______ sin (x) cos (x) + _______ | cos (x) dx
  253. # | |
  254. # | n n |
  255. #/ /
  256. #
  257. return (Rational(1, n) * sin(x) * cos(x)**(n - 1) +
  258. Rational(n - 1, n) * _cos_pow_integrate(n - 2, x))
  259. if n < 0:
  260. if n == -1:
  261. ##Recursion break
  262. return trigintegrate(1/cos(x), x)
  263. # n < 0
  264. # / /
  265. # | |
  266. # | n -1 n+1 n + 2 | n+2
  267. # | cos (x) dx = _______ sin (x) cos (x) + _______ | cos (x) dx
  268. # | |
  269. # | n + 1 n + 1 |
  270. #/ /
  271. #
  272. return (Rational(-1, n + 1) * sin(x) * cos(x)**(n + 1) +
  273. Rational(n + 2, n + 1) * _cos_pow_integrate(n + 2, x))
  274. else:
  275. # n == 0
  276. #Recursion Break.
  277. return x