Piecewise-defined Functions¶
This module implement piecewise functions in a single variable. See
sage.sets.real_set
for more information about how to construct
subsets of the real line for the domains.
EXAMPLES:
sage: f = piecewise([((0,1), x^3), ([-1,0], -x^2)]); f
piecewise(x|-->x^3 on (0, 1), x|-->-x^2 on [-1, 0]; x)
sage: 2*f
2*piecewise(x|-->x^3 on (0, 1), x|-->-x^2 on [-1, 0]; x)
sage: f(x=1/2)
1/8
sage: plot(f) # not tested
TODO:
- Implement max/min location and values,
AUTHORS:
- David Joyner (2006-04): initial version
- David Joyner (2006-09): added __eq__, extend_by_zero_to, unextend, convolution, trapezoid, trapezoid_integral_approximation, riemann_sum, riemann_sum_integral_approximation, tangent_line fixed bugs in __mul__, __add__
- David Joyner (2007-03): adding Hann filter for FS, added general FS filter methods for computing and plotting, added options to plotting of FS (eg, specifying rgb values are now allowed). Fixed bug in documentation reported by Pablo De Napoli.
- David Joyner (2007-09): bug fixes due to behaviour of SymbolicArithmetic
- David Joyner (2008-04): fixed docstring bugs reported by J Morrow; added support for Laplace transform of functions with infinite support.
- David Joyner (2008-07): fixed a left multiplication bug reported by C. Boncelet (by defining __rmul__ = __mul__).
- Paul Butler (2009-01): added indefinite integration and default_variable
- Volker Braun (2013): Complete rewrite
- Ralf Stephan (2015): Rewrite of convolution() and other calculus functions; many doctest adaptations
TESTS:
sage: fast_callable(f, vars=[x])(0.5)
0.125000000000...
-
class
sage.functions.piecewise.
PiecewiseFunction
¶ Bases:
sage.symbolic.function.BuiltinFunction
Piecewise function
EXAMPLES:
sage: var('x, y') (x, y) sage: f = piecewise([((0,1), x^2*y), ([-1,0], -x*y^2)], var=x); f piecewise(x|-->x^2*y on (0, 1), x|-->-x*y^2 on [-1, 0]; x) sage: f(1/2) 1/4*y sage: f(-1/2) 1/2*y^2
-
class
EvaluationMethods
¶ Bases:
object
-
convolution
(self, parameters, variable, other)¶ Return the convolution function, \(f*g(t)=\int_{-\infty}^\infty f(u)g(t-u)du\), for compactly supported \(f,g\).
EXAMPLES:
sage: x = PolynomialRing(QQ,'x').gen() sage: f = piecewise([[[0,1],1]]) ## example 0 sage: g = f.convolution(f); g piecewise(x|-->x on (0, 1], x|-->-x + 2 on (1, 2]; x) sage: h = f.convolution(g); h piecewise(x|-->1/2*x^2 on (0, 1], x|-->-x^2 + 3*x - 3/2 on (1, 2], x|-->1/2*x^2 - 3*x + 9/2 on (2, 3]; x) sage: f = piecewise([[(0,1),1],[(1,2),2],[(2,3),1]]) ## example 1 sage: g = f.convolution(f) sage: h = f.convolution(g); h piecewise(x|-->1/2*x^2 on (0, 1], x|-->2*x^2 - 3*x + 3/2 on (1, 3], x|-->-2*x^2 + 21*x - 69/2 on (3, 4], x|-->-5*x^2 + 45*x - 165/2 on (4, 5], x|-->-2*x^2 + 15*x - 15/2 on (5, 6], x|-->2*x^2 - 33*x + 273/2 on (6, 8], x|-->1/2*x^2 - 9*x + 81/2 on (8, 9]; x) sage: f = piecewise([[(-1,1),1]]) ## example 2 sage: g = piecewise([[(0,3),x]]) sage: f.convolution(g) piecewise(x|-->1/2*x^2 + x + 1/2 on (-1, 1], x|-->2*x on (1, 2], x|-->-1/2*x^2 + x + 4 on (2, 4]; x) sage: g = piecewise([[(0,3),1],[(3,4),2]]) sage: f.convolution(g) piecewise(x|-->x + 1 on (-1, 1], x|-->2 on (1, 2], x|-->x on (2, 3], x|-->-x + 6 on (3, 4], x|-->-2*x + 10 on (4, 5]; x)
Check that the bugs raised in trac ticket #12123 are fixed:
sage: f = piecewise([[(-2, 2), 2]]) sage: g = piecewise([[(0, 2), 3/4]]) sage: f.convolution(g) piecewise(x|-->3/2*x + 3 on (-2, 0], x|-->3 on (0, 2], x|-->-3/2*x + 6 on (2, 4]; x) sage: f = piecewise([[(-1, 1), 1]]) sage: g = piecewise([[(0, 1), x], [(1, 2), -x + 2]]) sage: f.convolution(g) piecewise(x|-->1/2*x^2 + x + 1/2 on (-1, 0], x|-->-1/2*x^2 + x + 1/2 on (0, 2], x|-->1/2*x^2 - 3*x + 9/2 on (2, 3]; x)
-
critical_points
(self, parameters, variable)¶ Return the critical points of this piecewise function.
EXAMPLES:
sage: R.<x> = QQ[] sage: f1 = x^0 sage: f2 = 10*x - x^2 sage: f3 = 3*x^4 - 156*x^3 + 3036*x^2 - 26208*x sage: f = piecewise([[(0,3),f1],[(3,10),f2],[(10,20),f3]]) sage: expected = [5, 12, 13, 14] sage: all(abs(e-a) < 0.001 for e,a in zip(expected, f.critical_points())) True
TESTS:
Use variables other than x (trac ticket #13836):
sage: R.<y> = QQ[] sage: f1 = y^0 sage: f2 = 10*y - y^2 sage: f3 = 3*y^4 - 156*y^3 + 3036*y^2 - 26208*y sage: f = piecewise([[(0,3),f1],[(3,10),f2],[(10,20),f3]]) sage: expected = [5, 12, 13, 14] sage: all(abs(e-a) < 0.001 for e,a in zip(expected, f.critical_points())) True
-
domain
(self, parameters, variable)¶ Return the domain
OUTPUT:
The union of the domains of the individual pieces as a
RealSet
.EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: f.domain() [0, 2)
-
domains
(self, parameters, variable)¶ Return the individual domains
See also
expressions()
.OUTPUT:
The collection of domains of the component functions as a tuple of
RealSet
.EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: f.domains() ({0}, (0, 2))
-
end_points
(self, parameters, variable)¶ Return a list of all interval endpoints for this function.
EXAMPLES:
sage: f1(x) = 1 sage: f2(x) = 1-x sage: f3(x) = x^2-5 sage: f = piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3]]) sage: f.end_points() [0, 1, 2, 3] sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: f.end_points() [0, 2]
-
expression_at
(self, parameters, variable, point)¶ Return the expression defining the piecewise function at
value
INPUT:
point
– a real number.
OUTPUT:
The symbolic expression defining the function value at the given
point
.EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: f.expression_at(0) sin(x) sage: f.expression_at(1) cos(x) sage: f.expression_at(2) Traceback (most recent call last): ... ValueError: point is not in the domain
-
expressions
(self, parameters, variable)¶ Return the individual domains
See also
domains()
.OUTPUT:
The collection of expressions of the component functions.
EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: f.expressions() (sin(x), cos(x))
-
extension
(self, parameters, variable, extension, extension_domain=None)¶ Extend the function
INPUT:
extension
– a symbolic expressionextension_domain
– aRealSet
orNone
(default). The domain of the extension. By default, the entire complement of the current domain.
EXAMPLES:
sage: f = piecewise([((-1,1), x)]); f piecewise(x|-->x on (-1, 1); x) sage: f(3) Traceback (most recent call last): ... ValueError: point 3 is not in the domain sage: g = f.extension(0); g piecewise(x|-->x on (-1, 1), x|-->0 on (-oo, -1] + [1, +oo); x) sage: g(3) 0 sage: h = f.extension(1, RealSet.unbounded_above_closed(1)); h piecewise(x|-->x on (-1, 1), x|-->1 on [1, +oo); x) sage: h(3) 1
-
fourier_series_cosine_coefficient
(self, parameters, variable, n, L)¶ Returns the n-th Fourier series coefficient of \(\cos(n\pi x/L)\), \(a_n\).
INPUT:
self
- the function f(x), defined over -L x Ln
- an integer n=0L
- (the period)/2
OUTPUT: \(a_n = \frac{1}{L}\int_{-L}^L f(x)\cos(n\pi x/L)dx\)
EXAMPLES:
sage: f(x) = x^2 sage: f = piecewise([[(-1,1),f]]) sage: f.fourier_series_cosine_coefficient(2,1) pi^(-2) sage: f(x) = x^2 sage: f = piecewise([[(-pi,pi),f]]) sage: f.fourier_series_cosine_coefficient(2,pi) 1 sage: f1(x) = -1 sage: f2(x) = 2 sage: f = piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_cosine_coefficient(5,pi) -3/5/pi
-
fourier_series_partial_sum
(self, parameters, variable, N, L)¶ Returns the partial sum
\[f(x) \sim \frac{a_0}{2} + \sum_{n=1}^N [a_n\cos(\frac{n\pi x}{L}) + b_n\sin(\frac{n\pi x}{L})],\]as a string.
EXAMPLE:
sage: f(x) = x^2 sage: f = piecewise([[(-1,1),f]]) sage: f.fourier_series_partial_sum(3,1) cos(2*pi*x)/pi^2 - 4*cos(pi*x)/pi^2 + 1/3 sage: f1(x) = -1 sage: f2(x) = 2 sage: f = piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_partial_sum(3,pi) -3*cos(x)/pi - 3*sin(2*x)/pi + 3*sin(x)/pi - 1/4
-
fourier_series_sine_coefficient
(self, parameters, variable, n, L)¶ Returns the n-th Fourier series coefficient of \(\sin(n\pi x/L)\), \(b_n\).
INPUT:
self
- the function f(x), defined over -L x Ln
- an integer n0L
- (the period)/2
OUTPUT: \(b_n = \frac{1}{L}\int_{-L}^L f(x)\sin(n\pi x/L)dx\)
EXAMPLES:
sage: f(x) = x^2 sage: f = piecewise([[(-1,1),f]]) sage: f.fourier_series_sine_coefficient(2,1) # L=1, n=2 0
-
integral
(self, parameters, variable, x=None, a=None, b=None, definite=False)¶ By default, return the indefinite integral of the function. If definite=True is given, returns the definite integral.
AUTHOR:
- Paul Butler
EXAMPLES:
sage: f1(x) = 1-x sage: f = piecewise([((0,1),1), ((1,2),f1)]) sage: f.integral(definite=True) 1/2
sage: f1(x) = -1 sage: f2(x) = 2 sage: f = piecewise([((0,pi/2),f1), ((pi/2,pi),f2)]) sage: f.integral(definite=True) 1/2*pi sage: f1(x) = 2 sage: f2(x) = 3 - x sage: f = piecewise([[(-2, 0), f1], [(0, 3), f2]]) sage: f.integral() piecewise(x|-->2*x + 4 on (-2, 0), x|-->-1/2*x^2 + 3*x + 4 on (0, 3); x) sage: f1(y) = -1 sage: f2(y) = y + 3 sage: f3(y) = -y - 1 sage: f4(y) = y^2 - 1 sage: f5(y) = 3 sage: f = piecewise([[[-4,-3],f1],[(-3,-2),f2],[[-2,0],f3],[(0,2),f4],[[2,3],f5]]) sage: F = f.integral(y) sage: F piecewise(y|-->-y - 4 on [-4, -3], y|-->1/2*y^2 + 3*y + 7/2 on (-3, -2), y|-->-1/2*y^2 - y - 1/2 on [-2, 0], y|-->1/3*y^3 - y - 1/2 on (0, 2), y|-->3*y - 35/6 on [2, 3]; y)
Ensure results are consistent with FTC:
sage: F(-3) - F(-4) -1 sage: F(-1) - F(-3) 1 sage: F(2) - F(0) 2/3 sage: f.integral(y, 0, 2) 2/3 sage: F(3) - F(-4) 19/6 sage: f.integral(y, -4, 3) 19/6 sage: f.integral(definite=True) 19/6
sage: f1(y) = (y+3)^2 sage: f2(y) = y+3 sage: f3(y) = 3 sage: f = piecewise([[(-infinity, -3), f1], [(-3, 0), f2], [(0, infinity), f3]]) sage: f.integral() piecewise(y|-->1/3*y^3 + 3*y^2 + 9*y + 9 on (-oo, -3), y|-->1/2*y^2 + 3*y + 9/2 on (-3, 0), y|-->3*y + 9/2 on (0, +oo); y)
sage: f1(x) = e^(-abs(x)) sage: f = piecewise([[(-infinity, infinity), f1]]) sage: f.integral(definite=True) 2 sage: f.integral() piecewise(x|-->-1/2*((sgn(x) - 1)*e^(2*x) - 2*e^x*sgn(x) + sgn(x) + 1)*e^(-x) - 1 on (-oo, +oo); x)
sage: f = piecewise([((0, 5), cos(x))]) sage: f.integral() piecewise(x|-->sin(x) on (0, 5); x)
TESTS:
Verify that piecewise integrals of zero work (trac ticket #10841):
sage: f0(x) = 0 sage: f = piecewise([[[0,1],f0]]) sage: f.integral(x,0,1) 0 sage: f = piecewise([[[0,1], 0]]) sage: f.integral(x,0,1) 0 sage: f = piecewise([[[0,1], SR(0)]]) sage: f.integral(x,0,1) 0
-
items
(self, parameters, variable)¶ Iterate over the pieces of the piecewise function
Note
You should probably use
pieces()
instead, which offers a nicer interface.OUTPUT:
This method iterates over pieces of the piecewise function, each represented by a pair. The first element is the support, and the second the function over that support.
EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]) sage: for support, function in f.items(): ....: print('support is {0}, function is {1}'.format(support, function)) support is {0}, function is sin(x) support is (0, 2), function is cos(x)
-
laplace
(self, parameters, variable, x='x', s='t')¶ Returns the Laplace transform of self with respect to the variable var.
INPUT:
x
- variable of selfs
- variable of Laplace transform.
We assume that a piecewise function is 0 outside of its domain and that the left-most endpoint of the domain is 0.
EXAMPLES:
sage: x, s, w = var('x, s, w') sage: f = piecewise([[(0,1),1],[[1,2], 1-x]]) sage: f.laplace(x, s) -e^(-s)/s + (s + 1)*e^(-2*s)/s^2 + 1/s - e^(-s)/s^2 sage: f.laplace(x, w) -e^(-w)/w + (w + 1)*e^(-2*w)/w^2 + 1/w - e^(-w)/w^2
sage: y, t = var('y, t') sage: f = piecewise([[[1,2], 1-y]]) sage: f.laplace(y, t) (t + 1)*e^(-2*t)/t^2 - e^(-t)/t^2
sage: s = var('s') sage: t = var('t') sage: f1(t) = -t sage: f2(t) = 2 sage: f = piecewise([[[0,1],f1],[(1,infinity),f2]]) sage: f.laplace(t,s) (s + 1)*e^(-s)/s^2 + 2*e^(-s)/s - 1/s^2
-
pieces
(self, parameters, variable)¶ Return the “pieces”.
OUTPUT:
A tuple of piecewise functions, each having only a single expression.
EXAMPLES:
sage: p = piecewise([((-1, 0), -x), ([0, 1], x)], var=x) sage: p.pieces() (piecewise(x|-->-x on (-1, 0); x), piecewise(x|-->x on [0, 1]; x))
-
piecewise_add
(self, parameters, variable, other)¶ Return a new piecewise function with domain the union of the original domains and functions summed. Undefined intervals in the union domain get function value \(0\).
EXAMPLES:
sage: f = piecewise([([0,1], 1), ((2,3), x)]) sage: g = piecewise([((1/2, 2), x)]) sage: f.piecewise_add(g).unextend_zero() piecewise(x|-->1 on (0, 1/2], x|-->x + 1 on (1/2, 1], x|-->x on (1, 2) + (2, 3); x)
-
restriction
(self, parameters, variable, restricted_domain)¶ Restrict the domain
INPUT:
restricted_domain
– aRealSet
or something that defines one.
OUTPUT:
A new piecewise function obtained by restricting the domain.
EXAMPLES:
sage: f = piecewise([((-oo, oo), x)]); f piecewise(x|-->x on (-oo, +oo); x) sage: f.restriction([[-1,1], [3,3]]) piecewise(x|-->x on [-1, 1] + {3}; x)
-
trapezoid
(self, parameters, variable, N)¶ Return the piecewise line function defined by the trapezoid rule for numerical integration based on a subdivision of each domain interval into N subintervals.
EXAMPLES:
sage: f = piecewise([[[0,1], x^2], [RealSet.open_closed(1,2), 5-x^2]]) sage: f.trapezoid(2) piecewise(x|-->1/2*x on (0, 1/2), x|-->3/2*x - 1/2 on (1/2, 1), x|-->7/2*x - 5/2 on (1, 3/2), x|-->-7/2*x + 8 on (3/2, 2); x) sage: f = piecewise([[[-1,1], 1-x^2]]) sage: f.trapezoid(4).integral(definite=True) 5/4 sage: f = piecewise([[[-1,1], 1/2+x-x^3]]) ## example 3 sage: f.trapezoid(6).integral(definite=True) 1
TESTS:
Use variables or rings other than x (trac ticket #13836):
sage: R.<y> = QQ[] sage: f1 = y^2 sage: f2 = 5-y^2 sage: f = piecewise([[[0,1],f1], [RealSet.open_closed(1,2),f2]]) sage: f.trapezoid(2) piecewise(y|-->1/2*y on (0, 1/2), y|-->3/2*y - 1/2 on (1/2, 1), y|-->7/2*y - 5/2 on (1, 3/2), y|-->-7/2*y + 8 on (3/2, 2); y)
-
unextend_zero
(self, parameters, variable)¶ Remove zero pieces.
EXAMPLES:
sage: f = piecewise([((-1,1), x)]); f piecewise(x|-->x on (-1, 1); x) sage: g = f.extension(0); g piecewise(x|-->x on (-1, 1), x|-->0 on (-oo, -1] + [1, +oo); x) sage: g(3) 0 sage: h = g.unextend_zero() sage: bool(h == f) True
-
which_function
(self, parameters, variable, point)¶ Return the expression defining the piecewise function at
value
INPUT:
point
– a real number.
OUTPUT:
The symbolic expression defining the function value at the given
point
.EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: f.expression_at(0) sin(x) sage: f.expression_at(1) cos(x) sage: f.expression_at(2) Traceback (most recent call last): ... ValueError: point is not in the domain
-
-
static
PiecewiseFunction.
in_operands
(ex)¶ Return whether a symbolic expression contains a piecewise function as operand
INPUT:
ex
– a symbolic expression.
OUTPUT:
Boolean
EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]); f piecewise(x|-->sin(x) on {0}, x|-->cos(x) on (0, 2); x) sage: piecewise.in_operands(f) True sage: piecewise.in_operands(1+sin(f)) True sage: piecewise.in_operands(1+sin(0*f)) False
-
static
PiecewiseFunction.
simplify
(ex)¶ Combine piecewise operands into single piecewise function
OUTPUT:
A piecewise function whose operands are not piecewiese if possible, that is, as long as the piecewise variable is the same.
EXAMPLES:
sage: f = piecewise([([0,0], sin(x)), ((0,2), cos(x))]) sage: piecewise.simplify(f) Traceback (most recent call last): ... NotImplementedError
-
class