Coverage for /usr/lib/python3/dist-packages/sympy/integrals/integrals.py: 9%
635 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-14 15:55 +0200
1from typing import Tuple as tTuple
3from sympy.concrete.expr_with_limits import AddWithLimits
4from sympy.core.add import Add
5from sympy.core.basic import Basic
6from sympy.core.containers import Tuple
7from sympy.core.expr import Expr
8from sympy.core.exprtools import factor_terms
9from sympy.core.function import diff
10from sympy.core.logic import fuzzy_bool
11from sympy.core.mul import Mul
12from sympy.core.numbers import oo, pi
13from sympy.core.relational import Ne
14from sympy.core.singleton import S
15from sympy.core.symbol import (Dummy, Symbol, Wild)
16from sympy.core.sympify import sympify
17from sympy.functions import Piecewise, sqrt, piecewise_fold, tan, cot, atan
18from sympy.functions.elementary.exponential import log
19from sympy.functions.elementary.integers import floor
20from sympy.functions.elementary.complexes import Abs, sign
21from sympy.functions.elementary.miscellaneous import Min, Max
22from .rationaltools import ratint
23from sympy.matrices import MatrixBase
24from sympy.polys import Poly, PolynomialError
25from sympy.series.formal import FormalPowerSeries
26from sympy.series.limits import limit
27from sympy.series.order import Order
28from sympy.tensor.functions import shape
29from sympy.utilities.exceptions import sympy_deprecation_warning
30from sympy.utilities.iterables import is_sequence
31from sympy.utilities.misc import filldedent
34class Integral(AddWithLimits):
35 """Represents unevaluated integral."""
37 __slots__ = ()
39 args: tTuple[Expr, Tuple]
41 def __new__(cls, function, *symbols, **assumptions):
42 """Create an unevaluated integral.
44 Explanation
45 ===========
47 Arguments are an integrand followed by one or more limits.
49 If no limits are given and there is only one free symbol in the
50 expression, that symbol will be used, otherwise an error will be
51 raised.
53 >>> from sympy import Integral
54 >>> from sympy.abc import x, y
55 >>> Integral(x)
56 Integral(x, x)
57 >>> Integral(y)
58 Integral(y, y)
60 When limits are provided, they are interpreted as follows (using
61 ``x`` as though it were the variable of integration):
63 (x,) or x - indefinite integral
64 (x, a) - "evaluate at" integral is an abstract antiderivative
65 (x, a, b) - definite integral
67 The ``as_dummy`` method can be used to see which symbols cannot be
68 targeted by subs: those with a prepended underscore cannot be
69 changed with ``subs``. (Also, the integration variables themselves --
70 the first element of a limit -- can never be changed by subs.)
72 >>> i = Integral(x, x)
73 >>> at = Integral(x, (x, x))
74 >>> i.as_dummy()
75 Integral(x, x)
76 >>> at.as_dummy()
77 Integral(_0, (_0, x))
79 """
81 #This will help other classes define their own definitions
82 #of behaviour with Integral.
83 if hasattr(function, '_eval_Integral'):
84 return function._eval_Integral(*symbols, **assumptions)
86 if isinstance(function, Poly):
87 sympy_deprecation_warning(
88 """
89 integrate(Poly) and Integral(Poly) are deprecated. Instead,
90 use the Poly.integrate() method, or convert the Poly to an
91 Expr first with the Poly.as_expr() method.
92 """,
93 deprecated_since_version="1.6",
94 active_deprecations_target="deprecated-integrate-poly")
96 obj = AddWithLimits.__new__(cls, function, *symbols, **assumptions)
97 return obj
99 def __getnewargs__(self):
100 return (self.function,) + tuple([tuple(xab) for xab in self.limits])
102 @property
103 def free_symbols(self):
104 """
105 This method returns the symbols that will exist when the
106 integral is evaluated. This is useful if one is trying to
107 determine whether an integral depends on a certain
108 symbol or not.
110 Examples
111 ========
113 >>> from sympy import Integral
114 >>> from sympy.abc import x, y
115 >>> Integral(x, (x, y, 1)).free_symbols
116 {y}
118 See Also
119 ========
121 sympy.concrete.expr_with_limits.ExprWithLimits.function
122 sympy.concrete.expr_with_limits.ExprWithLimits.limits
123 sympy.concrete.expr_with_limits.ExprWithLimits.variables
124 """
125 return super().free_symbols
127 def _eval_is_zero(self):
128 # This is a very naive and quick test, not intended to do the integral to
129 # answer whether it is zero or not, e.g. Integral(sin(x), (x, 0, 2*pi))
130 # is zero but this routine should return None for that case. But, like
131 # Mul, there are trivial situations for which the integral will be
132 # zero so we check for those.
133 if self.function.is_zero:
134 return True
135 got_none = False
136 for l in self.limits:
137 if len(l) == 3:
138 z = (l[1] == l[2]) or (l[1] - l[2]).is_zero
139 if z:
140 return True
141 elif z is None:
142 got_none = True
143 free = self.function.free_symbols
144 for xab in self.limits:
145 if len(xab) == 1:
146 free.add(xab[0])
147 continue
148 if len(xab) == 2 and xab[0] not in free:
149 if xab[1].is_zero:
150 return True
151 elif xab[1].is_zero is None:
152 got_none = True
153 # take integration symbol out of free since it will be replaced
154 # with the free symbols in the limits
155 free.discard(xab[0])
156 # add in the new symbols
157 for i in xab[1:]:
158 free.update(i.free_symbols)
159 if self.function.is_zero is False and got_none is False:
160 return False
162 def transform(self, x, u):
163 r"""
164 Performs a change of variables from `x` to `u` using the relationship
165 given by `x` and `u` which will define the transformations `f` and `F`
166 (which are inverses of each other) as follows:
168 1) If `x` is a Symbol (which is a variable of integration) then `u`
169 will be interpreted as some function, f(u), with inverse F(u).
170 This, in effect, just makes the substitution of x with f(x).
172 2) If `u` is a Symbol then `x` will be interpreted as some function,
173 F(x), with inverse f(u). This is commonly referred to as
174 u-substitution.
176 Once f and F have been identified, the transformation is made as
177 follows:
179 .. math:: \int_a^b x \mathrm{d}x \rightarrow \int_{F(a)}^{F(b)} f(x)
180 \frac{\mathrm{d}}{\mathrm{d}x}
182 where `F(x)` is the inverse of `f(x)` and the limits and integrand have
183 been corrected so as to retain the same value after integration.
185 Notes
186 =====
188 The mappings, F(x) or f(u), must lead to a unique integral. Linear
189 or rational linear expression, ``2*x``, ``1/x`` and ``sqrt(x)``, will
190 always work; quadratic expressions like ``x**2 - 1`` are acceptable
191 as long as the resulting integrand does not depend on the sign of
192 the solutions (see examples).
194 The integral will be returned unchanged if ``x`` is not a variable of
195 integration.
197 ``x`` must be (or contain) only one of of the integration variables. If
198 ``u`` has more than one free symbol then it should be sent as a tuple
199 (``u``, ``uvar``) where ``uvar`` identifies which variable is replacing
200 the integration variable.
201 XXX can it contain another integration variable?
203 Examples
204 ========
206 >>> from sympy.abc import a, x, u
207 >>> from sympy import Integral, cos, sqrt
209 >>> i = Integral(x*cos(x**2 - 1), (x, 0, 1))
211 transform can change the variable of integration
213 >>> i.transform(x, u)
214 Integral(u*cos(u**2 - 1), (u, 0, 1))
216 transform can perform u-substitution as long as a unique
217 integrand is obtained:
219 >>> i.transform(x**2 - 1, u)
220 Integral(cos(u)/2, (u, -1, 0))
222 This attempt fails because x = +/-sqrt(u + 1) and the
223 sign does not cancel out of the integrand:
225 >>> Integral(cos(x**2 - 1), (x, 0, 1)).transform(x**2 - 1, u)
226 Traceback (most recent call last):
227 ...
228 ValueError:
229 The mapping between F(x) and f(u) did not give a unique integrand.
231 transform can do a substitution. Here, the previous
232 result is transformed back into the original expression
233 using "u-substitution":
235 >>> ui = _
236 >>> _.transform(sqrt(u + 1), x) == i
237 True
239 We can accomplish the same with a regular substitution:
241 >>> ui.transform(u, x**2 - 1) == i
242 True
244 If the `x` does not contain a symbol of integration then
245 the integral will be returned unchanged. Integral `i` does
246 not have an integration variable `a` so no change is made:
248 >>> i.transform(a, x) == i
249 True
251 When `u` has more than one free symbol the symbol that is
252 replacing `x` must be identified by passing `u` as a tuple:
254 >>> Integral(x, (x, 0, 1)).transform(x, (u + a, u))
255 Integral(a + u, (u, -a, 1 - a))
256 >>> Integral(x, (x, 0, 1)).transform(x, (u + a, a))
257 Integral(a + u, (a, -u, 1 - u))
259 See Also
260 ========
262 sympy.concrete.expr_with_limits.ExprWithLimits.variables : Lists the integration variables
263 as_dummy : Replace integration variables with dummy ones
264 """
265 d = Dummy('d')
267 xfree = x.free_symbols.intersection(self.variables)
268 if len(xfree) > 1:
269 raise ValueError(
270 'F(x) can only contain one of: %s' % self.variables)
271 xvar = xfree.pop() if xfree else d
273 if xvar not in self.variables:
274 return self
276 u = sympify(u)
277 if isinstance(u, Expr):
278 ufree = u.free_symbols
279 if len(ufree) == 0:
280 raise ValueError(filldedent('''
281 f(u) cannot be a constant'''))
282 if len(ufree) > 1:
283 raise ValueError(filldedent('''
284 When f(u) has more than one free symbol, the one replacing x
285 must be identified: pass f(u) as (f(u), u)'''))
286 uvar = ufree.pop()
287 else:
288 u, uvar = u
289 if uvar not in u.free_symbols:
290 raise ValueError(filldedent('''
291 Expecting a tuple (expr, symbol) where symbol identified
292 a free symbol in expr, but symbol is not in expr's free
293 symbols.'''))
294 if not isinstance(uvar, Symbol):
295 # This probably never evaluates to True
296 raise ValueError(filldedent('''
297 Expecting a tuple (expr, symbol) but didn't get
298 a symbol; got %s''' % uvar))
300 if x.is_Symbol and u.is_Symbol:
301 return self.xreplace({x: u})
303 if not x.is_Symbol and not u.is_Symbol:
304 raise ValueError('either x or u must be a symbol')
306 if uvar == xvar:
307 return self.transform(x, (u.subs(uvar, d), d)).xreplace({d: uvar})
309 if uvar in self.limits:
310 raise ValueError(filldedent('''
311 u must contain the same variable as in x
312 or a variable that is not already an integration variable'''))
314 from sympy.solvers.solvers import solve
315 if not x.is_Symbol:
316 F = [x.subs(xvar, d)]
317 soln = solve(u - x, xvar, check=False)
318 if not soln:
319 raise ValueError('no solution for solve(F(x) - f(u), x)')
320 f = [fi.subs(uvar, d) for fi in soln]
321 else:
322 f = [u.subs(uvar, d)]
323 from sympy.simplify.simplify import posify
324 pdiff, reps = posify(u - x)
325 puvar = uvar.subs([(v, k) for k, v in reps.items()])
326 soln = [s.subs(reps) for s in solve(pdiff, puvar)]
327 if not soln:
328 raise ValueError('no solution for solve(F(x) - f(u), u)')
329 F = [fi.subs(xvar, d) for fi in soln]
331 newfuncs = {(self.function.subs(xvar, fi)*fi.diff(d)
332 ).subs(d, uvar) for fi in f}
333 if len(newfuncs) > 1:
334 raise ValueError(filldedent('''
335 The mapping between F(x) and f(u) did not give
336 a unique integrand.'''))
337 newfunc = newfuncs.pop()
339 def _calc_limit_1(F, a, b):
340 """
341 replace d with a, using subs if possible, otherwise limit
342 where sign of b is considered
343 """
344 wok = F.subs(d, a)
345 if wok is S.NaN or wok.is_finite is False and a.is_finite:
346 return limit(sign(b)*F, d, a)
347 return wok
349 def _calc_limit(a, b):
350 """
351 replace d with a, using subs if possible, otherwise limit
352 where sign of b is considered
353 """
354 avals = list({_calc_limit_1(Fi, a, b) for Fi in F})
355 if len(avals) > 1:
356 raise ValueError(filldedent('''
357 The mapping between F(x) and f(u) did not
358 give a unique limit.'''))
359 return avals[0]
361 newlimits = []
362 for xab in self.limits:
363 sym = xab[0]
364 if sym == xvar:
365 if len(xab) == 3:
366 a, b = xab[1:]
367 a, b = _calc_limit(a, b), _calc_limit(b, a)
368 if fuzzy_bool(a - b > 0):
369 a, b = b, a
370 newfunc = -newfunc
371 newlimits.append((uvar, a, b))
372 elif len(xab) == 2:
373 a = _calc_limit(xab[1], 1)
374 newlimits.append((uvar, a))
375 else:
376 newlimits.append(uvar)
377 else:
378 newlimits.append(xab)
380 return self.func(newfunc, *newlimits)
382 def doit(self, **hints):
383 """
384 Perform the integration using any hints given.
386 Examples
387 ========
389 >>> from sympy import Piecewise, S
390 >>> from sympy.abc import x, t
391 >>> p = x**2 + Piecewise((0, x/t < 0), (1, True))
392 >>> p.integrate((t, S(4)/5, 1), (x, -1, 1))
393 1/3
395 See Also
396 ========
398 sympy.integrals.trigonometry.trigintegrate
399 sympy.integrals.heurisch.heurisch
400 sympy.integrals.rationaltools.ratint
401 as_sum : Approximate the integral using a sum
402 """
403 if not hints.get('integrals', True):
404 return self
406 deep = hints.get('deep', True)
407 meijerg = hints.get('meijerg', None)
408 conds = hints.get('conds', 'piecewise')
409 risch = hints.get('risch', None)
410 heurisch = hints.get('heurisch', None)
411 manual = hints.get('manual', None)
412 if len(list(filter(None, (manual, meijerg, risch, heurisch)))) > 1:
413 raise ValueError("At most one of manual, meijerg, risch, heurisch can be True")
414 elif manual:
415 meijerg = risch = heurisch = False
416 elif meijerg:
417 manual = risch = heurisch = False
418 elif risch:
419 manual = meijerg = heurisch = False
420 elif heurisch:
421 manual = meijerg = risch = False
422 eval_kwargs = {"meijerg": meijerg, "risch": risch, "manual": manual, "heurisch": heurisch,
423 "conds": conds}
425 if conds not in ('separate', 'piecewise', 'none'):
426 raise ValueError('conds must be one of "separate", "piecewise", '
427 '"none", got: %s' % conds)
429 if risch and any(len(xab) > 1 for xab in self.limits):
430 raise ValueError('risch=True is only allowed for indefinite integrals.')
432 # check for the trivial zero
433 if self.is_zero:
434 return S.Zero
436 # hacks to handle integrals of
437 # nested summations
438 from sympy.concrete.summations import Sum
439 if isinstance(self.function, Sum):
440 if any(v in self.function.limits[0] for v in self.variables):
441 raise ValueError('Limit of the sum cannot be an integration variable.')
442 if any(l.is_infinite for l in self.function.limits[0][1:]):
443 return self
444 _i = self
445 _sum = self.function
446 return _sum.func(_i.func(_sum.function, *_i.limits).doit(), *_sum.limits).doit()
448 # now compute and check the function
449 function = self.function
450 if deep:
451 function = function.doit(**hints)
452 if function.is_zero:
453 return S.Zero
455 # hacks to handle special cases
456 if isinstance(function, MatrixBase):
457 return function.applyfunc(
458 lambda f: self.func(f, *self.limits).doit(**hints))
460 if isinstance(function, FormalPowerSeries):
461 if len(self.limits) > 1:
462 raise NotImplementedError
463 xab = self.limits[0]
464 if len(xab) > 1:
465 return function.integrate(xab, **eval_kwargs)
466 else:
467 return function.integrate(xab[0], **eval_kwargs)
469 # There is no trivial answer and special handling
470 # is done so continue
472 # first make sure any definite limits have integration
473 # variables with matching assumptions
474 reps = {}
475 for xab in self.limits:
476 if len(xab) != 3:
477 # it makes sense to just make
478 # all x real but in practice with the
479 # current state of integration...this
480 # doesn't work out well
481 # x = xab[0]
482 # if x not in reps and not x.is_real:
483 # reps[x] = Dummy(real=True)
484 continue
485 x, a, b = xab
486 l = (a, b)
487 if all(i.is_nonnegative for i in l) and not x.is_nonnegative:
488 d = Dummy(positive=True)
489 elif all(i.is_nonpositive for i in l) and not x.is_nonpositive:
490 d = Dummy(negative=True)
491 elif all(i.is_real for i in l) and not x.is_real:
492 d = Dummy(real=True)
493 else:
494 d = None
495 if d:
496 reps[x] = d
497 if reps:
498 undo = {v: k for k, v in reps.items()}
499 did = self.xreplace(reps).doit(**hints)
500 if isinstance(did, tuple): # when separate=True
501 did = tuple([i.xreplace(undo) for i in did])
502 else:
503 did = did.xreplace(undo)
504 return did
506 # continue with existing assumptions
507 undone_limits = []
508 # ulj = free symbols of any undone limits' upper and lower limits
509 ulj = set()
510 for xab in self.limits:
511 # compute uli, the free symbols in the
512 # Upper and Lower limits of limit I
513 if len(xab) == 1:
514 uli = set(xab[:1])
515 elif len(xab) == 2:
516 uli = xab[1].free_symbols
517 elif len(xab) == 3:
518 uli = xab[1].free_symbols.union(xab[2].free_symbols)
519 # this integral can be done as long as there is no blocking
520 # limit that has been undone. An undone limit is blocking if
521 # it contains an integration variable that is in this limit's
522 # upper or lower free symbols or vice versa
523 if xab[0] in ulj or any(v[0] in uli for v in undone_limits):
524 undone_limits.append(xab)
525 ulj.update(uli)
526 function = self.func(*([function] + [xab]))
527 factored_function = function.factor()
528 if not isinstance(factored_function, Integral):
529 function = factored_function
530 continue
532 if function.has(Abs, sign) and (
533 (len(xab) < 3 and all(x.is_extended_real for x in xab)) or
534 (len(xab) == 3 and all(x.is_extended_real and not x.is_infinite for
535 x in xab[1:]))):
536 # some improper integrals are better off with Abs
537 xr = Dummy("xr", real=True)
538 function = (function.xreplace({xab[0]: xr})
539 .rewrite(Piecewise).xreplace({xr: xab[0]}))
540 elif function.has(Min, Max):
541 function = function.rewrite(Piecewise)
542 if (function.has(Piecewise) and
543 not isinstance(function, Piecewise)):
544 function = piecewise_fold(function)
545 if isinstance(function, Piecewise):
546 if len(xab) == 1:
547 antideriv = function._eval_integral(xab[0],
548 **eval_kwargs)
549 else:
550 antideriv = self._eval_integral(
551 function, xab[0], **eval_kwargs)
552 else:
553 # There are a number of tradeoffs in using the
554 # Meijer G method. It can sometimes be a lot faster
555 # than other methods, and sometimes slower. And
556 # there are certain types of integrals for which it
557 # is more likely to work than others. These
558 # heuristics are incorporated in deciding what
559 # integration methods to try, in what order. See the
560 # integrate() docstring for details.
561 def try_meijerg(function, xab):
562 ret = None
563 if len(xab) == 3 and meijerg is not False:
564 x, a, b = xab
565 try:
566 res = meijerint_definite(function, x, a, b)
567 except NotImplementedError:
568 _debug('NotImplementedError '
569 'from meijerint_definite')
570 res = None
571 if res is not None:
572 f, cond = res
573 if conds == 'piecewise':
574 u = self.func(function, (x, a, b))
575 # if Piecewise modifies cond too
576 # much it may not be recognized by
577 # _condsimp pattern matching so just
578 # turn off all evaluation
579 return Piecewise((f, cond), (u, True),
580 evaluate=False)
581 elif conds == 'separate':
582 if len(self.limits) != 1:
583 raise ValueError(filldedent('''
584 conds=separate not supported in
585 multiple integrals'''))
586 ret = f, cond
587 else:
588 ret = f
589 return ret
591 meijerg1 = meijerg
592 if (meijerg is not False and
593 len(xab) == 3 and xab[1].is_extended_real and xab[2].is_extended_real
594 and not function.is_Poly and
595 (xab[1].has(oo, -oo) or xab[2].has(oo, -oo))):
596 ret = try_meijerg(function, xab)
597 if ret is not None:
598 function = ret
599 continue
600 meijerg1 = False
601 # If the special meijerg code did not succeed in
602 # finding a definite integral, then the code using
603 # meijerint_indefinite will not either (it might
604 # find an antiderivative, but the answer is likely
605 # to be nonsensical). Thus if we are requested to
606 # only use Meijer G-function methods, we give up at
607 # this stage. Otherwise we just disable G-function
608 # methods.
609 if meijerg1 is False and meijerg is True:
610 antideriv = None
611 else:
612 antideriv = self._eval_integral(
613 function, xab[0], **eval_kwargs)
614 if antideriv is None and meijerg is True:
615 ret = try_meijerg(function, xab)
616 if ret is not None:
617 function = ret
618 continue
620 final = hints.get('final', True)
621 # dotit may be iterated but floor terms making atan and acot
622 # continuous should only be added in the final round
623 if (final and not isinstance(antideriv, Integral) and
624 antideriv is not None):
625 for atan_term in antideriv.atoms(atan):
626 atan_arg = atan_term.args[0]
627 # Checking `atan_arg` to be linear combination of `tan` or `cot`
628 for tan_part in atan_arg.atoms(tan):
629 x1 = Dummy('x1')
630 tan_exp1 = atan_arg.subs(tan_part, x1)
631 # The coefficient of `tan` should be constant
632 coeff = tan_exp1.diff(x1)
633 if x1 not in coeff.free_symbols:
634 a = tan_part.args[0]
635 antideriv = antideriv.subs(atan_term, Add(atan_term,
636 sign(coeff)*pi*floor((a-pi/2)/pi)))
637 for cot_part in atan_arg.atoms(cot):
638 x1 = Dummy('x1')
639 cot_exp1 = atan_arg.subs(cot_part, x1)
640 # The coefficient of `cot` should be constant
641 coeff = cot_exp1.diff(x1)
642 if x1 not in coeff.free_symbols:
643 a = cot_part.args[0]
644 antideriv = antideriv.subs(atan_term, Add(atan_term,
645 sign(coeff)*pi*floor((a)/pi)))
647 if antideriv is None:
648 undone_limits.append(xab)
649 function = self.func(*([function] + [xab])).factor()
650 factored_function = function.factor()
651 if not isinstance(factored_function, Integral):
652 function = factored_function
653 continue
654 else:
655 if len(xab) == 1:
656 function = antideriv
657 else:
658 if len(xab) == 3:
659 x, a, b = xab
660 elif len(xab) == 2:
661 x, b = xab
662 a = None
663 else:
664 raise NotImplementedError
666 if deep:
667 if isinstance(a, Basic):
668 a = a.doit(**hints)
669 if isinstance(b, Basic):
670 b = b.doit(**hints)
672 if antideriv.is_Poly:
673 gens = list(antideriv.gens)
674 gens.remove(x)
676 antideriv = antideriv.as_expr()
678 function = antideriv._eval_interval(x, a, b)
679 function = Poly(function, *gens)
680 else:
681 def is_indef_int(g, x):
682 return (isinstance(g, Integral) and
683 any(i == (x,) for i in g.limits))
685 def eval_factored(f, x, a, b):
686 # _eval_interval for integrals with
687 # (constant) factors
688 # a single indefinite integral is assumed
689 args = []
690 for g in Mul.make_args(f):
691 if is_indef_int(g, x):
692 args.append(g._eval_interval(x, a, b))
693 else:
694 args.append(g)
695 return Mul(*args)
697 integrals, others, piecewises = [], [], []
698 for f in Add.make_args(antideriv):
699 if any(is_indef_int(g, x)
700 for g in Mul.make_args(f)):
701 integrals.append(f)
702 elif any(isinstance(g, Piecewise)
703 for g in Mul.make_args(f)):
704 piecewises.append(piecewise_fold(f))
705 else:
706 others.append(f)
707 uneval = Add(*[eval_factored(f, x, a, b)
708 for f in integrals])
709 try:
710 evalued = Add(*others)._eval_interval(x, a, b)
711 evalued_pw = piecewise_fold(Add(*piecewises))._eval_interval(x, a, b)
712 function = uneval + evalued + evalued_pw
713 except NotImplementedError:
714 # This can happen if _eval_interval depends in a
715 # complicated way on limits that cannot be computed
716 undone_limits.append(xab)
717 function = self.func(*([function] + [xab]))
718 factored_function = function.factor()
719 if not isinstance(factored_function, Integral):
720 function = factored_function
721 return function
723 def _eval_derivative(self, sym):
724 """Evaluate the derivative of the current Integral object by
725 differentiating under the integral sign [1], using the Fundamental
726 Theorem of Calculus [2] when possible.
728 Explanation
729 ===========
731 Whenever an Integral is encountered that is equivalent to zero or
732 has an integrand that is independent of the variable of integration
733 those integrals are performed. All others are returned as Integral
734 instances which can be resolved with doit() (provided they are integrable).
736 References
737 ==========
739 .. [1] https://en.wikipedia.org/wiki/Differentiation_under_the_integral_sign
740 .. [2] https://en.wikipedia.org/wiki/Fundamental_theorem_of_calculus
742 Examples
743 ========
745 >>> from sympy import Integral
746 >>> from sympy.abc import x, y
747 >>> i = Integral(x + y, y, (y, 1, x))
748 >>> i.diff(x)
749 Integral(x + y, (y, x)) + Integral(1, y, (y, 1, x))
750 >>> i.doit().diff(x) == i.diff(x).doit()
751 True
752 >>> i.diff(y)
753 0
755 The previous must be true since there is no y in the evaluated integral:
757 >>> i.free_symbols
758 {x}
759 >>> i.doit()
760 2*x**3/3 - x/2 - 1/6
762 """
764 # differentiate under the integral sign; we do not
765 # check for regularity conditions (TODO), see issue 4215
767 # get limits and the function
768 f, limits = self.function, list(self.limits)
770 # the order matters if variables of integration appear in the limits
771 # so work our way in from the outside to the inside.
772 limit = limits.pop(-1)
773 if len(limit) == 3:
774 x, a, b = limit
775 elif len(limit) == 2:
776 x, b = limit
777 a = None
778 else:
779 a = b = None
780 x = limit[0]
782 if limits: # f is the argument to an integral
783 f = self.func(f, *tuple(limits))
785 # assemble the pieces
786 def _do(f, ab):
787 dab_dsym = diff(ab, sym)
788 if not dab_dsym:
789 return S.Zero
790 if isinstance(f, Integral):
791 limits = [(x, x) if (len(l) == 1 and l[0] == x) else l
792 for l in f.limits]
793 f = self.func(f.function, *limits)
794 return f.subs(x, ab)*dab_dsym
796 rv = S.Zero
797 if b is not None:
798 rv += _do(f, b)
799 if a is not None:
800 rv -= _do(f, a)
801 if len(limit) == 1 and sym == x:
802 # the dummy variable *is* also the real-world variable
803 arg = f
804 rv += arg
805 else:
806 # the dummy variable might match sym but it's
807 # only a dummy and the actual variable is determined
808 # by the limits, so mask off the variable of integration
809 # while differentiating
810 u = Dummy('u')
811 arg = f.subs(x, u).diff(sym).subs(u, x)
812 if arg:
813 rv += self.func(arg, (x, a, b))
814 return rv
816 def _eval_integral(self, f, x, meijerg=None, risch=None, manual=None,
817 heurisch=None, conds='piecewise',final=None):
818 """
819 Calculate the anti-derivative to the function f(x).
821 Explanation
822 ===========
824 The following algorithms are applied (roughly in this order):
826 1. Simple heuristics (based on pattern matching and integral table):
828 - most frequently used functions (e.g. polynomials, products of
829 trig functions)
831 2. Integration of rational functions:
833 - A complete algorithm for integrating rational functions is
834 implemented (the Lazard-Rioboo-Trager algorithm). The algorithm
835 also uses the partial fraction decomposition algorithm
836 implemented in apart() as a preprocessor to make this process
837 faster. Note that the integral of a rational function is always
838 elementary, but in general, it may include a RootSum.
840 3. Full Risch algorithm:
842 - The Risch algorithm is a complete decision
843 procedure for integrating elementary functions, which means that
844 given any elementary function, it will either compute an
845 elementary antiderivative, or else prove that none exists.
846 Currently, part of transcendental case is implemented, meaning
847 elementary integrals containing exponentials, logarithms, and
848 (soon!) trigonometric functions can be computed. The algebraic
849 case, e.g., functions containing roots, is much more difficult
850 and is not implemented yet.
852 - If the routine fails (because the integrand is not elementary, or
853 because a case is not implemented yet), it continues on to the
854 next algorithms below. If the routine proves that the integrals
855 is nonelementary, it still moves on to the algorithms below,
856 because we might be able to find a closed-form solution in terms
857 of special functions. If risch=True, however, it will stop here.
859 4. The Meijer G-Function algorithm:
861 - This algorithm works by first rewriting the integrand in terms of
862 very general Meijer G-Function (meijerg in SymPy), integrating
863 it, and then rewriting the result back, if possible. This
864 algorithm is particularly powerful for definite integrals (which
865 is actually part of a different method of Integral), since it can
866 compute closed-form solutions of definite integrals even when no
867 closed-form indefinite integral exists. But it also is capable
868 of computing many indefinite integrals as well.
870 - Another advantage of this method is that it can use some results
871 about the Meijer G-Function to give a result in terms of a
872 Piecewise expression, which allows to express conditionally
873 convergent integrals.
875 - Setting meijerg=True will cause integrate() to use only this
876 method.
878 5. The "manual integration" algorithm:
880 - This algorithm tries to mimic how a person would find an
881 antiderivative by hand, for example by looking for a
882 substitution or applying integration by parts. This algorithm
883 does not handle as many integrands but can return results in a
884 more familiar form.
886 - Sometimes this algorithm can evaluate parts of an integral; in
887 this case integrate() will try to evaluate the rest of the
888 integrand using the other methods here.
890 - Setting manual=True will cause integrate() to use only this
891 method.
893 6. The Heuristic Risch algorithm:
895 - This is a heuristic version of the Risch algorithm, meaning that
896 it is not deterministic. This is tried as a last resort because
897 it can be very slow. It is still used because not enough of the
898 full Risch algorithm is implemented, so that there are still some
899 integrals that can only be computed using this method. The goal
900 is to implement enough of the Risch and Meijer G-function methods
901 so that this can be deleted.
903 Setting heurisch=True will cause integrate() to use only this
904 method. Set heurisch=False to not use it.
906 """
908 from sympy.integrals.risch import risch_integrate, NonElementaryIntegral
909 from sympy.integrals.manualintegrate import manualintegrate
911 if risch:
912 try:
913 return risch_integrate(f, x, conds=conds)
914 except NotImplementedError:
915 return None
917 if manual:
918 try:
919 result = manualintegrate(f, x)
920 if result is not None and result.func != Integral:
921 return result
922 except (ValueError, PolynomialError):
923 pass
925 eval_kwargs = {"meijerg": meijerg, "risch": risch, "manual": manual,
926 "heurisch": heurisch, "conds": conds}
928 # if it is a poly(x) then let the polynomial integrate itself (fast)
929 #
930 # It is important to make this check first, otherwise the other code
931 # will return a SymPy expression instead of a Polynomial.
932 #
933 # see Polynomial for details.
934 if isinstance(f, Poly) and not (manual or meijerg or risch):
935 # Note: this is deprecated, but the deprecation warning is already
936 # issued in the Integral constructor.
937 return f.integrate(x)
939 # Piecewise antiderivatives need to call special integrate.
940 if isinstance(f, Piecewise):
941 return f.piecewise_integrate(x, **eval_kwargs)
943 # let's cut it short if `f` does not depend on `x`; if
944 # x is only a dummy, that will be handled below
945 if not f.has(x):
946 return f*x
948 # try to convert to poly(x) and then integrate if successful (fast)
949 poly = f.as_poly(x)
950 if poly is not None and not (manual or meijerg or risch):
951 return poly.integrate().as_expr()
953 if risch is not False:
954 try:
955 result, i = risch_integrate(f, x, separate_integral=True,
956 conds=conds)
957 except NotImplementedError:
958 pass
959 else:
960 if i:
961 # There was a nonelementary integral. Try integrating it.
963 # if no part of the NonElementaryIntegral is integrated by
964 # the Risch algorithm, then use the original function to
965 # integrate, instead of re-written one
966 if result == 0:
967 return NonElementaryIntegral(f, x).doit(risch=False)
968 else:
969 return result + i.doit(risch=False)
970 else:
971 return result
973 # since Integral(f=g1+g2+...) == Integral(g1) + Integral(g2) + ...
974 # we are going to handle Add terms separately,
975 # if `f` is not Add -- we only have one term
977 # Note that in general, this is a bad idea, because Integral(g1) +
978 # Integral(g2) might not be computable, even if Integral(g1 + g2) is.
979 # For example, Integral(x**x + x**x*log(x)). But many heuristics only
980 # work term-wise. So we compute this step last, after trying
981 # risch_integrate. We also try risch_integrate again in this loop,
982 # because maybe the integral is a sum of an elementary part and a
983 # nonelementary part (like erf(x) + exp(x)). risch_integrate() is
984 # quite fast, so this is acceptable.
985 from sympy.simplify.fu import sincos_to_sum
986 parts = []
987 args = Add.make_args(f)
988 for g in args:
989 coeff, g = g.as_independent(x)
991 # g(x) = const
992 if g is S.One and not meijerg:
993 parts.append(coeff*x)
994 continue
996 # g(x) = expr + O(x**n)
997 order_term = g.getO()
999 if order_term is not None:
1000 h = self._eval_integral(g.removeO(), x, **eval_kwargs)
1002 if h is not None:
1003 h_order_expr = self._eval_integral(order_term.expr, x, **eval_kwargs)
1005 if h_order_expr is not None:
1006 h_order_term = order_term.func(
1007 h_order_expr, *order_term.variables)
1008 parts.append(coeff*(h + h_order_term))
1009 continue
1011 # NOTE: if there is O(x**n) and we fail to integrate then
1012 # there is no point in trying other methods because they
1013 # will fail, too.
1014 return None
1016 # c
1017 # g(x) = (a*x+b)
1018 if g.is_Pow and not g.exp.has(x) and not meijerg:
1019 a = Wild('a', exclude=[x])
1020 b = Wild('b', exclude=[x])
1022 M = g.base.match(a*x + b)
1024 if M is not None:
1025 if g.exp == -1:
1026 h = log(g.base)
1027 elif conds != 'piecewise':
1028 h = g.base**(g.exp + 1) / (g.exp + 1)
1029 else:
1030 h1 = log(g.base)
1031 h2 = g.base**(g.exp + 1) / (g.exp + 1)
1032 h = Piecewise((h2, Ne(g.exp, -1)), (h1, True))
1034 parts.append(coeff * h / M[a])
1035 continue
1037 # poly(x)
1038 # g(x) = -------
1039 # poly(x)
1040 if g.is_rational_function(x) and not (manual or meijerg or risch):
1041 parts.append(coeff * ratint(g, x))
1042 continue
1044 if not (manual or meijerg or risch):
1045 # g(x) = Mul(trig)
1046 h = trigintegrate(g, x, conds=conds)
1047 if h is not None:
1048 parts.append(coeff * h)
1049 continue
1051 # g(x) has at least a DiracDelta term
1052 h = deltaintegrate(g, x)
1053 if h is not None:
1054 parts.append(coeff * h)
1055 continue
1057 from .singularityfunctions import singularityintegrate
1058 # g(x) has at least a Singularity Function term
1059 h = singularityintegrate(g, x)
1060 if h is not None:
1061 parts.append(coeff * h)
1062 continue
1064 # Try risch again.
1065 if risch is not False:
1066 try:
1067 h, i = risch_integrate(g, x,
1068 separate_integral=True, conds=conds)
1069 except NotImplementedError:
1070 h = None
1071 else:
1072 if i:
1073 h = h + i.doit(risch=False)
1075 parts.append(coeff*h)
1076 continue
1078 # fall back to heurisch
1079 if heurisch is not False:
1080 from sympy.integrals.heurisch import (heurisch as heurisch_,
1081 heurisch_wrapper)
1082 try:
1083 if conds == 'piecewise':
1084 h = heurisch_wrapper(g, x, hints=[])
1085 else:
1086 h = heurisch_(g, x, hints=[])
1087 except PolynomialError:
1088 # XXX: this exception means there is a bug in the
1089 # implementation of heuristic Risch integration
1090 # algorithm.
1091 h = None
1092 else:
1093 h = None
1095 if meijerg is not False and h is None:
1096 # rewrite using G functions
1097 try:
1098 h = meijerint_indefinite(g, x)
1099 except NotImplementedError:
1100 _debug('NotImplementedError from meijerint_definite')
1101 if h is not None:
1102 parts.append(coeff * h)
1103 continue
1105 if h is None and manual is not False:
1106 try:
1107 result = manualintegrate(g, x)
1108 if result is not None and not isinstance(result, Integral):
1109 if result.has(Integral) and not manual:
1110 # Try to have other algorithms do the integrals
1111 # manualintegrate can't handle,
1112 # unless we were asked to use manual only.
1113 # Keep the rest of eval_kwargs in case another
1114 # method was set to False already
1115 new_eval_kwargs = eval_kwargs
1116 new_eval_kwargs["manual"] = False
1117 new_eval_kwargs["final"] = False
1118 result = result.func(*[
1119 arg.doit(**new_eval_kwargs) if
1120 arg.has(Integral) else arg
1121 for arg in result.args
1122 ]).expand(multinomial=False,
1123 log=False,
1124 power_exp=False,
1125 power_base=False)
1126 if not result.has(Integral):
1127 parts.append(coeff * result)
1128 continue
1129 except (ValueError, PolynomialError):
1130 # can't handle some SymPy expressions
1131 pass
1133 # if we failed maybe it was because we had
1134 # a product that could have been expanded,
1135 # so let's try an expansion of the whole
1136 # thing before giving up; we don't try this
1137 # at the outset because there are things
1138 # that cannot be solved unless they are
1139 # NOT expanded e.g., x**x*(1+log(x)). There
1140 # should probably be a checker somewhere in this
1141 # routine to look for such cases and try to do
1142 # collection on the expressions if they are already
1143 # in an expanded form
1144 if not h and len(args) == 1:
1145 f = sincos_to_sum(f).expand(mul=True, deep=False)
1146 if f.is_Add:
1147 # Note: risch will be identical on the expanded
1148 # expression, but maybe it will be able to pick out parts,
1149 # like x*(exp(x) + erf(x)).
1150 return self._eval_integral(f, x, **eval_kwargs)
1152 if h is not None:
1153 parts.append(coeff * h)
1154 else:
1155 return None
1157 return Add(*parts)
1159 def _eval_lseries(self, x, logx=None, cdir=0):
1160 expr = self.as_dummy()
1161 symb = x
1162 for l in expr.limits:
1163 if x in l[1:]:
1164 symb = l[0]
1165 break
1166 for term in expr.function.lseries(symb, logx):
1167 yield integrate(term, *expr.limits)
1169 def _eval_nseries(self, x, n, logx=None, cdir=0):
1170 expr = self.as_dummy()
1171 symb = x
1172 for l in expr.limits:
1173 if x in l[1:]:
1174 symb = l[0]
1175 break
1176 terms, order = expr.function.nseries(
1177 x=symb, n=n, logx=logx).as_coeff_add(Order)
1178 order = [o.subs(symb, x) for o in order]
1179 return integrate(terms, *expr.limits) + Add(*order)*x
1181 def _eval_as_leading_term(self, x, logx=None, cdir=0):
1182 series_gen = self.args[0].lseries(x)
1183 for leading_term in series_gen:
1184 if leading_term != 0:
1185 break
1186 return integrate(leading_term, *self.args[1:])
1188 def _eval_simplify(self, **kwargs):
1189 expr = factor_terms(self)
1190 if isinstance(expr, Integral):
1191 from sympy.simplify.simplify import simplify
1192 return expr.func(*[simplify(i, **kwargs) for i in expr.args])
1193 return expr.simplify(**kwargs)
1195 def as_sum(self, n=None, method="midpoint", evaluate=True):
1196 """
1197 Approximates a definite integral by a sum.
1199 Parameters
1200 ==========
1202 n :
1203 The number of subintervals to use, optional.
1204 method :
1205 One of: 'left', 'right', 'midpoint', 'trapezoid'.
1206 evaluate : bool
1207 If False, returns an unevaluated Sum expression. The default
1208 is True, evaluate the sum.
1210 Notes
1211 =====
1213 These methods of approximate integration are described in [1].
1215 Examples
1216 ========
1218 >>> from sympy import Integral, sin, sqrt
1219 >>> from sympy.abc import x, n
1220 >>> e = Integral(sin(x), (x, 3, 7))
1221 >>> e
1222 Integral(sin(x), (x, 3, 7))
1224 For demonstration purposes, this interval will only be split into 2
1225 regions, bounded by [3, 5] and [5, 7].
1227 The left-hand rule uses function evaluations at the left of each
1228 interval:
1230 >>> e.as_sum(2, 'left')
1231 2*sin(5) + 2*sin(3)
1233 The midpoint rule uses evaluations at the center of each interval:
1235 >>> e.as_sum(2, 'midpoint')
1236 2*sin(4) + 2*sin(6)
1238 The right-hand rule uses function evaluations at the right of each
1239 interval:
1241 >>> e.as_sum(2, 'right')
1242 2*sin(5) + 2*sin(7)
1244 The trapezoid rule uses function evaluations on both sides of the
1245 intervals. This is equivalent to taking the average of the left and
1246 right hand rule results:
1248 >>> e.as_sum(2, 'trapezoid')
1249 2*sin(5) + sin(3) + sin(7)
1250 >>> (e.as_sum(2, 'left') + e.as_sum(2, 'right'))/2 == _
1251 True
1253 Here, the discontinuity at x = 0 can be avoided by using the
1254 midpoint or right-hand method:
1256 >>> e = Integral(1/sqrt(x), (x, 0, 1))
1257 >>> e.as_sum(5).n(4)
1258 1.730
1259 >>> e.as_sum(10).n(4)
1260 1.809
1261 >>> e.doit().n(4) # the actual value is 2
1262 2.000
1264 The left- or trapezoid method will encounter the discontinuity and
1265 return infinity:
1267 >>> e.as_sum(5, 'left')
1268 zoo
1270 The number of intervals can be symbolic. If omitted, a dummy symbol
1271 will be used for it.
1273 >>> e = Integral(x**2, (x, 0, 2))
1274 >>> e.as_sum(n, 'right').expand()
1275 8/3 + 4/n + 4/(3*n**2)
1277 This shows that the midpoint rule is more accurate, as its error
1278 term decays as the square of n:
1280 >>> e.as_sum(method='midpoint').expand()
1281 8/3 - 2/(3*_n**2)
1283 A symbolic sum is returned with evaluate=False:
1285 >>> e.as_sum(n, 'midpoint', evaluate=False)
1286 2*Sum((2*_k/n - 1/n)**2, (_k, 1, n))/n
1288 See Also
1289 ========
1291 Integral.doit : Perform the integration using any hints
1293 References
1294 ==========
1296 .. [1] https://en.wikipedia.org/wiki/Riemann_sum#Riemann_summation_methods
1297 """
1299 from sympy.concrete.summations import Sum
1300 limits = self.limits
1301 if len(limits) > 1:
1302 raise NotImplementedError(
1303 "Multidimensional midpoint rule not implemented yet")
1304 else:
1305 limit = limits[0]
1306 if (len(limit) != 3 or limit[1].is_finite is False or
1307 limit[2].is_finite is False):
1308 raise ValueError("Expecting a definite integral over "
1309 "a finite interval.")
1310 if n is None:
1311 n = Dummy('n', integer=True, positive=True)
1312 else:
1313 n = sympify(n)
1314 if (n.is_positive is False or n.is_integer is False or
1315 n.is_finite is False):
1316 raise ValueError("n must be a positive integer, got %s" % n)
1317 x, a, b = limit
1318 dx = (b - a)/n
1319 k = Dummy('k', integer=True, positive=True)
1320 f = self.function
1322 if method == "left":
1323 result = dx*Sum(f.subs(x, a + (k-1)*dx), (k, 1, n))
1324 elif method == "right":
1325 result = dx*Sum(f.subs(x, a + k*dx), (k, 1, n))
1326 elif method == "midpoint":
1327 result = dx*Sum(f.subs(x, a + k*dx - dx/2), (k, 1, n))
1328 elif method == "trapezoid":
1329 result = dx*((f.subs(x, a) + f.subs(x, b))/2 +
1330 Sum(f.subs(x, a + k*dx), (k, 1, n - 1)))
1331 else:
1332 raise ValueError("Unknown method %s" % method)
1333 return result.doit() if evaluate else result
1335 def principal_value(self, **kwargs):
1336 """
1337 Compute the Cauchy Principal Value of the definite integral of a real function in the given interval
1338 on the real axis.
1340 Explanation
1341 ===========
1343 In mathematics, the Cauchy principal value, is a method for assigning values to certain improper
1344 integrals which would otherwise be undefined.
1346 Examples
1347 ========
1349 >>> from sympy import Integral, oo
1350 >>> from sympy.abc import x
1351 >>> Integral(x+1, (x, -oo, oo)).principal_value()
1352 oo
1353 >>> f = 1 / (x**3)
1354 >>> Integral(f, (x, -oo, oo)).principal_value()
1355 0
1356 >>> Integral(f, (x, -10, 10)).principal_value()
1357 0
1358 >>> Integral(f, (x, -10, oo)).principal_value() + Integral(f, (x, -oo, 10)).principal_value()
1359 0
1361 References
1362 ==========
1364 .. [1] https://en.wikipedia.org/wiki/Cauchy_principal_value
1365 .. [2] https://mathworld.wolfram.com/CauchyPrincipalValue.html
1366 """
1367 if len(self.limits) != 1 or len(list(self.limits[0])) != 3:
1368 raise ValueError("You need to insert a variable, lower_limit, and upper_limit correctly to calculate "
1369 "cauchy's principal value")
1370 x, a, b = self.limits[0]
1371 if not (a.is_comparable and b.is_comparable and a <= b):
1372 raise ValueError("The lower_limit must be smaller than or equal to the upper_limit to calculate "
1373 "cauchy's principal value. Also, a and b need to be comparable.")
1374 if a == b:
1375 return S.Zero
1377 from sympy.calculus.singularities import singularities
1379 r = Dummy('r')
1380 f = self.function
1381 singularities_list = [s for s in singularities(f, x) if s.is_comparable and a <= s <= b]
1382 for i in singularities_list:
1383 if i in (a, b):
1384 raise ValueError(
1385 'The principal value is not defined in the given interval due to singularity at %d.' % (i))
1386 F = integrate(f, x, **kwargs)
1387 if F.has(Integral):
1388 return self
1389 if a is -oo and b is oo:
1390 I = limit(F - F.subs(x, -x), x, oo)
1391 else:
1392 I = limit(F, x, b, '-') - limit(F, x, a, '+')
1393 for s in singularities_list:
1394 I += limit(((F.subs(x, s - r)) - F.subs(x, s + r)), r, 0, '+')
1395 return I
1399def integrate(*args, meijerg=None, conds='piecewise', risch=None, heurisch=None, manual=None, **kwargs):
1400 """integrate(f, var, ...)
1402 .. deprecated:: 1.6
1404 Using ``integrate()`` with :class:`~.Poly` is deprecated. Use
1405 :meth:`.Poly.integrate` instead. See :ref:`deprecated-integrate-poly`.
1407 Explanation
1408 ===========
1410 Compute definite or indefinite integral of one or more variables
1411 using Risch-Norman algorithm and table lookup. This procedure is
1412 able to handle elementary algebraic and transcendental functions
1413 and also a huge class of special functions, including Airy,
1414 Bessel, Whittaker and Lambert.
1416 var can be:
1418 - a symbol -- indefinite integration
1419 - a tuple (symbol, a) -- indefinite integration with result
1420 given with ``a`` replacing ``symbol``
1421 - a tuple (symbol, a, b) -- definite integration
1423 Several variables can be specified, in which case the result is
1424 multiple integration. (If var is omitted and the integrand is
1425 univariate, the indefinite integral in that variable will be performed.)
1427 Indefinite integrals are returned without terms that are independent
1428 of the integration variables. (see examples)
1430 Definite improper integrals often entail delicate convergence
1431 conditions. Pass conds='piecewise', 'separate' or 'none' to have
1432 these returned, respectively, as a Piecewise function, as a separate
1433 result (i.e. result will be a tuple), or not at all (default is
1434 'piecewise').
1436 **Strategy**
1438 SymPy uses various approaches to definite integration. One method is to
1439 find an antiderivative for the integrand, and then use the fundamental
1440 theorem of calculus. Various functions are implemented to integrate
1441 polynomial, rational and trigonometric functions, and integrands
1442 containing DiracDelta terms.
1444 SymPy also implements the part of the Risch algorithm, which is a decision
1445 procedure for integrating elementary functions, i.e., the algorithm can
1446 either find an elementary antiderivative, or prove that one does not
1447 exist. There is also a (very successful, albeit somewhat slow) general
1448 implementation of the heuristic Risch algorithm. This algorithm will
1449 eventually be phased out as more of the full Risch algorithm is
1450 implemented. See the docstring of Integral._eval_integral() for more
1451 details on computing the antiderivative using algebraic methods.
1453 The option risch=True can be used to use only the (full) Risch algorithm.
1454 This is useful if you want to know if an elementary function has an
1455 elementary antiderivative. If the indefinite Integral returned by this
1456 function is an instance of NonElementaryIntegral, that means that the
1457 Risch algorithm has proven that integral to be non-elementary. Note that
1458 by default, additional methods (such as the Meijer G method outlined
1459 below) are tried on these integrals, as they may be expressible in terms
1460 of special functions, so if you only care about elementary answers, use
1461 risch=True. Also note that an unevaluated Integral returned by this
1462 function is not necessarily a NonElementaryIntegral, even with risch=True,
1463 as it may just be an indication that the particular part of the Risch
1464 algorithm needed to integrate that function is not yet implemented.
1466 Another family of strategies comes from re-writing the integrand in
1467 terms of so-called Meijer G-functions. Indefinite integrals of a
1468 single G-function can always be computed, and the definite integral
1469 of a product of two G-functions can be computed from zero to
1470 infinity. Various strategies are implemented to rewrite integrands
1471 as G-functions, and use this information to compute integrals (see
1472 the ``meijerint`` module).
1474 The option manual=True can be used to use only an algorithm that tries
1475 to mimic integration by hand. This algorithm does not handle as many
1476 integrands as the other algorithms implemented but may return results in
1477 a more familiar form. The ``manualintegrate`` module has functions that
1478 return the steps used (see the module docstring for more information).
1480 In general, the algebraic methods work best for computing
1481 antiderivatives of (possibly complicated) combinations of elementary
1482 functions. The G-function methods work best for computing definite
1483 integrals from zero to infinity of moderately complicated
1484 combinations of special functions, or indefinite integrals of very
1485 simple combinations of special functions.
1487 The strategy employed by the integration code is as follows:
1489 - If computing a definite integral, and both limits are real,
1490 and at least one limit is +- oo, try the G-function method of
1491 definite integration first.
1493 - Try to find an antiderivative, using all available methods, ordered
1494 by performance (that is try fastest method first, slowest last; in
1495 particular polynomial integration is tried first, Meijer
1496 G-functions second to last, and heuristic Risch last).
1498 - If still not successful, try G-functions irrespective of the
1499 limits.
1501 The option meijerg=True, False, None can be used to, respectively:
1502 always use G-function methods and no others, never use G-function
1503 methods, or use all available methods (in order as described above).
1504 It defaults to None.
1506 Examples
1507 ========
1509 >>> from sympy import integrate, log, exp, oo
1510 >>> from sympy.abc import a, x, y
1512 >>> integrate(x*y, x)
1513 x**2*y/2
1515 >>> integrate(log(x), x)
1516 x*log(x) - x
1518 >>> integrate(log(x), (x, 1, a))
1519 a*log(a) - a + 1
1521 >>> integrate(x)
1522 x**2/2
1524 Terms that are independent of x are dropped by indefinite integration:
1526 >>> from sympy import sqrt
1527 >>> integrate(sqrt(1 + x), (x, 0, x))
1528 2*(x + 1)**(3/2)/3 - 2/3
1529 >>> integrate(sqrt(1 + x), x)
1530 2*(x + 1)**(3/2)/3
1532 >>> integrate(x*y)
1533 Traceback (most recent call last):
1534 ...
1535 ValueError: specify integration variables to integrate x*y
1537 Note that ``integrate(x)`` syntax is meant only for convenience
1538 in interactive sessions and should be avoided in library code.
1540 >>> integrate(x**a*exp(-x), (x, 0, oo)) # same as conds='piecewise'
1541 Piecewise((gamma(a + 1), re(a) > -1),
1542 (Integral(x**a*exp(-x), (x, 0, oo)), True))
1544 >>> integrate(x**a*exp(-x), (x, 0, oo), conds='none')
1545 gamma(a + 1)
1547 >>> integrate(x**a*exp(-x), (x, 0, oo), conds='separate')
1548 (gamma(a + 1), re(a) > -1)
1550 See Also
1551 ========
1553 Integral, Integral.doit
1555 """
1556 doit_flags = {
1557 'deep': False,
1558 'meijerg': meijerg,
1559 'conds': conds,
1560 'risch': risch,
1561 'heurisch': heurisch,
1562 'manual': manual
1563 }
1564 integral = Integral(*args, **kwargs)
1566 if isinstance(integral, Integral):
1567 return integral.doit(**doit_flags)
1568 else:
1569 new_args = [a.doit(**doit_flags) if isinstance(a, Integral) else a
1570 for a in integral.args]
1571 return integral.func(*new_args)
1574def line_integrate(field, curve, vars):
1575 """line_integrate(field, Curve, variables)
1577 Compute the line integral.
1579 Examples
1580 ========
1582 >>> from sympy import Curve, line_integrate, E, ln
1583 >>> from sympy.abc import x, y, t
1584 >>> C = Curve([E**t + 1, E**t - 1], (t, 0, ln(2)))
1585 >>> line_integrate(x + y, C, [x, y])
1586 3*sqrt(2)
1588 See Also
1589 ========
1591 sympy.integrals.integrals.integrate, Integral
1592 """
1593 from sympy.geometry import Curve
1594 F = sympify(field)
1595 if not F:
1596 raise ValueError(
1597 "Expecting function specifying field as first argument.")
1598 if not isinstance(curve, Curve):
1599 raise ValueError("Expecting Curve entity as second argument.")
1600 if not is_sequence(vars):
1601 raise ValueError("Expecting ordered iterable for variables.")
1602 if len(curve.functions) != len(vars):
1603 raise ValueError("Field variable size does not match curve dimension.")
1605 if curve.parameter in vars:
1606 raise ValueError("Curve parameter clashes with field parameters.")
1608 # Calculate derivatives for line parameter functions
1609 # F(r) -> F(r(t)) and finally F(r(t)*r'(t))
1610 Ft = F
1611 dldt = 0
1612 for i, var in enumerate(vars):
1613 _f = curve.functions[i]
1614 _dn = diff(_f, curve.parameter)
1615 # ...arc length
1616 dldt = dldt + (_dn * _dn)
1617 Ft = Ft.subs(var, _f)
1618 Ft = Ft * sqrt(dldt)
1620 integral = Integral(Ft, curve.limits).doit(deep=False)
1621 return integral
1624### Property function dispatching ###
1626@shape.register(Integral)
1627def _(expr):
1628 return shape(expr.function)
1630# Delayed imports
1631from .deltafunctions import deltaintegrate
1632from .meijerint import meijerint_definite, meijerint_indefinite, _debug
1633from .trigonometry import trigintegrate