Revision history [back]

First, when asking questions, it might speed up a potential answer if you would care to give enough detail about your problem ... so, I only guess that this is what you want $$F(\mu,y) = \int_0^\infty \frac{x^{y}}{e^{x-\mu}+1} dx$$

Second, I'd be very careful about www-ing statements like

"...which can only be done numerically..."

since, to my recollection $F(\mu,y)$ is just another way of writing the polylog for all $y>-1$, i.e down to those values of $y$ where the singularity of the integrand remains harmless. Therefore, at least I do not understand your claim, since in 1D the DOS in a semiconductor is $\propto x^{-1/2}$ (which you also did not specify that in your question). So it can be done analytically ...

Third, if still you really want to do this numerically in sage ... just look here http://www.sagemath.org/doc/reference/calculus/sage/gsl/integration.html and e.g. do this

def F(y,mu):
    ul = 100+mu # increase the 100 if you really need
    return numerical_integral(x^y/(e^(x+mu)+1),0,ul)

F(-1/2,0)
(1.072154920510084, 8.286087609477639e-07)

Finally, and only to my very humble understanding: sage is a chain of existing tools, i.e. it is kind of as strong as it's weakest link. Since in plain sage things like Maxima are the link for many of those calculations which 90%-99% of non-mathematician might call 'symbolic', plain sage frequently is of little help if compared to MM or Maple once you come to such 'symbolic' problems.

Applied to your case this means:

sage: var('y,mu')
(y, mu)
sage: assume(y+1>0)
sage: assume(y,'integer')
sage: integral(x^y/(e^(x+mu)+1),x,0,oo)
integrate(x^y/(e^(mu + x) + 1), x, 0, +Infinity)

so, nothing ... or if you try special values

sage: integral(x^(-1/2)/(e^(x+mu)+1),x,0,oo)
integrate(1/((e^(mu + x) + 1)*sqrt(x)), x, 0, +Infinity)

sage: integral(x^1/(e^(x+mu)+1),x,0,oo)
polylog(2, -e^mu) + limit(1/2*x^2 - x*log(e^(mu + x) + 1) - polylog(2, -e^(mu + x)), x, +Infinity, minus)

again, not too helpful. Maybe there are 'tricks' to get further but I don't know.

If however you use MM then, out of the box:

sage: mathematica('Integrate[x^a/ (Exp[x + b] + 1), {x, 0, Infinity}, Assumptions -> {Re[a] > -1}]')
-(a*Gamma[a]*PolyLog[1 + a, -E^(-b)])

so MM 'knows' that your integral can be done analytically and provides the answer, and

sage: mathematica('-(a*Gamma[a]*PolyLog[1 + a, -E^(-b)])/. {a -> -1/2, b -> 0} // N')
1.072154929940191

which is the numerical value we got previously.

In conclusion: To answer if that implies you should better be using MM, there are more powerful experts around here.

First, when asking questions, it might speed up a potential answer if you would care to give enough detail about your problem ... so, I only guess that this is what you want $$F(\mu,y) $$F(y,\mu) = \int_0^\infty \frac{x^{y}}{e^{x-\mu}+1} \frac{x^{y}}{e^{x+\mu}+1} dx$$

Second, I'd be very careful about www-ing statements like

"...which can only be done numerically..."

since, to my recollection $F(\mu,y)$ is just another way of writing the polylog for all $y>-1$, i.e down to those values of $y$ where the singularity of the integrand remains harmless. Therefore, at least I do not understand your claim, since in 1D the DOS in a semiconductor is $\propto x^{-1/2}$ (which you also did not specify that in your question). So it can be done analytically ...

Third, if still you really want to do this numerically in sage ... just look here http://www.sagemath.org/doc/reference/calculus/sage/gsl/integration.html and e.g. do this

def F(y,mu):
    ul = 100+mu # increase the 100 if you really need
    return numerical_integral(x^y/(e^(x+mu)+1),0,ul)

F(-1/2,0)
(1.072154920510084, 8.286087609477639e-07)

Finally, and only to my very humble understanding: sage is a chain of existing tools, i.e. it is kind of as strong as it's weakest link. Since in plain sage things like Maxima are the link for many of those calculations which 90%-99% of non-mathematician might call 'symbolic', plain sage frequently is of little help if compared to MM or Maple once you come to such 'symbolic' problems.

Applied to your case this means:

sage: var('y,mu')
(y, mu)
sage: assume(y+1>0)
sage: assume(y,'integer')
sage: integral(x^y/(e^(x+mu)+1),x,0,oo)
integrate(x^y/(e^(mu + x) + 1), x, 0, +Infinity)

so, nothing ... or if you try special values

sage: integral(x^(-1/2)/(e^(x+mu)+1),x,0,oo)
integrate(1/((e^(mu + x) + 1)*sqrt(x)), x, 0, +Infinity)

sage: integral(x^1/(e^(x+mu)+1),x,0,oo)
polylog(2, -e^mu) + limit(1/2*x^2 - x*log(e^(mu + x) + 1) - polylog(2, -e^(mu + x)), x, +Infinity, minus)

again, not too helpful. Maybe there are 'tricks' to get further but I don't know.

If however you use MM then, out of the box:

sage: mathematica('Integrate[x^a/ (Exp[x + b] + 1), {x, 0, Infinity}, Assumptions -> {Re[a] > -1}]')
-(a*Gamma[a]*PolyLog[1 + a, -E^(-b)])

so MM 'knows' that your integral can be done analytically and provides the answer, and

sage: mathematica('-(a*Gamma[a]*PolyLog[1 + a, -E^(-b)])/. {a -> -1/2, b -> 0} // N')
1.072154929940191

which is the numerical value we got previously.

In conclusion: To answer if that implies you should better be using MM, there are more powerful experts around here.

First, when asking questions, it might speed up a potential answer if you would care to give enough detail about your problem ... so, I only guess that this is what you want $$F(y,\mu) = \int_0^\infty \frac{x^{y}}{e^{x+\mu}+1} dx$$

Second, I'd be very careful about www-ing statements like

"...which can only be done numerically..."

since, to my recollection $F(\mu,y)$ $F(y,\mu)$ is just another way of writing the polylog for all $y>-1$, i.e down to those values of $y$ where the singularity of the integrand remains harmless. Therefore, at least I do not understand your claim, since in 1D the DOS in a semiconductor is $\propto x^{-1/2}$ (which you also did not specify that in your question). So it can be done analytically ...

Third, if still you really want to do this numerically in sage ... just look here http://www.sagemath.org/doc/reference/calculus/sage/gsl/integration.html and e.g. do this

def F(y,mu):
    ul = 100+mu # increase the 100 if you really need
    return numerical_integral(x^y/(e^(x+mu)+1),0,ul)

F(-1/2,0)
(1.072154920510084, 8.286087609477639e-07)

Finally, and only to my very humble understanding: sage is a chain of existing tools, i.e. it is kind of as strong as it's weakest link. Since in plain sage things like Maxima are the link for many of those calculations which 90%-99% of non-mathematician might call 'symbolic', plain sage frequently is of little help if compared to MM or Maple once you come to such 'symbolic' problems.

Applied to your case this means:

sage: var('y,mu')
(y, mu)
sage: assume(y+1>0)
sage: assume(y,'integer')
sage: integral(x^y/(e^(x+mu)+1),x,0,oo)
integrate(x^y/(e^(mu + x) + 1), x, 0, +Infinity)

so, nothing ... or if you try special values

sage: integral(x^(-1/2)/(e^(x+mu)+1),x,0,oo)
integrate(1/((e^(mu + x) + 1)*sqrt(x)), x, 0, +Infinity)

sage: integral(x^1/(e^(x+mu)+1),x,0,oo)
polylog(2, -e^mu) + limit(1/2*x^2 - x*log(e^(mu + x) + 1) - polylog(2, -e^(mu + x)), x, +Infinity, minus)

again, not too helpful. Maybe there are 'tricks' to get further but I don't know.

If however you use MM then, out of the box:

sage: mathematica('Integrate[x^a/ (Exp[x + b] + 1), {x, 0, Infinity}, Assumptions -> {Re[a] > -1}]')
-(a*Gamma[a]*PolyLog[1 + a, -E^(-b)])

so MM 'knows' that your integral can be done analytically and provides the answer, and

sage: mathematica('-(a*Gamma[a]*PolyLog[1 + a, -E^(-b)])/. {a -> -1/2, b -> 0} // N')
1.072154929940191

which is the numerical value we got previously.

In conclusion: To answer if that implies you should better be using MM, there are more powerful experts around here.