Module: sage.functions.piecewise
Piecewise-defined Functions.
Sage implements a very simple class of piecewise-defined functions. Functions may be any type of symbolic expression. Infinite intervals are not supported. The endpoints of each interval must line up.
Implemented methods: latex outout __call__ __eq__ plot (using matplotlib) integral convolution derivative critical points tangent_line trapezoid trapezoid_integral_approximation riemann_sum (right and left end points) riemann_sum_integral_approximation (right and left end points) fourier series fourier series value coefficients (also sine series and cosine series, and Cesaro sum) partial sum (in string format) plot of partial sum plot_fourier_series_partial_sum_cesaro plot_fourier_series_partial_sum_hann plot_fourier_series_partial_sum_filter laplace transform latex output option domain range list __add__ - addition (of functions) __mul__ - multiplication (of functions, or fcn*scalar - ie, *right* multiplication by QQ) extend_by_zero_to unextend
TODO: [] Implement (a) max/min location and values, [] (b) multiplication by a scalar. [] (c) Extend the implementation of the trick to pass Sage's pi back and forth with Maxima's [[Passing the constants to maxima is already implemented; maybe need passing them back?]] [] Need: parent object - ring of piecewise functions [] This class should derive from an element-type class, and should define _add_, _mul_, etc. That will automatically take care of left multiplication and proper coercion. The coercion mentioned below for scalar mult on right is bad, since it only allows ints and rationals. The right way is to use an element class and only define _mul_, and have a parent, so anything gets coerced properly.
Author Log:
Module-level Functions
list_of_pairs) |
Returns a piecewise function from a list of (interval, function) pairs.
list_of_pairs
is a list of pairs (I, fcn), where fcn is
a SAGE function (such as a polynomial over RR, or functions
using the lambda notation), and I is an interval such as I = (1,3).
Two consecutive intervals must share a common endpoint.
We assume that these definitions are consistent (ie, no checking is done).
sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f(1) -1 sage: f(3) 2
x) |
list_of_pairs) |
Returns a piecewise function from a list of (interval, function) pairs.
list_of_pairs
is a list of pairs (I, fcn), where fcn is
a SAGE function (such as a polynomial over RR, or functions
using the lambda notation), and I is an interval such as I = (1,3).
Two consecutive intervals must share a common endpoint.
We assume that these definitions are consistent (ie, no checking is done).
sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f(1) -1 sage: f(3) 2
Class: PiecewisePolynomial
sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f(1) -1 sage: f(3) 2
self, list_of_pairs) |
list_of_pairs
is a list of pairs (I, fcn), where fcn is
a SAGE function (such as a polynomial over RR, or functions
using the lambda notation), and I is an interval such as I = (1,3).
Two consecutive intervals must share a common endpoint.
We assume that these definitions are consistent (ie, no checking is done).
Functions: base_ring,
convolution,
cosine_series_coefficient,
critical_points,
derivative,
domain,
end_points,
extend_by_zero_to,
fourier_series_cosine_coefficient,
fourier_series_partial_sum,
fourier_series_partial_sum_cesaro,
fourier_series_partial_sum_filtered,
fourier_series_partial_sum_hann,
fourier_series_sine_coefficient,
fourier_series_value,
functions,
integral,
intervals,
laplace,
length,
list,
plot,
plot_fourier_series_partial_sum,
plot_fourier_series_partial_sum_cesaro,
plot_fourier_series_partial_sum_filtered,
plot_fourier_series_partial_sum_hann,
riemann_sum,
riemann_sum_integral_approximation,
sine_series_coefficient,
tangent_line,
trapezoid,
trapezoid_integral_approximation,
unextend,
which_function
self) |
Returns the base-ring (ie, QQ[x]) - useful when this class is extended.
self, other) |
Returns the convolution function,
,
for compactly supported
.
sage: x = PolynomialRing(QQ,'x').gen() sage: f = Piecewise([[(0,1),1*x^0]]) ## example 0 sage: g = f.convolution(f) sage: h = f.convolution(g) sage: P = f.plot(); Q = g.plot(rgbcolor=(1,1,0)); R = h.plot(rgbcolor=(0,1,1)); sage: # Type show(P+Q+R) to view sage: f = Piecewise([[(0,1),1*x^0],[(1,2),2*x^0],[(2,3),1*x^0]]) ## example 1 sage: g = f.convolution(f) sage: h = f.convolution(g) sage: P = f.plot(); Q = g.plot(rgbcolor=(1,1,0)); R = h.plot(rgbcolor=(0,1,1)); sage: # Type show(P+Q+R) to view sage: f = Piecewise([[(-1,1),1]]) ## example 2 sage: g = Piecewise([[(0,3),x]]) sage: f.convolution(g) Piecewise defined function with 3 parts, [[(-1, 1), 0], [(1, 2), -3/2*x], [(2, 4), -3/2*x]] sage: g = Piecewise([[(0,3),1*x^0],[(3,4),2*x^0]]) sage: f.convolution(g) Piecewise defined function with 5 parts, [[(-1, 1), x + 1], [(1, 2), 3], [(2, 3), x], [(3, 4), -x + 8], [(4, 5), -2*x + 10]]
self, n, L) |
Returns the n-th cosine series coefficient of
,
.
Input:
sage: f = lambda x:x sage: f = Piecewise([[(0,1),f]]) sage: f.cosine_series_coefficient(2,1) 0 sage: f.cosine_series_coefficient(3,1) -4/(9*pi^2) sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f.cosine_series_coefficient(2,pi) 0 sage: f.cosine_series_coefficient(3,pi) 2/pi sage: f.cosine_series_coefficient(111,pi) 2/(37*pi)
self) |
Return the critical points of this piecewise function.
WARNINGS: Uses maxima, which prints the warning to use results with caution. Only works for piecewise functions whose parts are polynomials with real critical not occurring on the interval endpoints.
sage: x = PolynomialRing(QQ, 'x').0 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: f.critical_points() [5.0, 12.000000000000171, 12.9999999999996, 14.000000000000229]
self) |
Returns the derivative (as computed by maxima)
Piecewise(I,
), as I runs over the intervals
belonging to self. self must be piecewise polynomial.
sage: f1 = lambda x:1 sage: f2 = lambda x:1-x sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: f.derivative() Piecewise defined function with 2 parts, [[(0, 1), 0], [(1, 2), -1]] sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f.derivative() Piecewise defined function with 2 parts, [[(0, pi/2), 0], [(pi/2, pi), 0]]
self) |
Returns the domain of the function.
sage: f1(x) = 1 sage: f2(x) = 1-x sage: f3(x) = exp(x) sage: f4(x) = sin(2*x) sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: f.domain() (0, 10)
self, [xmin=-1000], [xmax=1000]) |
This function simply returns the piecewise defined function which is extended by 0 so it is defined on all of (xmin,xmax). This is needed to add two piecewise functions in a reasonable way.
self, n, L) |
Returns the n-th Fourier series coefficient of
,
.
Input:
sage: f = lambda x:x^2 sage: f = Piecewise([[(-1,1),f]]) sage: f.fourier_series_cosine_coefficient(2,1) 1/pi^2 sage: f = lambda x:x^2 sage: f = Piecewise([[(-pi,pi),f]]) sage: float(f.fourier_series_cosine_coefficient(2,pi)) 1.0 sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_cosine_coefficient(5,pi) -3/(5*pi)
self, N, L) |
Returns the partial sum
as a string.
sage: f = lambda 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 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_partial_sum(3,pi) -3*sin(2*x)/pi + 3*sin(x)/pi - 3*cos(x)/pi - 1/4
self, N, L) |
Returns the Cesaro partial sum
as a string. This is a "smoother" partial sum - the Gibbs phenomenon is mollified.
sage: f = lambda x:x^2 sage: f = Piecewise([[(-1,1),f]]) sage: f.fourier_series_partial_sum_cesaro(3,1) cos(2*pi*x)/(3*pi^2) - 8*cos(pi*x)/(3*pi^2) + 1/3 sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_partial_sum_cesaro(3,pi) -sin(2*x)/pi + 2*sin(x)/pi - 2*cos(x)/pi - 1/4
self, N, L, F) |
Returns the "filtered" partial sum
as a string, where
sage: f = lambda x:x^2 sage: f = Piecewise([[(-1,1),f]]) sage: f.fourier_series_partial_sum_filtered(3,1,[1,1,1]) cos(2*pi*x)/pi^2 - 4*cos(pi*x)/pi^2 + 1/3 sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_partial_sum_filtered(3,pi,[1,1,1]) -3*sin(2*x)/pi + 3*sin(x)/pi - 3*cos(x)/pi - 1/4
self, N, L) |
Returns the Hann-filtered partial sum (named after von Hann, not Hamming)
as a string, where
sage: f = lambda x:x^2 sage: f = Piecewise([[(-1,1),f]]) sage: f.fourier_series_partial_sum_hann(3,1) 0.500000000000000*(cos(2*pi/3) + 1)*cos(2*pi*x)/pi^2 - 2.00000000000000*(cos(pi/3) + 1)*cos(pi*x)/pi^2 + 1/3 sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_partial_sum_hann(3,pi) -1.50000000000000*(cos(2*pi/3) + 1)*sin(2*x)/pi + 1.50000000000000*(cos(pi/3) + 1)*sin(x)/pi - 1.50000000000000*(cos(pi/3) + 1)*cos(x)/pi - 1/4
self, n, L) |
Returns the n-th Fourier series coefficient of
,
.
Input:
sage: f = lambda x:x^2 sage: f = Piecewise([[(-1,1),f]]) sage: f.fourier_series_sine_coefficient(2,1) # L=1, n=2 0
self, x, L) |
Returns the value of the Fourier series coefficient of self at
,
This method applies to piecewise non-polynomial functions as well.
Input:
sage: f1 = lambda x:1 sage: f2 = lambda x:1-x sage: f3 = lambda x:exp(x) sage: f4 = lambda x:sin(2*x) sage: f = Piecewise([[(-10,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: f.fourier_series_value(101,10) 1/2 sage: f.fourier_series_value(100,10) 1 sage: f.fourier_series_value(10,10) sin(20) sage: f.fourier_series_value(20,10) 1 sage: f.fourier_series_value(30,10) sin(20) sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,0),lambda x:0],[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f.fourier_series_value(-1,pi) 0 sage: f.fourier_series_value(20,pi) -1 sage: f.fourier_series_value(pi/2,pi) 1/2
self) |
Returns the list of functions (the "pieces").
sage: f1(x) = 1 sage: f2(x) = 1-x sage: f3(x) = exp(x) sage: f4(x) = sin(2*x) sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: f.functions() [x |--> 1, x |--> 1 - x, x |--> e^x, x |--> sin(2*x)]
self, [x=None]) |
Returns the definite integral (as computed by maxima)
, as I runs over the intervals
belonging to self.
sage: f1(x) = 1 sage: f2(x) = 1-x sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: f.integral() 1/2 sage: f1(x) = -1 sage: f2(x) = 2 sage: f = Piecewise([[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: f.integral() pi/2
self) |
A piecewise non-polynomial example.
sage: f1(x) = 1 sage: f2(x) = 1-x sage: f3(x) = exp(x) sage: f4(x) = sin(2*x) sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: f.intervals() [(0, 1), (1, 2), (2, 3), (3, 10)]
self, x, s) |
Returns the Laplace transform of self with respect to the variable var.
Input:
We assume that a piecewise function is 0 outside of its domain and that the left-most endpoint of the domain is 0.
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 - e^(-s)/s^2 + (s + 1)*e^(-(2*s))/s^2 + 1/s sage: f.laplace(x, w) -e^(-w)/w - e^(-w)/w^2 + (w + 1)*e^(-(2*w))/w^2 + 1/w
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
self) |
Returns the plot of self.
Keyword arguments are passed onto the plot command for each piece of the function. E.g., the plot_points keyword affects each segment of the plot.
sage: f1 = lambda x:1 sage: f2 = lambda x:1-x sage: f3 = lambda x:exp(x) sage: f4 = lambda x:sin(2*x) sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: P = f.plot(rgbcolor=(0.7,0.1,0), plot_points=40)
Remember: to view this, type show(P) or P.save("<path>/myplot.png") and then open it in a graphics viewer such as GIMP.
self, N, L, xmin, xmax) |
Plots the partial sum
over xmin < x < xmin.
sage: f1 = lambda x:-2 sage: f2 = lambda x:1 sage: f3 = lambda x:-1 sage: f4 = lambda x:2 sage: f = Piecewise([[(-pi,-pi/2),f1],[(-pi/2,0),f2],[(0,pi/2),f3],[(pi/2,pi),f4]]) sage: P = f.plot_fourier_series_partial_sum(3,pi,-5,5) # long time sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: P = f.plot_fourier_series_partial_sum(15,pi,-5,5) # long time
Remember, to view this type show(P) or P.save("<path>/myplot.png") and then open it in a graphics viewer such as GIMP.
self, N, L, xmin, xmax) |
Plots the partial sum
over xmin < x < xmin. This is a "smoother" partial sum - the Gibbs phenomenon is mollified.
sage: f1 = lambda x:-2 sage: f2 = lambda x:1 sage: f3 = lambda x:-1 sage: f4 = lambda x:2 sage: f = Piecewise([[(-pi,-pi/2),f1],[(-pi/2,0),f2],[(0,pi/2),f3],[(pi/2,pi),f4]]) sage: P = f.plot_fourier_series_partial_sum_cesaro(3,pi,-5,5) # long time sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: P = f.plot_fourier_series_partial_sum_cesaro(15,pi,-5,5) # long time
Remember, to view this type show(P) or P.save("<path>/myplot.png") and then open it in a graphics viewer such as GIMP.
self, N, L, F, xmin, xmax) |
Plots the partial sum
over xmin < x < xmin, where
sage: f1 = lambda x:-2 sage: f2 = lambda x:1 sage: f3 = lambda x:-1 sage: f4 = lambda x:2 sage: f = Piecewise([[(-pi,-pi/2),f1],[(-pi/2,0),f2],[(0,pi/2),f3],[(pi/2,pi),f4]]) sage: P = f.plot_fourier_series_partial_sum_filtered(3,pi,[1]*3,-5,5) # long time sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,-pi/2),f1],[(-pi/2,0),f2],[(0,pi/2),f1],[(pi/2,pi),f2]]) sage: P = f.plot_fourier_series_partial_sum_filtered(15,pi,[1]*15,-5,5) # long time
Remember, to view this type show(P) or P.save("<path>/myplot.png") and then open it in a graphics viewer such as GIMP.
self, N, L, xmin, xmax) |
Plots the partial sum
over xmin < x < xmin, where H_N(x) = (0.5)+(0.5)*cos(x*pi/N) is the N-th Hann filter.
sage: f1 = lambda x:-2 sage: f2 = lambda x:1 sage: f3 = lambda x:-1 sage: f4 = lambda x:2 sage: f = Piecewise([[(-pi,-pi/2),f1],[(-pi/2,0),f2],[(0,pi/2),f3],[(pi/2,pi),f4]]) sage: P = f.plot_fourier_series_partial_sum_hann(3,pi,-5,5) # long time sage: f1 = lambda x:-1 sage: f2 = lambda x:2 sage: f = Piecewise([[(-pi,pi/2),f1],[(pi/2,pi),f2]]) sage: P = f.plot_fourier_series_partial_sum_hann(15,pi,-5,5) # long time
Remember, to view this type show(P) or P.save("<path>/myplot.png") and then open it in a graphics viewer such as GIMP.
self, N, [mode=None]) |
Returns the piecewise line function defined by the Riemann sums in numerical integration based on a subdivision into N subintervals. Set mode="midpoint" for the height of the rectangles to be determined by the midpoint of the subinterval; set mode="right" for the height of the rectangles to be determined by the right-hand endpoint of the subinterval; the default is mode="left" (the height of the rectangles to be determined by the leftt-hand endpoint of the subinterval).
sage: f1(x) = x^2 sage: f2(x) = 5-x^2 sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: f.riemann_sum(6,mode="midpoint") Piecewise defined function with 6 parts, [[(0, 1/3), 1/36], [(1/3, 2/3), 1/4], [(2/3, 1), 25/36], [(1, 4/3), 131/36], [(4/3, 5/3), 11/4], [(5/3, 2), 59/36]] sage: f = Piecewise([[(-1,1),1-x^2]]) sage: rsf = f.riemann_sum(7) sage: P = f.plot(rgbcolor=(0.7,0.1,0.5), plot_points=40) sage: Q = rsf.plot(rgbcolor=(0.7,0.6,0.6), plot_points=40) sage: L = add([line([[pf[0][0],0],[pf[0][0],pf[1](x=pf[0][0])]],rgbcolor=(0.7,0.6,0.6)) for pf in rsf.list()]) sage: ## To view this, type show(P+Q+L). sage: f = Piecewise([[(-1,1),1/2+x-x^3]]) ## example 3 sage: rsf = f.riemann_sum(8) sage: P = f.plot(rgbcolor=(0.7,0.1,0.5), plot_points=40) sage: Q = rsf.plot(rgbcolor=(0.7,0.6,0.6), plot_points=40) sage: L = add([line([[pf[0][0],0],[pf[0][0],pf[1](x=pf[0][0])]],rgbcolor=(0.7,0.6,0.6)) for pf in rsf.list()]) sage: ## To view this, type show(P+Q+L).
self, N, [mode=None]) |
Returns the piecewise line function defined by the Riemann sums in numerical integration based on a subdivision into N subintervals. Set mode="midpoint" for the height of the rectangles to be determined by the midpoint of the subinterval; set mode="right" for the height of the rectangles to be determined by the right-hand endpoint of the subinterval; the default is mode="left" (the height of the rectangles to be determined by the leftt-hand endpoint of the subinterval).
sage: f1 = x^2 ## example 1 sage: f2 = 5-x^2 sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: f.riemann_sum_integral_approximation(6) 17/6 sage: f.riemann_sum_integral_approximation(6,mode="right") 19/6 sage: f.riemann_sum_integral_approximation(6,mode="midpoint") 3 sage: f.integral() 3
self, n, L) |
Returns the n-th sine series coefficient of
,
.
Input:
sage: f = lambda x:1 sage: f = Piecewise([[(0,1),f]]) sage: f.sine_series_coefficient(2,1) 0 sage: f.sine_series_coefficient(3,1) 4/(3*pi)
self, pt) |
Computes the linear function defining the tangent line of the piecewise function self.
sage: f1 = lambda x:x^2 sage: f2 = lambda x:5-x^3+x sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: tf = f.tangent_line(0.9) ## tangent line at x=0.9 sage: P = f.plot(rgbcolor=(0.7,0.1,0.5), plot_points=40) sage: Q = tf.plot(rgbcolor=(0.7,0.2,0.2), plot_points=40)
Type show(P+Q) to view the graph of the function and the tangent line.
self, N) |
Returns the piecewise line function defined by the trapezoid rule for numerical integration based on a subdivision into N subintervals.
sage: f1 = lambda x:x^2 ## example 1 sage: f2 = lambda x:5-x^2 sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: f.trapezoid(4) Piecewise defined function with 4 parts, [[(0, 1/2), 1/2*x], [(1/2, 1), 9/2*x - 2], [(1, 3/2), 1/2*x + 2], [(3/2, 2), -7/2*x + 8]] sage: x = PolynomialRing(QQ,'x').gen() ## example 2 sage: f = Piecewise([[(-1,1),1-x^2]]) sage: tf = f.trapezoid(4) sage: P = f.plot(rgbcolor=(0.7,0.1,0.5), plot_points=40) sage: Q = tf.plot(rgbcolor=(0.7,0.6,0.6), plot_points=40) sage: L = add([line([[pf[0][0],0],[pf[0][0],pf[1](pf[0][0])]],rgbcolor=(0.7,0.6,0.6)) for pf in tf.list()]) sage: ## To view this, type show(P+Q+L). sage: f = Piecewise([[(-1,1),1/2+x-x^3]]) ## example 3 sage: tf = f.trapezoid(6) sage: P = f.plot(rgbcolor=(0.7,0.1,0.5), plot_points=40) sage: Q = tf.plot(rgbcolor=(0.7,0.6,0.6), plot_points=40) sage: L = add([line([[pf[0][0],0],[pf[0][0],pf[1](pf[0][0])]],rgbcolor=(0.7,0.6,0.6)) for pf in tf.list()]) sage: ## To view this, type show(P+Q+L).
self, N) |
Returns the approximation given by the trapezoid rule for numerical integration based on a subdivision into N subintervals.
sage: f1(x) = x^2 ## example 1 sage: f2(x) = 1-(1-x)^2 sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: P = f.plot(rgbcolor=(0.7,0.1,0.5), plot_points=40) sage: tf = f.trapezoid(6) sage: Q = tf.plot(rgbcolor=(0.7,0.6,0.6), plot_points=40) sage: ta = f.trapezoid_integral_approximation(6) sage: t = text('trapezoid approximation = %s'%ta, (1.5, 0.25)) sage: a = f.integral() sage: tt = text('area under curve = %s'%a, (1.5, -0.5)) sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) ## example 2 sage: tf = f.trapezoid(4) sage: ta = f.trapezoid_integral_approximation(4) sage: Q = tf.plot(rgbcolor=(0.7,0.6,0.6), plot_points=40) sage: t = text('trapezoid approximation = %s'%ta, (1.5, 0.25)) sage: a = f.integral() sage: tt = text('area under curve = %s'%a, (1.5, -0.5)) sage: ## To view this, type show(P+Q+L).
self) |
This removes any parts in the front or back of the function which is zero (the inverse to extend_by_zero_to).
sage: x = PolynomialRing(QQ,'x').gen() sage: f = Piecewise([[(-3,-1),1+2+x],[(-1,1),1-x^2]]) sage: e = f.extend_by_zero_to(-10,10); e Piecewise defined function with 4 parts, [[(-10, -3), 0], [(-3, -1), x + 3], [(-1, 1), -x^2 + 1], [(1, 10), 0]] sage: d = e.unextend(); d Piecewise defined function with 2 parts, [[(-3, -1), x + 3], [(-1, 1), -x^2 + 1]] sage: d==f True
self, x0) |
Returns the function piece used to evaluate self at x0.
sage: f1(z) = z sage: f2(x) = 1-x sage: f3(y) = exp(y) sage: f4(t) = sin(2*t) sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: f.which_function(3/2) x |--> 1 - x
Special Functions: __add__,
__call__,
__eq__,
__init__,
__mul__,
__repr__,
_latex_
self, other) |
Returns the piecewise defined function which is the sum of self and other. Doesnot require both domains be the same.
sage: x = PolynomialRing(QQ,'x').gen() sage: f1 = x^0 sage: f2 = 1-x sage: f3 = 2*x sage: f4 = 10-x sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: g1 = x-2 sage: g2 = x-5 sage: g = Piecewise([[(0,5),g1],[(5,10),g2]]) sage: h = f+g sage: h Piecewise defined function with 5 parts, [[(0, 1), x - 1], [(1, 2), -1], [(2, 3), 3*x - 2], [(3, 5), 8], [(5, 10), 5]]
Note that in this case the functions must be defined using polynomial expressions *not* using the lambda notation.
self, x0) |
Evaluates self at x0. Returns the average value of the jump if x0 is an interior endpoint of one of the intervals of self and the usual value otherwise.
sage: f1(x) = 1 sage: f2(x) = 1-x sage: f3(x) = exp(x) sage: f4(x) = sin(2*x) sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: f(0.5) 1 sage: f(2.5) 12.1824939607034... sage: f(1) 1/2
self, other) |
Implements Boolean == operator.
self, other) |
Returns the piecewise defined function which is the product of one piecewise function (self) with another one (other).
sage: x = PolynomialRing(QQ,'x').gen() sage: f1 = x^0 sage: f2 = 1-x sage: f3 = 2*x sage: f4 = 10-x sage: f = Piecewise([[(0,1),f1],[(1,2),f2],[(2,3),f3],[(3,10),f4]]) sage: g1 = x-2 sage: g2 = x-5 sage: g = Piecewise([[(0,5),g1],[(5,10),g2]]) sage: h = f*g sage: h Piecewise defined function with 5 parts, [[(0, 1), x - 2], [(1, 2), -x^2 + 3*x - 2], [(2, 3), 2*x^2 - 4*x], [(3, 5), -x^2 + 12*x - 20], [(5, 10), -x^2 + 15*x - 50]] sage: g*(11/2) Piecewise defined function with 2 parts, [[(0, 5), 11/2*x - 11], [(5, 10), 11/2*x - 55/2]]
Note that in this method the functions must be defined using polynomial expressions *not* using the lambda notation.
self) |
sage: f1(x) = 1 sage: f2(x) = 1 - x sage: f = Piecewise([[(0,1),f1],[(1,2),f2]]) sage: latex(f) \begin{cases} x \ {\mapsto}\ 1 \&\text{on $(0, 1)$}\cr x \ {\mapsto}\ 1 - x \&\text{on $(1, 2)$} \end{cases}
sage: f(x) = sin(x*pi/2) sage: g(x) = 1-(x-1)^2 sage: h(x) = -x sage: P = Piecewise([[(0,1), f], [(1,3),g], [(3,5), h]]) sage: latex(P) \begin{cases} x \ {\mapsto}\ \sin \left( \frac{{\pi x}}{2} \right) \&\text{on $(0, 1)$}\cr x \ {\mapsto}\ 1 - {\left( x - 1 \right)}^{2} \&\text{on $(1, 3)$}\cr x \ {\mapsto}\ -x \&\text{on $(3, 5)$} \end{cases}
See About this document... for information on suggesting changes.