Ask Your Question

numerical_integral (GSL C library) does not work with complex integrands?

asked 2018-03-12 10:30:59 -0500

joaoff gravatar image

I am unable to use numerical_integral (GSL C library) with a complex integrand. I haven't found anything about this issue in the documentation. Is this correct? nintegral works correctly.


f(x) = (exp(i*x) + exp(-i*x))/2

the above code outputs

(0.8414709848078965, 9.34220461887732e-15, 21, 0)
Error in lines 3-3
TypeError: unable to coerce to a real number

Tested with SageMath version 7.5.1 and Cocalc.


edit retag flag offensive close merge delete

1 answer

Sort by ยป oldest newest most voted

answered 2018-03-12 13:14:58 -0500

dan_fulea gravatar image

In the documentation of numerical_integral there is no explicit statement that the function should be (defined on some proper or improper interval) with values in $\mathbb R$, but the following sentences in the doc string strongly suggest this is the case:

Note the coercion to the real field RR, which prevents underflow:

      sage: f = exp(-x**2)
      sage: numerical_integral(f, -Infinity, +Infinity)  # random output
      (1.7724538509060035, 3.4295192165889879e-08)


One can integrate any real-valued callable function:

      sage: numerical_integral(lambda x: abs(zeta(x)), [1.1,1.5])  # random output
      (1.8488570602160455, 2.052643677492633e-14)

In our case the solution / a work-around is simple, we just isolate the imaginary and the real parts of the function to be integrated.

sage: f(x) = (exp(i*x) + exp(-i*x))/2
sage: numerical_integral( real(f), (0,1) )
(0.8414709848078964, 9.342204618877319e-15)
sage: numerical_integral( imag(f), (0,1) )
(0.0, 0.0)

which matches

sage: numerical_integral( cos(x), (0,1) )
(0.8414709848078964, 9.342204618877319e-15)

Further examples with different, but similar data:

sage: numerical_integral( lambda x: exp(i*x).real(), (0, pi/6) )
(0.4999999999999999, 5.551115123125781e-15)
sage: numerical_integral( lambda x: exp(i*x).imag(), (0, pi/6) )
(0.13397459621556138, 1.487416814333746e-15)
edit flag offensive delete link more



In general, the integral of a complex function depends on the integration path if the function has any singularities ; see this for a smalll introduction to a vast subject.

This answer assumes that the integration path is the segment $[0~1]$ of $\mathbb{R}$. What does justify this assumption ?

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2018-03-12 13:46:33 -0500 )edit

Thank you for the replies. In fact, I was trying to optimize my code finding a replacer for nintegral(), that works for me, but seems to be much slower than the equivalent function in Mathematica.

joaoff gravatar imagejoaoff ( 2018-03-13 04:45:17 -0500 )edit

@Emmanuel Charpentier The question also assumes the integration path is the segment $[0,1]$ in $\mathbb R$. Moreover, the integrand is entire, so given this example and still the question, the only point of the question could have been the functionality of numerical_integral. If we really need to understand monodromy, we can still reduce the numerical calculus to the above:

def myint(f, P=0, eps=0.1):
    var( 't' );
    path(t)  = P + eps*exp(i*t)
    g(t) = f( path(t) ) * diff( path(t), t )
    return ( numerical_integral( lambda u: real(g(u)), (0, 2*pi) )[0]
             numerical_integral( lambda u: imag(g(u)), (0, 2*pi) )[0] * i  )

Then myint(log), myint(log,0,1) shows "it".

dan_fulea gravatar imagedan_fulea ( 2018-03-13 18:47:37 -0500 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower


Asked: 2018-03-12 10:30:59 -0500

Seen: 35 times

Last updated: Mar 12