Reputation:
I'm trying to integrate a piecewise function using Sagemath, and finding it to be impossible. My original code is below, but it's wrong due to accidental evaluation described here.
def f(x):
if(x < 0):
return 3 * x + 3
else:
return -3 * x + 3
g(x) = integrate(f(t), t, 0, x)
The fix for plotting mentioned on the website is to use f
instead of f(t)
, but this is apparently unsupported for the integrate()
function since a TypeError
is raised.
Is there a fix for this that I'm unaware of?
Upvotes: 5
Views: 1821
Reputation: 3965
In newer SageMath versions, use piecewise
(lowercase) instead.
sage: f = piecewise([((-oo, 0), 1), ((0, 1), 2), ((1, oo), 3)], var=x)
sage: f
piecewise(x|-->1 on (-oo, 0), x|-->2 on (0, 1), x|-->3 on (1, +oo); x)
sage: g = integrate(f, x)
sage: g
piecewise(x|-->x on (-oo, 0), x|-->2*x on (0, 1), x|-->3*x - 1 on (1, +oo); x)
sage: plot(g, (x, -2, 2))
Launched png viewer for Graphics object consisting of 1 graphics primitive
Apart from that, you also have sgn
, abs
and heaviside
which you can use to build up logical expressions.
sage: f = heaviside(x)*1 + (1-heaviside(x))*2
sage: f
-heaviside(x) + 2
sage: g = integrate(f, x)
Warning: piecewise indefinite integration does not return a continuous antiderivative
sage: g
-x*heaviside(x) + 2*x
sage: plot(g, (x, -2, 2))
Launched png viewer for Graphics object consisting of 1 graphics primitive
Ideally SageMath could expose Maxima's charfun
(https://stackoverflow.com/a/24915794), but it currently doesn't. integrate
mostly relies on Maxima anyway.
Upvotes: 0
Reputation:
In addition to using the Piecewise class, this can easily be fixed by defining g(x)
as a Python function as well:
def f(x):
if(x < 0):
return 3 * x + 3
else:
return -3 * x + 3
def g(x):
(y, e) = integral_numerical(f, 0, x)
return y
Then plot(g)
works just fine.
Upvotes: 1
Reputation:
Instead of defining a piecewise function via def
, use the built-in piecewise class:
f = Piecewise([[(-infinity, 0), 3*x+3],[(0, infinity), -3*x+3]])
f.integral()
Output:
Piecewise defined function with 2 parts, [[(-Infinity, 0), x |--> 3/2*x^2 + 3*x], [(0, +Infinity), x |--> -3/2*x^2 + 3*x]]
The piecewise functions have their own methods, such as .plot()
. Plotting does not support infinite intervals, though. A plot can be obtained with finite intervals
f = Piecewise([[(-5, 0), 3*x+3],[(0, 5), -3*x+3]])
g = f.integral()
g.plot()
But you also want to subtract g(0) from g. This is not as straightforward as g-g(0), but not too bad, either: get the list of pieces with g.list()
, subtract g(0) from each function, then recombine.
g0 = Piecewise([(piece[0], piece[1] - g(0)) for piece in g.list()])
g0.plot()
And there you have it:
By extending this approach, we don't even need to put finite intervals in f from the beginning. The following plots g - g(0) on a given interval [a,b], by modifying the domain:
a = -2
b = 3
g0 = Piecewise([((max(piece[0][0], a), min(piece[0][1], b)), piece[1] - g(0)) for piece in g.list()])
g.plot()
Upvotes: 3