ASKSAGE: Sage Q&A Forum - Latest question feedhttp://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Tue, 11 Dec 2018 09:47:56 -0600Integration on a region defined by an inequalityhttp://ask.sagemath.org/question/44636/integration-on-a-region-defined-by-an-inequality/ Does Sage or Maxima have a routine for -directly- calculating an integral over a region defined implicitly by an inequality. Suppose region R is defined by "f(x,y) less than c" and we want to find "integral g(x,y) over R". Does Sage have a procedure for this? The required output will be a function of c.maesumiTue, 11 Dec 2018 09:47:56 -0600http://ask.sagemath.org/question/44636/what instance is integral from?http://ask.sagemath.org/question/43244/what-instance-is-integral-from/ When the result has `integrate` in it, I check for it as follows
sage: anti=integrate(1/(sqrt(x + 1)*sqrt(-x + 1) + 5), x)
integrate(1/(sqrt(x + 1)*sqrt(-x + 1) + 5), x)
sage: isinstance(anti.operator(), sage.symbolic.integration.integral.IndefiniteIntegral)
True
But the above does not work when the result is `integral` instead of `integrate`.
My question is, what instance is `integral` coming from?
sage: anti=integrate(cos(b*x + a)*cos_integral(d*x + c)/x,x, algorithm="fricas")
integral(cos(b*x + a)*cos_integral(d*x + c)/x, x)
sage: anti.operator()
integral
sage: isinstance(anti.operator(), sage.symbolic.integration.integral.IndefiniteIntegral)
False
I looked at http://doc.sagemath.org/html/en/reference/calculus/sage/symbolic/integration/integral.html
but still do not know how to check for `intergal` vs. `integrate`
Any suggestions?
thanks
--Nasser
NasserThu, 02 Aug 2018 22:26:05 -0500http://ask.sagemath.org/question/43244/Problem evaluating a very tiny integralhttp://ask.sagemath.org/question/39514/problem-evaluating-a-very-tiny-integral/Hi everyone,
I'm pretty new with sage. I was forced to change from MATLAB to Sage, because I was told Sage does approximate very tiny numbers better as it can work with sqrt(2) as sqrt(2) and not as the rational number approximating it.
Approximations are very important for my problem. I need to evaluate this integral
$$\sum_{c=1}^{d}\int_{\min(256c-0.5,\frac{y(y+1)}{2})}^{\min(256c+0.5,\frac{y(y+1)}{2})}\frac{1}{\sqrt{2\pi}}e^{-\frac{x^2}{2}}\, dx$$
Suppose d = 1, then this is simply the integral
$$ \int_{255.5}^{256.5}\frac{1}{\sqrt{2\pi}}e^{-\frac{x^2}{2}}\, dx$$
To evaluate this integral I wrote the following code
> T = RealDistribution('gaussian,1)
> print T.cum_distribution_function(256.5)-T.cum_distribution_function(255.5)
because the integral above is the same as the difference of the distribution function of a standard distributed gaussian random variable between the boundaries of the integral.
However, and you can check yourselves if you don't believe it, the result I get with sage is 0.
I guess that this is due to some approximation, which sage does. Indeed the exact value of the integral is pretty close (and with pretty I mean a lot) to 0. My problem is that I need to be able to have the exact value, because the sum of integrals I'm working with and my whole work behind this integral requires me to be very careful with those little tiny numbers.
I tried to use the function
>integrate
to deal with this problem, but something funny and apparently inexplicable happened when I was trying to use it.
To be precise I defined this code:
def VarianceOfvy(y):
temp1 = 0
temp2 = 0
for r in range(0,y+1):
for x in range(0,r+1):
temp1 = temp1 + (255/256)^x * 1/256 * (r-x)^2
for r in range(0,y+1):
for x in range(0,r+1):
temp2 = temp2 + ((255/256)^x * 1/256 * (r-x))^2
sigma = temp1 - temp2
return sqrt(sigma)
def Integerd(y):
b = y*(y+1)/2
d = 1
c = 0
while min((c+1)*256-0.5,b) == (c+1)*256-0.5:
d = d+1
c = c+1
return d
def Probabilityvynequiv(y):
var('c')
b = (y*(y+1))/2
sigma = 2
mu = 1
d = Integerd(y)
factor = 1/(integrate(1/sqrt(2*pi)*e^(-(x/sigma)^2/2),x,-oo,(b-mu)/sigma) - integrate(1/sqrt(2*pi)*e^(-(x)^2/2),x,-oo,(-mu)/sigma))
p = sum(factor*1/sigma*integrate(1/sqrt(2*pi)*e^((-x^2)/(2)),x,c*256+0.5,min((c+1)*256-0.5,b)),c,0,d)
return p
And if I let it run, the result I get is
1/2*(erf(255.75*sqrt(2)) - erf(128.25*sqrt(2)) + erf(127.75*sqrt(2)) - erf(0.25*sqrt(2)))/(3*erf(1/4*sqrt(2)) + 1)
which I assume is correct, and at least it tells me that Sage is able to read my code and output a result.
If I call the function VarianceOfvy(2), the result I get is 3/65536*sqrt(11133895), which is also correct. Now, if I'm changing the command
sigma = 2
with
sigma = VarianceOfvy(2)
and try to let the whole program run again, Sage is not able anymore to output a result.
I'm really lost and I don't know what to do. Could someone advise me and give me some hints on how to evaluate those tiny integrals, in such a way that I don't loose any approximation?Alessandro VerzasconiMon, 13 Nov 2017 08:03:58 -0600http://ask.sagemath.org/question/39514/Definite Integral Fails due to Runtime Errorhttp://ask.sagemath.org/question/39354/definite-integral-fails-due-to-runtime-error/Here is the code that I have been writing; it's meant to predict the moment of inertia of a perfect cylinder tilted and offset from its axis of rotation, among other things. However, the definite_integral command fails with an error thrown by the gcd function. I suspect the issue is due to Sage's inner workings; however, I do not have the time or know-how to fully debug this issue. Any help would be appreciated.
from sage.symbolic.integration.integral import definite_integral
from sage.calculus.integration import numerical_integral
dSpringStretchedFully=10.625
dSpringUnstretchedLength=7.875
dArmatureInches=8
dSpringAnchor=4.75
thetaArmatureDefault=95
dSpringPitch=dSpringUnstretchedLength/sin(thetaArmatureDefault*pi/180)*sin(pi/180*(180-thetaArmatureDefault-(180/pi*arcsin(dSpringAnchor*sin(thetaArmatureDefault*pi/180)/dSpringUnstretchedLength))))
fSpringPoundsPerInch=65.69
kNewtonsPerPound=4.44822
kMetersPerInch=0.0254
fSpringNewtonsPerInch=fSpringPoundsPerInch*kNewtonsPerPound
dArmatureMeters=dArmatureInches*kMetersPerInch
dHalfRopeInches=dArmatureInches*cos(pi/180*(thetaArmatureDefault-90))
dHalfRopeMeters=dHalfRopeInches*kMetersPerInch
fSpringNewtonsPerMeter=fSpringNewtonsPerInch/kMetersPerInch
fSpringMinimumLoad=0
kSpringDirection=1
rRopeInches=0.5
rRopeMeters=rRopeInches*kMetersPerInch
var('x')
var('dHalfRope')
var('thetaArmatureToRope')
thetaPointToAxis(x,thetaArmatureToRope)=thetaArmatureToRope-(180/pi*arctan((x-rRopeMeters)/dArmatureMeters))
dPointFromFarAxis(x,dHalfRope)=((dHalfRope^2)+(x-rRopeMeters)^2)^(1/2)
dPointFromNearAxis(x,dHalfRope,thetaArmatureToRope)=((dArmatureMeters^2)+((dPointFromFarAxis(x,dHalfRope))^2)-2*dPointFromFarAxis(x,dHalfRope)*dArmatureMeters*cos(pi/180*thetaPointToAxis(x,thetaArmatureToRope)))^(1/2)
IMultiplier(x)=2*(((rRopeMeters^2)-(x-rRopeMeters)^2)^(1/2))
ILine(x,dHalfRope,thetaArmatureToRope)=((dPointFromNearAxis(x,dHalfRope, thetaArmatureToRope))^2)*IMultiplier(x)
ISlice(dHalfRope,thetaArmatureToRope)=definite_integral(ILine(x,dHalfRope,thetaArmatureToRope),x,0,2*rRopeMeters)
ISlice(3,4)Benji_ZhangTue, 31 Oct 2017 11:17:36 -0500http://ask.sagemath.org/question/39354/Trouble with an integralhttp://ask.sagemath.org/question/37395/trouble-with-an-integral/Hi,
Let f be the following function:
f(t, z0, z1) = (z0*z1 - 1)*(z0 - z1)/((t*z0 - 1)*(t*z1 - 1)*(t - z0)*(t - z1)*z0)
f is a nice holomorphic function and I would like to compute the contour integral
\int_{|z1|=1} \int_{|z0|=1} f(t, z0, z1) dz0/z0 dz1/z1
Assuming 0 < t < 1, the result should be -4 pi^2/(t^2-1).
In Sage:
**var("z0, z1, theta0, theta1")**
**assume(0 < t < 1)**
**f = (z0*z1 - 1)*(z0 - z1)/((t*z0 - 1)*(t*z1 - 1)*(t - z0)*(t - z1)*z0)**
**g = f.subs({z0: exp(I*theta0), z1: exp(I*theta1)})**
Now if I integrate with respect to theta0 first: **g.integrate(theta0, 0, 2*pi)** Sage answers that the integral is zero.
If I integrate with respect to theta1 first: **factor(g.integrate(theta1, 0, 2*pi).integrate(theta0, 0, 2*pi))** Sage answers that the integral is -4*pi^2/t^2 which is also clearly wrong...
Maple finds the right answer.
What can I do to make Sage compute it right? (This is a test I need to compute much more complicated integrals after this so I have to make sure Sage gives me the right answer).
RoroMon, 24 Apr 2017 10:21:44 -0500http://ask.sagemath.org/question/37395/Integrate piecewise function with change of variablehttp://ask.sagemath.org/question/37114/integrate-piecewise-function-with-change-of-variable/I would like to integrate a piecewise defined function while operating a change of variable. I start by defining the function and another variable involved in the change of variable:
phi(x) = piecewise([([-1,1], (1-abs(x))*(1-abs(x))*(1+2*abs(x)))]);
phi(x) = phi.extension(0);
h=pi/n;
h=h.n();
What I would like to do is integrate the function `phi(x/h-1)` between `0` and `pi` so I try it and results in
integral(phi(x/h-1),x,0,pi)
ValueError: substituting the piecewise variable must result in real number
So I then try to use another variable which I try to define to be 'real'
t=var('t')
assume(t,'real');
integral(phi(t/h-1),t,0,pi)
but it results in the same error... Now I try the "lambda" method since it worked when calling the `plot` function with the same change of variable; but fail again
integral(lambda t: phi(t/h-1),t,0,pi)
TypeError: unable to convert <function <lambda> at 0x16d71f140> to a symbolic expression
Now I try to use another integration method with `definite_integral` but get the same errors, only different for the "lambda" method
definite_integral(lambda x: phi(x/h-1),x,0,pi)
TypeError: cannot coerce arguments: no canonical coercion from <type 'function'> to Symbolic Ring
Is there any way around this? I really do not know what else to try...
jrojasquTue, 28 Mar 2017 18:28:56 -0500http://ask.sagemath.org/question/37114/Sage Interactive has Attribute errorhttp://ask.sagemath.org/question/30132/sage-interactive-has-attribute-error/ Herr is my code for a double integrator interactive I am making that both make a double integral and shows a graph with added rectangles. However in my code I get an attribute error I don't know.
from sage.plot.plot3d.shapes import Box
x,y = var('x,y')
html("<h1>triple integrater<h1>")
permutations = ["dx dy","dy dx"]
@interact
def interplay(order= permutations,function= input_box(sin(x*y)),lower_x_bound= input_box(0.1.1),upper_x_bound = input_box(1),lower_y_bound=input_box(0),upper_y_bound=input_box(1),showGraph = checkbox(default = False), numrecs= input_box(50)):
if permutations == "dx dy":
result = integral(integral(function,x,lower_x_bound,upper_x_bound),y,lower_y_bound,upper_y_bound)+function(lower_x_bound,lower_y_bound)+function(lower_x_bound,lower_y_bound)#function.integrate(x,lower_x_bound,upper_x_bound).integrate(y,lower_y_bound,upper_y_bound)
q="$\int_%s^{%s} \int_%s^%s %s $" % (lower_y_bound, upper_y_bound,lower_x_bound,upper_x_bound,function)
# puts expressions inside latex
else:
result =integral(integral(function,y,lower_y_bound,upper_y_bound),y,lower_x_bound,upper_x_bound)
#actually calculates integral
q="$\int_%s^{%s} \int_%s^%s %s $" % (lower_x_bound,upper_x_bound,lower_y_bound,upper_y_bound,function)
#if type(lower_x_bound)!=float or type(upper_x_bound)!=float or type(lower_y_bound)!= float or type(upper_y_bound)!= float :
if lower_x_bound in RR or upper_x_bound in RR or lower_y_bound in RR or lower_y_bound in RR:
#check if bounds are are real before graphing because I do nto know how to graph with x as a bound.
showGraph = True
html("Can only grah numerical bounds")
else:
showGraph = False
if showGraph == True :
graph = plot3d(function,(x,lower_x_bound,upper_x_bound),(y,lower_y_bound,upper_y_bound),fill=True,color = "orange",spin = 4)
html("<h3>Graph of Integrated Region</h3>")
delta = upper_x_bound - lower_x_bound
B = sum([Box([1/numrecs,.5,abs(function(lower_x_bound+(i-1)*delta/numrecs , .2))], color="orange").translate((lower_x_bound+delta*i/numrecs,0,function(lower_x_bound+i*delta/numrecs)/2)) for i in [0..numrecs]])
#makes a whole bunch of rectangels which approximate graph of function being integrated
show(graph+B)
html("<h3>Numerical Result</h3>")
#s = "$\int_{2}^{3} \int_{4}^{5} {0} \,dxdy = {1} $"
#p = s.format(function, result,lower_x_bound,upper_x_bound,lower_y_bound,upper_y_bound)
html("%s" %q)
#shows the result
I am getting this error. I do not know what it means.
Error in lines 5-29
Traceback (most recent call last):
File "/projects/454b81d2-ef23-4b95-bd57-5c718a468ea1/.sagemathcloud/sage_server.py", line 881, in execute
exec compile(block+'\n', '', 'single') in namespace, locals
File "", line 2, in <module>
File "sage/structure/element.pyx", line 418, in sage.structure.element.Element.__getattr__ (/projects/sage/sage-6.9/src/build/cythonized/sage/structure/element.c:4670)
return getattr_from_other_class(self, P._abstract_element_class, name)
File "sage/structure/misc.pyx", line 259, in sage.structure.misc.getattr_from_other_class (/projects/sage/sage-6.9/src/build/cythonized/sage/structure/misc.c:1771)
raise dummy_attribute_error
AttributeError: 'sage.rings.real_mpfr.RealLiteral' object has no attribute 'gen'collabmathTue, 20 Oct 2015 01:35:13 -0500http://ask.sagemath.org/question/30132/Is this a known bug with integral()http://ask.sagemath.org/question/30075/is-this-a-known-bug-with-integral/I've tried to compute the following integral wth integral() in a SageMathCloud worksheet: $\displaystyle \int_{-\pi/6}^{\pi/6}\frac{\cos x}{1+\sin x}dx$.
The output was an error message (saying the integral is divergent), just like the one I got in SageMathCell (see link):
https://sagecell.sagemath.org/?z=eJzzVLBVyMwrSU0vSszRSM4v1qjQ1Ncw1C7OzAOyNHUUKnQUdAsy9c10FECkJi9XcUZ-uYanJgDa5Q_i&lang=sage
So I tried with integrate() and with numerical_integral() as well. I was never able to obtain the value of this integral, which turns out to be $\ln(3)$ after an obvious substitution.
Is this a bug?
Note that replacing 1 by 1.1 yields this:
https://sagecell.sagemath.org/?z=eJzzVLBVyMwrSU0vSszRSM4v1qjQ1Ncw1DPULs7MA7I1dRQqdBR0CzL1zXQUQKQmL1dxRn65hqcmAPXGEEE=&lang=sage
while we get that when replacing 1 by 2:
https://sagecell.sagemath.org/?z=eJzzVLBVyMwrSU0vSszRSM4v1qjQ1Ncw0i7OzAOyNHUUKnQUdAsy9c10FECkJi9XcUZ-uYanJgDbCA_j&lang=sageJulienSat, 17 Oct 2015 10:01:02 -0500http://ask.sagemath.org/question/30075/Simplify result of this definite integralhttp://ask.sagemath.org/question/10503/simplify-result-of-this-definite-integral/It is well known that for $n\in\mathbb{N}$ and $n>0$ (an maybe even for more than these restrictions): $$I_n = \int_0^\infty\frac{x^n}{e^x-1}dx = \zeta(n+1)n!$$ which, analytically can be shown easily by expanding $1/(1-e^{-x})$ into a geometric series, which leads to trivial integrals, and by using $\zeta(n+1)=\sum_{l=1}^\infty l^{-(n+1)}$. So, eg.: $$I_1 = \pi^2/6$$ $$I_2 = 2\zeta(3)$$ a.s.o...
Now, if I try even the simplest case with sage, I get this 'nifty' little results
sage: integrate(x/(exp(x)-1),x,0,oo)
-1/6*pi^2 + limit(-1/2*x^2 + x*log(-e^x + 1) + polylog(2, e^x), x,
+Infinity, minus)
**Is there any trick to simplify this** down to the final result, or is this about as far as I can get with sage alone?
PS.: it is probably needless to say that (once again ... :( ...)
In[1]:= Integrate[x/(Exp[x] - 1), {x, 0, Infinity}]
Out[1]:= Pi^2/6
MarkWed, 04 Sep 2013 10:06:49 -0500http://ask.sagemath.org/question/10503/Solving this DE containing an integralhttp://ask.sagemath.org/question/28643/solving-this-de-containing-an-integral/I am trying to solve
$$ 0 = - \partial_a F(a)-e(a)F(a) + \int_0^{a} e(a')\partial_a F(a') d a' + n $$
optimally, without specifying e(a). But if necessary (as I guess), e(a) = k1*exp(k2*a). Here's my code:
var('a b k_1 k_2 n')
e(b) = k_1*exp(b*k_2)
F = function('F', a)
g(b) = e(b)*diff(F,b,1)
de = -diff(F, a, 1) - e(a)*F(a) + g(b).integral(b, 0, a) + n
y = desolve(de, F, ivar=a); y
And the output is
(_C + n*Ei(k_1*e^(a*k_2)/k_2)/k_2)*e^(-k_1*e^(a*k_2)/k_2)
However, I believe that something is wrong with my integral, it is not doing what I expect it to do. For example, if I change the integration boundary to 1:
de = -diff(F, a, 1) - e(a)*F(a) + g(b).integral(b, 0, a) + n
I will still get the same result:
(_C + n*Ei(k_1*e^(a*k_2)/k_2)/k_2)*e^(-k_1*e^(a*k_2)/k_2)
Where am I going wrong?FooBarSun, 19 Jul 2015 07:13:44 -0500http://ask.sagemath.org/question/28643/Problem with integralhttp://ask.sagemath.org/question/27235/problem-with-integral/Considering my last question ('Problem with hypergeometric') I tried to replace my Ansatz with a different formula. Here what happened:
F = lambda z: (2/pi)*integral((4*cos(x)^2-1)^z*sin(x)^2,x,0,pi)
print F(1/2)
RuntimeError: ECL says: Error executing code in Maxima:
With Maple:
F := z -> (2/Pi)*int((4*cos(x)^2-1)^z*sin(x)^2,x=0..Pi);
evalf(F(1/2)); .3697166864+.4838437554*I
Oh mei oh mei ...Peter LuschnyMon, 29 Jun 2015 07:59:58 -0500http://ask.sagemath.org/question/27235/Another problem with integralhttp://ask.sagemath.org/question/27237/another-problem-with-integral/This is a follow-up to 'Problem with integral' which is a follow-up of 'Problem with hypergeometric'. I am still trying to solve the same basic problem with Sage and I am just reporting what I am experiencing.
F = lambda z: (1/pi)*integral((x-1)^z*sqrt(1/x-1/4), x,0,4)
print F(1/2).n()
NaN
print F(3/2).n()
0.509025648974361 Exception ValueError: ValueError('negative number to a fractional power not real',) in 'sage.gsl.integration.c_ff' ignored [The error message some 20 times.]
With Maple:
F := z -> (1/Pi)*int((x-1)^z*sqrt(1/x-1/4), x=0..4);
evalf(F(1/2)); 0.3697166867 + 0.4838248688*I
evalf(F(3/2)); 0.5090256475 - 0.3669993270 IPeter LuschnyMon, 29 Jun 2015 08:21:52 -0500http://ask.sagemath.org/question/27237/integral should not be zerohttp://ask.sagemath.org/question/24532/integral-should-not-be-zero/
F = sqrt((cos(x) - 1)^2 + sin(x)^2)
F.integrate(x, 0, 2*pi)
yields 0, the expected answer is 8.
jllbSun, 19 Oct 2014 17:44:26 -0500http://ask.sagemath.org/question/24532/What does "Runtime Error: ECL says: ... is not of type FIXNUM" mean and how to fix it?http://ask.sagemath.org/question/10908/what-does-runtime-error-ecl-says-is-not-of-type-fixnum-mean-and-how-to-fix-it/I am trying to integrate $x(0.6x^{0.5}+0.6)$ from $0$ to 1.
Here's what I tried:
var('x')
(x*(0.6*x^0.5+0.6)).integrate(x,0,1)
But I received the following error:
RuntimeError: ECL says: 3.0 is not of type FIXNUM.
What does the error mean (I am especially curious where the 3 came from) and what can I do to correctly evaluate the definite integral? According to Wolfram Alpha, the answer is 0.54.ensabaMon, 13 Jan 2014 03:42:37 -0600http://ask.sagemath.org/question/10908/symbolic integrationhttp://ask.sagemath.org/question/10724/symbolic-integration/Hi,
I am a sage-noob. To calculate the vortex sheet of a flat plate in potential flow, I have to solve some integrales of such a king:
**Mathematica syntax:** Integrate[Exp[-I * t * 2 * m * k]*((1/t+1)^0.5-1),t,0, Infinity]
And this is the reult:
If[Im[k m] < 0,
((0. + 0.5 I) - (0. + 0.886227 I) HypergeometricU[-0.5, 0, (2 I) k m])/km
My approach in **Sagemath** is:
f(x) = exp(-I*x*2*m*k)*((1/x+1)^0.5-1)
integral(f, x, 0, infinity)
But Sagemath is asking for assumptions regarding m and k. Im[k m] < 0 would be the appropriate assumption. What would be the right syntax? Thank you very much for your answers.
Best wishes, Chris
chris42Mon, 11 Nov 2013 02:24:05 -0600http://ask.sagemath.org/question/10724/Fermi-Dirac integral of half orderhttp://ask.sagemath.org/question/10517/fermi-dirac-integral-of-half-order/I am trying to implement a 1D model for semiconductor pn-junctions. This involves evaluating the Fermi-Dirac integral of half order, which can only be done numerically. I saw that the GSL library has an algorithm to perform this integration, but I can't figure out how to call it in Sage, is this possible? If not, what are some alternatives to evaluate this integral? Thanks in advance.wlp2Fri, 06 Sep 2013 07:20:20 -0500http://ask.sagemath.org/question/10517/How to integrate in 2D, along the locus of a line.http://ask.sagemath.org/question/10348/how-to-integrate-in-2d-along-the-locus-of-a-line/If I have a 2D continuous and differentiable function (f(x,y)) and a general 2D line L (ax + by + c = 0), how do I evaluate the definite integral of f with respect to x and y along the line L between an arbitrary start and end point that both lie on the line?mmsood99Thu, 11 Jul 2013 12:56:46 -0500http://ask.sagemath.org/question/10348/Change of variable in an integrationhttp://ask.sagemath.org/question/8679/change-of-variable-in-an-integration/How to indicate a change of variable to Sage in an integration when Sage seems clueless?Green diodSat, 28 Jan 2012 11:59:27 -0600http://ask.sagemath.org/question/8679/