ZeroDivisionError in SymPy definite integral; variable integration bounds

1
votes

I'm using Julia and calling SymPy to integrate a bivariate normal PDF over its two variables, a and n. There is a third variable p that I want to include in the integration bounds, but SymPy throws a ZeroDivisionError when I try this.

Running the following,

using SymPy
σ = 0.2
ρ = -0.9
f(a,n) = (1/(2*π*(σ^2)*sqrt(1-ρ^2)))*exp(-(1/(2*(1-ρ^2)))*((a^2)/σ^2 + (n^2)/σ^2 -((2*ρ*a*n)/σ^2)))
@vars a n p
integrate(a*f(a,n), (n, -oo, p*a), (a, -oo, oo))

Yields

ERROR: PyError ($(Expr(:escape, :(ccall(#= C:\Users\Anshu\.julia\packages\PyCall\zqDXB\src\pyfncall.jl:43 =# @pysym(:PyObject_Call), PyPtr, (PyPtr, PyPtr, PyPtr), o, pyargsptr, kw))))) <class 'ZeroDivisionError'>

Weirdly, the following attempts run with no issues:

Just integrating f, not a*f:

integrate(f(a,n), (n, -oo, p*a), (a, -oo, oo))

Replacing p with a constant 1:

integrate(a*f(a,n), (n, -oo, a), (a, -oo, oo))
1
pythonjuliasympyintegral

1 Answers

1
votes

After defining sigma and rho, I get:

julia> integrate(a*f(a,n), (n, -oo, p*a), (a, -oo, oo)) |> string
"0.159154943091895*Integral(a*exp(-a^2/(-2*ρ^2*σ^2 + 2*σ^2))*Integral(exp(-n^2/(-2*ρ^2*σ^2 + 2*σ^2))*exp(2*a*n*ρ/(-2*ρ^2*σ^2 + 2*σ^2)), (n, -oo, a*p)), (a, -oo, oo))/(σ^2*sqrt(-(ρ - 1)*(ρ + 1)))"

Perhaps it is your sympy version, I have

julia> sympy.__version__
"1.5.1"