图片解析应用
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.

73 lines
2.2 KiB

  1. """
  2. This module implements the Residue function and related tools for working
  3. with residues.
  4. """
  5. from sympy.core.mul import Mul
  6. from sympy.core.singleton import S
  7. from sympy.core.sympify import sympify
  8. from sympy.utilities.timeutils import timethis
  9. @timethis('residue')
  10. def residue(expr, x, x0):
  11. """
  12. Finds the residue of ``expr`` at the point x=x0.
  13. The residue is defined as the coefficient of ``1/(x-x0)`` in the power series
  14. expansion about ``x=x0``.
  15. Examples
  16. ========
  17. >>> from sympy import Symbol, residue, sin
  18. >>> x = Symbol("x")
  19. >>> residue(1/x, x, 0)
  20. 1
  21. >>> residue(1/x**2, x, 0)
  22. 0
  23. >>> residue(2/sin(x), x, 0)
  24. 2
  25. This function is essential for the Residue Theorem [1].
  26. References
  27. ==========
  28. .. [1] https://en.wikipedia.org/wiki/Residue_theorem
  29. """
  30. # The current implementation uses series expansion to
  31. # calculate it. A more general implementation is explained in
  32. # the section 5.6 of the Bronstein's book {M. Bronstein:
  33. # Symbolic Integration I, Springer Verlag (2005)}. For purely
  34. # rational functions, the algorithm is much easier. See
  35. # sections 2.4, 2.5, and 2.7 (this section actually gives an
  36. # algorithm for computing any Laurent series coefficient for
  37. # a rational function). The theory in section 2.4 will help to
  38. # understand why the resultant works in the general algorithm.
  39. # For the definition of a resultant, see section 1.4 (and any
  40. # previous sections for more review).
  41. from sympy.series.order import Order
  42. from sympy.simplify.radsimp import collect
  43. expr = sympify(expr)
  44. if x0 != 0:
  45. expr = expr.subs(x, x + x0)
  46. for n in (0, 1, 2, 4, 8, 16, 32):
  47. s = expr.nseries(x, n=n)
  48. if not s.has(Order) or s.getn() >= 0:
  49. break
  50. s = collect(s.removeO(), x)
  51. if s.is_Add:
  52. args = s.args
  53. else:
  54. args = [s]
  55. res = S.Zero
  56. for arg in args:
  57. c, m = arg.as_coeff_mul(x)
  58. m = Mul(*m)
  59. if not (m in (S.One, x) or (m.is_Pow and m.exp.is_Integer)):
  60. raise NotImplementedError('term of unexpected form: %s' % m)
  61. if m == 1/x:
  62. res += c
  63. return res