Approximations to Piecewise Continuous Functions 

 

Univ.-Prof. Dr.-Ing. habil. Josef  BETTEN 

RWTH Aachen University 

Templergraben 55 

D-52056  A a c h e n ,  Germany 

betten@mmw.rwth-aachen.de 

 

 

Abstract 

 

This worksheet is concerned with approximations to piecewise continuous functions. It has been illustrated that Maplesoft furnishes powerful tools in finding best approximations to several functions. Some examples are discussed in more detail.  

 

Keywords:  step functions;  HEAVISIDE functions; FOURIER approximation; L-two norm  

 

                   

FOURIER Approximation 

 

 

>

restart:

 

 

>	FOURIER_series(x):=          a[0]/2+sum(a[k]*cos(k*x)+b[k]*sin(k*x),k=1..infinity);

 

 

(1)

 

 

 

>	a[k]:=(1/Pi)*Int(f(x)*cos(k*x),x=0..2*Pi);

 

 

(2)

 

 

>	a[0]:=simplify(subs(k=0,%));

 

 

(3)

 

 

>	b[k]:=(1/Pi)*Int(f(x)*sin(k*x),x=0..2*Pi);

 

 

(4)

 

 

>

 

 

Unit Step Function or HEAVISIDE Unit Function 

 

 

>	F(x):=piecewise(x>0 and x<Pi,1,x>Pi and x<2*Pi,-1,x>2*Pi,1);

 

 

(5)

 

 

>	plot(%,x=0..2.1*Pi,labels=[x,F],color=black);

 

 

 

 

Alternatively, the piecewise continuous function F(x) can be represented as a HEAVISIDE  

function: 

 

>	alias(H=Heaviside):

 

 

>	HEAVISIDE[F]:=convert(F(x),H);

 

 

(6)

 

 

>	plot(%,x=0..2.1*Pi,color=black,                                  title="HEAVISIDE unit step function");

 

 

 

 

>	A[0]:=value(subs(f(x)=F(x),a[0]));

 

 

(7)

 

 

>	A[k]:=simplify(value(subs(f(x)=F(x),a[k])));

 

 

(8)

 

 

>	A[k]:=subs(sin(k*Pi)=0,%);

 

 

(9)

 

 

>	B[k]:=simplify(value(subs(f(x)=F(x),b[k])));

 

 

(10)

 

 

>	B[k]:=simplify(subs(cos(k*Pi)=(-1)^k,%));

 

 

(11)

 

 

>	for i in [1,3,5] do y(x,n=i):=sum(B[k]*sin(k*x),k=1..i) od;

 

 

(12)

 

 

(12)

 

 

(12)

 

 

>

 

 

>	plot({y(x,n=1),y(x,n=3),y(x,n=5)},x=0..3.5*Pi,          labels=[x,y],color=black,                                    title="Approximations with n = [1, 3, 5]");

 

 

 

 

>	for i in [1,3,99] do y(x,n=i):=sum(B[k]*sin(k*x),k=1..i) od:

 

 

>	plot({y(x,n=1),y(x,n=3),y(x,n=99)},x=0..3.5*Pi,         labels=[x,y],color=black,                                   title="Approximations with n = [1, 3, 99]");

 

 

 

 

>	L_zwei[n]:=sqrt((1/Pi)*Int((F(x_)-y(x,n))^2,x=0..Pi));

 

 

(13)

 

 

>	for i in [1,3,5,99] do     L_zwei[n=i]:=            evalf(sqrt((1/Pi)*int((F(x)-y(x,n=i))^2,x=0.05..0.95*Pi)),4) od;

 

 

(14)

 

 

(14)

 

 

(14)

 

 

(14)

 

 

 

>

 

Example:  h(x) 

 

 

>	h(x):=piecewise(x>0 and x<Pi,Pi,x>Pi and x<2*Pi,x-Pi,                      x>2*Pi and x<3*Pi,Pi);

 

 

(15)

 

 

 

 

 

>	plot(%,x=0..3*Pi,color=black,labels=[x,h]);

 

 

 

 

>	alias(H=Heaviside):

 

 

>	HEAVISIDE[h]:=convert(h(x),H);

 

 

(16)

 

 

>	plot(%,x=0..3*Pi,color=black,                                       title="HEAVISIDE  Function for  h(x)");

 

 

 

 

>	A[0]:=value(subs(f(x)=h(x),a[0]));

 

 

(17)

 

 

>	A[k]:=simplify(value(subs(f(x)=h(x),a[k]))):

 

 

>	A[k]:=subs({sin(k*Pi)=0,cos(k*Pi)=(-1)^k,(cos(k*Pi))^2=1},%);

 

 

(18)

 

 

>	B[k]:=simplify(value(subs(f(x)=h(x),b[k]))):

 

 

>	B[k]:=subs({sin(k*Pi)=0,cos(k*Pi)=(-1)^k,(cos(k*Pi))^2=1},%);

 

 

(19)

 

 

>	y(x,n):=A[0]/2+sum(A[k]*cos(k*x)+B[k]*sin(k*x),k=1..n);

 

 

(20)

 

 

>	for i in [1,2,3] do                           y(x,n=i):=A[0]/2+sum(A[k]*cos(k*x)+B[k]*sin(k*x),k=1..i)                      od;

 

 

(21)

 

 

(21)

 

 

(21)

 

 

>

 

 

>	alias(H=Heaviside,co=color):

 

 

>	p[1]:=plot({y(x,n=1),y(x,n=2),y(x,n=3)},x=0..3*Pi,co=black):

 

 

>	p[2]:=plot(1.72*H(x-3*Pi),x=0..3.01*Pi,co=black):

 

 

>	p[3]:=plot(h(x),x=0..3*Pi,                                  title="Approximations with n = [1,  2,  3]" ):

 

 

>	plots[display]({seq(p[k],k=1..3)});

 

 

 

 

>	for i in [1,2,50] do                          y(x,n=i):=A[0]/2+sum(A[k]*cos(k*x)+B[k]*sin(k*x),k=1..i)                      od:

 

 

>	p[4]:=plot({y(x,n=1),y(x,n=2),y(x,n=50)},x=0..3*Pi,co=black):

 

 

>	p[5]:=plot(1.72*H(x-3*Pi),x=0..3.01*Pi,co=black,            title="Approximations with n = [1, 2, 50]"):

 

 

>	plots[display]({p[4],p[5]});

 

 

 

Smoothing of  the FOURIER-Series  y(x, n=3): 

 

 

>	g(k,N):=sin(Pi*k/N)/(Pi*k/N);  # smoothing factor

 

 

(22)

 

 

>	G(x,n,N):=A[0]/2+  sum(g(kappa,Nu)*(A[kappa]*cos(kappa*x)+B[kappa]*sin(kappa*x)),       kappa=1..n);  # smoothing function

 

 

(23)

 

 

>	g(k,4):=subs(N=4,g(k,N));       # N := n + 1

 

 

(24)

 

 

>	G(x,n=3,N=4):=subs({n=3,Nu=4,kappa=k,g(kappa,Nu)=g(k,4),     A[kappa]=A[k],B[kappa]=B[k]},G(x,n,N));

 

 

(25)

 

 

>	alias(H=Heaviside,co=color):

 

 

>	p[6]:=plot({y(x,n=3),G(x,n=3,N=4)},x=0..3*Pi,co=black):

 

 

>	p[7]:=plot(1.72*H(x-3*Pi),x=0..3.01*Pi,co=black):

 

 

>	p[8]:=plot(h(x),x=0..3*Pi,                                   title="Approximation with n = 3 and its Smoothing"):

 

 

>	plots[display](seq(p[k],k=6..8));

 

 

 

Smoothing of  the FOURIER-Series  y(x, n = 2): 

 

g(k,3):=subs(N=3,g(k,N)); 

 

(26)

 

 

>	G(x,n=2,N=3):=evalf(A[0]/2+           sum(g(k,3)*(A[k]*cos(k*x)+B[k]*sin(k*x)),k=1..2),4);

 

 

(27)

 

 

>	alias(H=Heaviside,co=color):

 

 

>	p[9]:=plot({h(x),y(x,n=2),G(x,n=2,N=3)},x=0..3*Pi,co=black,       title="Approximation with n = 2 and its Smoothing"):

 

 

>	plots[display]({p[7],p[9]});

 

 

 

The above Figures illustrate the good approximations  to the given function h(x). 

 

Example:  K(x) 

 

 

>	K(x):=piecewise(x>-Pi and x<0,3,x>0 and x<1/2,0,                                            x>1/2 and x<Pi,1/x);

 

 

(28)

 

 

>

 

 

>	plot(%,x=-Pi..Pi,color=black,labels=[x,K]);

 

 

 

 

>	alpha[0]:=(1/Pi)*Int(f(x),x=-Pi..Pi);

 

 

(29)

 

 

>	alpha[k]:=(1/Pi)*Int(f(x)*cos(k*x),x=-Pi..Pi);

 

 

(30)

 

 

>	beta[k]:=(1/Pi)*Int(f(x)*sin(k*x),x=-Pi..Pi);

 

 

(31)

 

 

>	Alpha[0]:=value(subs(f(x)=K(x),alpha[0]));

 

 

(32)

 

 

>	Alpha[k]:=simplify(value(subs(f(x)=K(x),alpha[k])));

 

 

(33)

 

 

>	A[k]:=subs(sin(k*Pi)=0,%);

 

 

(34)

 

 

>	BETA[k]:=simplify(value(subs(f(x)=K(x),beta[k])));

 

 

(35)

 

 

>	BETA[k]:=subs(cos(k*Pi)=(-1)^k,%);

 

 

(36)

 

 

>	y(x,n):=Alpha[0]/2+               sum(Alpha[k]*cos(k*x)+BETA[k]*sin(k*x),k=1..n):

 

 

>	y(x,n):=subs(sin(k*Pi)=0,%);

 

 

(37)

 

 

>	y(x,1):=subs(n=1,y(x,n)):

 

 

>	Y(x,1):=evalf(%,4);

 

 

(38)

 

 

>	y(x,3):=subs(n=3,y(x,n)):

 

 

>	Y(x,3):=evalf(%,4);

 

 

(39)

 

 

>	y(x,149):=subs(n=149,y(x,n)):

 

 

>	Y(x,149):=evalf(%):

 

 

>	plot({Y(x,1),Y(x,3),Y(x,149)},x=-Pi..Pi,co=black,            title="Approximations with n = [1, 3, 149]");

 

 

 

 

>	L_two[n]:=sqrt((1/2/Pi)*Int((K(x_)-Y(x,n))^2,x=-Pi..Pi));

 

 

(40)

 

 

>	for i in [1,3,149] do                L_two[n=i]:=evalf(sqrt((1/2/Pi)*int((K(x)-Y(x,i))^2,x=-3..3)),4) od;

 

 

(41)

 

 

(41)

 

 

(41)

 

 

>

 

 

Example:  M(x) 

 

 

>

restart:

 

 

>	M(x):=piecewise(-Pi<x and x<0,1,                                              x>0 and x<1/2,2, x>1/2 and x<Pi,1/x);

 

  

 

(42)

 

 

>	plot(%,x=-Pi..Pi,color=black,labels=[x,M]);

 

 

 

 

>

 

 

>	a[0]:=(1/Pi)*Int(f(x),x=-Pi..Pi);

 

 

(43)

 

 

>	a[k]:=(1/Pi)*Int(f(x)*cos(k*x),x=-Pi..Pi);

 

 

(44)

 

 

>	b[k]:=(1/Pi)*Int(f(x)*sin(k*x),x=-Pi..Pi);

 

 

(45)

 

 

>	A[0]:=value(subs(f(x)=M(x),a[0]));

 

 

(46)

 

 

>	A[k]:=simplify(value(subs(f(x)=M(x),a[k])));

 

 

(47)

 

 

>	A[k]:=subs(sin(Pi*k)=0,%);

 

 

(48)

 

 

>	B[k]:=simplify(value(subs(f(x)=M(x),b[k])));

 

 

(49)

 

 

>	B[k]:=subs(cos(Pi*k)=(-1)^k,%);

 

 

(50)

 

 

>	for i in [1,2,3,10] do                                       y(x,n=i):=evalf(A[0]/2+                sum(A[k]*cos(k*x)+B[k]*sin(k*x),k=1..i),4) od;

 

 

(51)

 

 

(51)

 

 

(51)

 

 

(51)

 

 

>	alias(H=Heaviside,th=thickness,co=color): p[1]:=plot({seq(y(x,n=i),i=1..3)},x=-Pi..Pi,co=black):

 

 

>	p[2]:=plot(M(x),x=-Pi..Pi,th=2,co=black):

 

 

>	p[3]:=plot(0.7191*H(x+Pi),x=-1.001*Pi..-0.999*Pi,co=black):

 

 

>	p[4]:=plot(0.7196*H(x-Pi),x=0.999*Pi..1.001*Pi,co=black,    title="Approximations with  n = [1, 2, 3]"):

 

 

>	plots[display]({seq(p[k],k=1..4)});

 

 

 

 

>

 

 

>	y(x,n=99):=evalf(A[0]/2+             sum(A[k]*cos(k*x)+B[k]*sin(k*x),k=1..99),4):

 

 

>	plot({M(x),y(x,n=99)},x=-Pi..Pi,color=black,                          title="M(x)  and  y(x, n = 99)");

 

 

 

 

 

>

 

 

>	L_two[n]:=sqrt((1/2/Pi)*Int((M(x_)-y(x,n))^2,x=-Pi..Pi));

 

 

(52)

 

 

>	for i in [1,2,3,99] do     L_two[n=i]:=evalf(sqrt((1/2/Pi)*int((M(x)-y(x,n=i))^2,               x=-3..3)),4) od;

 

 

(53)

 

 

(53)

 

 

(53)

 

 

(53)

 

 

>

 

Smoothing of  the FOURIER-Series  y(x, n = 3): 

 

>	g(k,N):=N*sin(k*Pi/N)/k/Pi;

 

 

(54)

 

 

>	g(k,4):=subs(N=4,%);

 

 

(55)

 

 

>	G(x,n=3,N=4):=evalf(A[0]/2+  sum(g(k,4)*(A[k]*cos(k*x)+B[k]*sin(k*x)),k=1..3),4);

 

 

(56)

 

 

>	plot({M(x),y(x,n=3),G(x,n=3,N=4)},x=-Pi..Pi,color=black,          title="Smoothing of the FOURIER-Series  y(x, n = 3)");

 

 

 

 

>

 

 

>	g(k,11):=subs(N=11,g(k,N));

 

 

(57)

 

 

>	G(x,n=10,N=11):=evalf(A[0]/2+ sum(g(k,11)*(A[k]*cos(k*x)+B[k]*sin(k*x)),k=1..10),4);

 

 

(58)

 

 

>	plot({M(x),y(x,n=10),G(x,n=10,N=11)},x=-Pi..Pi,color=black, title="FOURIER-Series  y(x, n = 10)  and  its Smoothing");

 

 

 

 

>

 

R?sum? 

 

In this worksheet it has been illustrated that Maplesoft furnishes powerful tools in finding best approximations to several piecewise continuous functions.