scipy.integrate.tplquad¶
- 
scipy.integrate.tplquad(func, a, b, gfun, hfun, qfun, rfun, args=(), epsabs=1.49e-08, epsrel=1.49e-08)[source]¶
- Compute a triple (definite) integral. - Return the triple integral of - func(z, y, x)from- x = a..b,- y = gfun(x)..hfun(x), and- z = qfun(x,y)..rfun(x,y).- Parameters
- funcfunction
- A Python function or method of at least three variables in the order (z, y, x). 
- a, bfloat
- The limits of integration in x: a < b 
- gfunfunction or float
- The lower boundary curve in y which is a function taking a single floating point argument (x) and returning a floating point result or a float indicating a constant boundary curve. 
- hfunfunction or float
- The upper boundary curve in y (same requirements as gfun). 
- qfunfunction or float
- The lower boundary surface in z. It must be a function that takes two floats in the order (x, y) and returns a float or a float indicating a constant boundary surface. 
- rfunfunction or float
- The upper boundary surface in z. (Same requirements as qfun.) 
- argstuple, optional
- Extra arguments to pass to func. 
- epsabsfloat, optional
- Absolute tolerance passed directly to the innermost 1-D quadrature integration. Default is 1.49e-8. 
- epsrelfloat, optional
- Relative tolerance of the innermost 1-D integrals. Default is 1.49e-8. 
 
- Returns
- yfloat
- The resultant integral. 
- abserrfloat
- An estimate of the error. 
 
 - See also - quad
- Adaptive quadrature using QUADPACK 
- quadrature
- Adaptive Gaussian quadrature 
- fixed_quad
- Fixed-order Gaussian quadrature 
- dblquad
- Double integrals 
- nquad
- N-dimensional integrals 
- romb
- Integrators for sampled data 
- simps
- Integrators for sampled data 
- ode
- ODE integrators 
- odeint
- ODE integrators 
- scipy.special
- For coefficients and roots of orthogonal polynomials 
 - Examples - Compute the triple integral of - x * y * z, over- xranging from 1 to 2,- yranging from 2 to 3,- zranging from 0 to 1.- >>> from scipy import integrate >>> f = lambda z, y, x: x*y*z >>> integrate.tplquad(f, 1, 2, lambda x: 2, lambda x: 3, ... lambda x, y: 0, lambda x, y: 1) (1.8750000000000002, 3.324644794257407e-14) 
