Opened 8 years ago

# Numerical approximation of a divergent integral

Reported by: Owned by: eviatarbach burcin major sage-6.4 calculus N/A

Sage is numerically approximating this integral, even though it's divergent:

```sage: integrate(x^3/sqrt(x^7+1), x, 1, oo).n()
-2.0585298599983344
```

It seems that if we don't allow Maxima to detect its divergence (`numerical_integral` passes it directly to GSL), GSL will also fail on simpler divergent integrals:

```sage: numerical_integral(1/x, 1, oo)
(65.94931131932763, 8.156214940519742)
sage: numerical_integral(x,1,oo)
(-0.4999999993521234, 1.3615531480049015e-09)
```

```sage: numerical_integral(e^(-x)/x,0,oo)
(37.191280375549404, 6.239196965189217)
```

### comment:1 Changed 8 years ago by eviatarbach

• Description modified (diff)

### comment:2 Changed 8 years ago by kcrisman

Regarding the first one; basically, when we return a noun form

```sage: integrate(x^3/sqrt(x^7+1), x, 1, oo)
integrate(x^3/sqrt(x^7 + 1), x, 1, +Infinity)
```

and call `n` we do

```sage: A._convert(RR)
-2.0585298599985333
```

because we do

```sage: sage.symbolic.integration.integral.DefiniteIntegral._evalf_??

from sage.gsl.integration import numerical_integral
# The gsl routine returns a tuple, which also contains the error.
# We only return the result.
return numerical_integral(f, a, b)[0]
```

So these are both manifestations of the same thing.

So... is it user error to numerically integrate a divergent integral? I certainly don't know that we should be checking every numerical integral for divergence, particularly since Maxima apparently can't (yet) do the first one in any case!

### comment:3 Changed 8 years ago by eviatarbach

Maybe the options being passed to GSL could be changed? It seems absurd that it should give a numerical answer for `numerical_integral(x^3,1,oo)`, for example. It apparently has an error for divergence (http://www.gnu.org/software/gsl/manual/html_node/Numerical-integration-error-codes.html).

### comment:4 Changed 8 years ago by eviatarbach

In fact, with the code at http://www.physics.ohio-state.edu/~ntg/780/gsl_examples/qagiu_test.cpp, it returns errors for all these integrals.

### comment:5 Changed 8 years ago by kcrisman

If you present a patch, I think we'd be very interested in reviewing it. Silly to return nonsense in these cases.

### comment:6 Changed 8 years ago by jdemeyer

• Milestone changed from sage-5.11 to sage-5.12

### comment:7 Changed 8 years ago by vbraun_spam

• Milestone changed from sage-6.1 to sage-6.2

### comment:8 Changed 7 years ago by vbraun_spam

• Milestone changed from sage-6.2 to sage-6.3

### comment:9 Changed 7 years ago by vbraun_spam

• Milestone changed from sage-6.3 to sage-6.4

### comment:10 Changed 7 years ago by rws

Edit: sorry wrong ticket

Last edited 7 years ago by rws (previous) (diff)

### comment:11 Changed 6 years ago by tmonteil

• Description modified (diff)
Note: See TracTickets for help on using tickets.