mp de sousa - 8 months ago 57

Python Question

I try to integrate this :

integrate(integrate(integrate(2*sin(z)*cos(atan((2*cos(y)-0.5+x)/(2*sin(y)))),y,0,pi/2),x,0,1),z,0,pi/2);

Wolfram find the solution but I would like to control the accuracy. I try with tplquad but there is some error.

`def f(x,y,z):`

return 2*sin(z)*cos(atan((2*cos(y)-0.5+x)/(2*sin(y))))

tplquad(f,0,1,0,pi/2,0,pi/2)

The errors are:

Blockquote

File "", line 3, in

File "/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.py", line 526, in tplquad

return dblquad(_infunc2,a,b,gfun,hfun,(func,qfun,rfun,args),epsabs=epsabs,epsrel=epsrel)

File "/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.py", line 461, in dblquad

return quad(_infunc,a,b,(func,gfun,hfun,args),epsabs=epsabs,epsrel=epsrel)

File "/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.py", line 281, in quad

retval = _quad(func,a,b,args,full_output,epsabs,epsrel,limit,points)

File "/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.py", line 345, in _quad

return _quadpack._qagse(func,a,b,args,full_output,epsabs,epsrel,limit)

File "/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.py", line 406, in _infunc

a = gfun(x)

Have you an idea where the error could come from ?

Answer

The documentation of `tplquad`

states that the integration limits of the inner integrals should be provided as functions of the outer integration variables (even if they happen to be constants as in your case).

The correct usage of `tplquad`

in your case is shown below. Note that the ordering of the arguments in the definition of `f`

should correspond to the ordering of integrations. The first (last) argument of `f`

is the last (first) to be integrated. In this case the ordering is irrelevant due to fixed integration limits.

```
import numpy as np
from scipy.integrate import tplquad
def f(y,x,z):
return 2*np.sin(z)*np.cos(np.arctan((2*np.cos(y)-0.5+x)/(2*np.sin(y))))
tplquad(f,0,np.pi/2, lambda z: 0, lambda z:1, lambda z, x: 0, lambda z, x: np.pi/2)
```

`(1.9999999999999998, 2.492629060475153e-14)`

Source (Stackoverflow)