Posted: 08/24/22

Due: 09/01/22, 11:55PM

Track 1: 1,2,3,6. Track 2: 1-6.

∴ 

using Printf;

∴ 

1. Construct a convergence plot in logarithmic coordinates for the following continued fraction approximation of $e$

   ${\displaystyle e=2+\frac{1}{1+\frac{1}{2+\frac{2}{3+\frac{3}{\ddots }}}}}$

   (1)

   Identify the terms in the general expression of a continued fraction

   ${\displaystyle F_n=b_0+\stackrel{n}{\underset{k=1}{\text{K}}}\frac{a_k}{b_k}.}$

   Compare with the additive approximation from a McLaurin series

   ${\displaystyle e^x=1+\frac{x}{1!}+\frac{x^2}{2!}+\cdots }$

   Estimate the rate and order of convergence for both approximations.

   Solution. Rewrite (1) as

   ${\displaystyle e=2+\frac{1}{f},f=1+\frac{1}{2+\frac{2}{3+\frac{3}{\ddots }}},}$

   and introduce the continued fraction sequence ${\left\{F_n\right\}}_{n\in \mathbb{N}}$,

   ${\displaystyle F_n=b_0+\stackrel{n}{\underset{k=1}{\text{K}}}\frac{a_k}{b_k}=b_0+\frac{a_1}{b_1+},}$ ${\displaystyle F_0=b_0,F_1=b_0+\frac{a_1}{b_1},F_2=b_0+\frac{a_1}{b_1+},F_3=b_0+\frac{a_1}{b_1+}}$

   with assumed limit ${lim}_{n\to \infty }=f$. Writing out the first few terms

   ${\displaystyle F_0=1,F_1=1+\frac{1}{2},F_2=1+\frac{1}{2+\frac{2}{3}},F_3=1+\frac{1}{2+\frac{2}{3+\frac{3}{4}}},\dots ,}$

   leads to identification of coefficients as

   ${\displaystyle b_k=k+1,a_k=k.}$

   $\circ$Define a function $f$ to compute $F_n$ and return $E_n=2+1/F_n\to e$.

   ∴	function f(n,dbg=false) b=n+1; F=1.0*b; E=2+1/F for k=n:-1:1 b=k; a=k; F = b + a/F; E=2+1/F if dbg @printf("k=%d b_k=%d a_k=%d F_k=%f E_k=%f\n",k,b,a,F,E) end end return E end;

   ∴

   Test the function

   ∴

   f(5,true)

   k=5 b_k=5 a_k=5 F_k=5.833333 E_k=2.171429

   k=4 b_k=4 a_k=4 F_k=4.685714 E_k=2.213415

   k=3 b_k=3 a_k=3 F_k=3.640244 E_k=2.274707

   k=2 b_k=2 a_k=2 F_k=2.549414 E_k=2.392247

   k=1 b_k=1 a_k=1 F_k=1.392247 E_k=2.718263

   $2.718263331760264$

   ∴

   exp(1.0)

   $2.718281828459045$

   Comparison with the builtin $e^1$ shows 5-digit accuracy for $n=5$.

   ∴

   $\circ$Define a function $g$ to compute

   ${\displaystyle M_n=1+\frac{1}{1!}+\frac{1}{2!}+\cdots +\frac{1}{n!}.}$

   ∴	function g(n) M=1.0; fact=1.0 for k=1:n fact=k*fact M=M+1/fact end return M end

   g

   ∴

   Test the function

   ∴

   g.(1:5)

   ${\displaystyle \left[\begin{array}{c}2.0\\ 2.5\\ 2.6666666666666665\\ 2.708333333333333\\ 2.7166666666666663\end{array}\right]}$

   (2)

   Note that only 3 accurate digits are obtained for $n=5$.

   ∴

   The convergence behavior of the two approximations $E_n,M_n$ is shown in Fig. 1.

   $\circ$

   Figure 1. Convergence of continued fraction and additive approximations of $e$.

   ∴	N=15; n=1:N; e=exp(1.0); x=log.(n);

   ∴	errE=log.(abs.(f.(n) .- e)); errM=log.(abs.(g.(n) .- e));

   ∴	clf(); plot(x,errE,"-o",x,errM,"-x"); xlabel("log(n)"); ylabel("log(err)"); title("e Sequence convergence"); grid("on"); plot([0,1],[-10,-11],"b-"); plot([0,1],[-10,-12],"g-"); legend(["errE","errM","1st","2nd"]);

   ∴	savefig(homedir()*"/courses/MATH661/images/H01Fig01.eps")

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   ∴

   The definition of order $p$ and rate $r$ of convergence

   ${\displaystyle {lim}_{n\to \infty }\frac{|x_{n+1}-x|}{{|x_n-x|}^p}=r,}$

   is based upon the assumption of power-law decrease of the error $e_n=|x_n-x|$,

   ${\displaystyle e_{n+1}\sim e_n^p\iff e_{n+1}=r{e}_n^p.}$

   In log-coordinates this assumption leads to a straight line representation

   ${\displaystyle \log e_{n+1}=p\log e_n+\log r.}$

   The $(\log n,\log e_n)$ representation in Fig. 1 does not allow direct identification of an order of convergence. However a plot of $(\log e_n,\log e_{n+1})$ is easily constructed (Fig. 2), and shows $p\cong 1$, $\ln r\cong -2\Rightarrow r=e^{-2}\cong 0.135$.

   $\circ$

   Figure 2. Order-of-convergence plot for approximations of $e$. Visual estimation indicates $p=1$, linear sequence convergence.

   ∴	N=15; n=1:N; e=exp(1.0); x=log.(n);

   ∴	errE=log.(abs.(f.(n) .- e)); errM=log.(abs.(g.(n) .- e));

   ∴	clf(); plot(errE[1:N-1],errE[2:N],"-o",errM[1:N-1],errM[2:N],"-x"); xlabel("log(e_n)"); ylabel("log(e_(n+1))"); title("e Sequence convergence"); grid("on"); plot([-25,-20],[-25,-20],"b-"); plot([-25,-20],[-25,-15],"g-"); legend(["errE","errM","1st","2nd"]);

   ∴	savefig(homedir()*"/courses/MATH661/images/H01Fig02.eps")

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   ∴

   exp(-2.0)

   $0.1353352832366127$

   ∴

2. Apply convergence acceleration to both the above approximations of $e$. Construct the convergence plot of the accelerated sequences, and estimate the new rate and order of convergence.

   Solution. Since both sequences exhibit linear convergence the Aitken formula

   ${\displaystyle a_n=x_n-\frac{{(x_n-x_{n-1})}^2}{x_n-2x_{n-1}+x_{n-2}}}$

   is applicable. The resulting accelerated convergence plot is shown in Fig. 3. Convergence acceleration to second order is observed for a small range of errors, $\ln e\in [-7,-4]$. For smaller errors, the floating point system cannot separate the small differences appearing in the Aitken correction.

   $\circ$

   Figure 3. Convergence of Aitken acceleration of $e$ approximation sequences

   ∴	function Aitken(x) a=copy(x); N=length(x) for n=3:N Δx = x[n] - x[n-1] Δ2x = x[n] - 2*x[n-1] + x[n-2] a[n] = x[n] - Δx^2/Δ2x end return a end;

   ∴	N=6; n=1:N; E=f.(n); M=g.(n); aE=Aitken(E); aM=Aitken(M);

   ∴	erraE=log.(abs.(aE .- e)); erraM=log.(abs.(aM .- e));

   ∴

   clf(); plot(errE[1:N-1],errE[2:N],"-ob"); plot(errM[1:N-1],errM[2:N],"-xb"); plot(erraE[1:N-1],erraE[2:N],"-og"); plot(erraM[1:N-1],erraM[2:N],"-xg"); xlabel("log(e_n)"); ylabel("log(e_(n+1))"); title("e Sequence convergence"); grid("on"); plot([-10,-5],[-10,-5],"b-"); plot([-10,-5],[-10,0],"g-"); legend(["errE","errM","erraE","erraM","1st","2nd"]);

   ∴	savefig(homedir()*"/courses/MATH661/images/H01Fig03.eps")

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   ∴

3. Completely state the mathematical problem of taking the $n^{\mathrm{th}}$ root of a positive real, $n\in \mathbb{N}$. Find the absolute and relative condition numbers.

   Solution. First, assume $n>0$ fixed leading to the problem $f:{ℝ}_{+}\to ℝ$, $f\left(x\right)=x^{1/n}$. The absolute condition number is

   ${\displaystyle \hat{\kappa }={lim}_{\epsilon \to 0}{sup}_{\left|\delta x\right|⩽\epsilon }\frac{||f(x+\delta x)-f\left(x\right)||}{||\delta x||}.}$

   The condition number furnishes a bound for the change in the solution upon a change in the input

   ${\displaystyle ||f(x+\delta x)-f\left(x\right)||⩽\hat{\kappa }||\delta x||.}$

   Using the absolute value norm $||x||=\left|x\right|$ for $x\in {ℝ}_{+}$ gives

   ${\displaystyle \hat{\kappa }={lim}_{\epsilon \to 0}{sup}_{\left|\delta x\right|⩽\epsilon }\frac{|f(x+\delta x)-f\left(x\right)|}{\left|\delta x\right|}=\left|\frac{df\left(x\right)}{dx}\right|=\frac{1}{n}{x}^{\frac{1}{n}-1}=\frac{1}{n}\frac{1}{x^{(n-1)/n}},\mathrm{for}x>0.}$

   Consider some simple cases:

   • $n=1$, $f\left(x\right)=x$, $\hat{\kappa }=1$, hence perturbations in the input are not amplified

     ${\displaystyle f\left(x\right)=x,f(x+\delta x)=x+\delta x,f(x+\delta x)-f\left(x\right)=\delta x.}$

     $\circ$The problem is well-conditioned.

     ∴	f(x)=x; δf(x,δx)=f(x+δx)-f(x);

     ∴	δx=0.1; δf.(0:0.2:1,δx)/δx

     ${\displaystyle \left[\begin{array}{c}1.0\\ 1.0000000000000002\\ 0.9999999999999998\\ 0.9999999999999998\\ 0.9999999999999998\\ 1.0000000000000009\end{array}\right]}$

     (3)

   • $n=2$, $f\left(x\right)=\sqrt{x}$,

     ${\displaystyle \hat{\kappa }=\frac{1}{2\sqrt{x}}.}$

     $\circ$at $x=1/4$, $\hat{\kappa }=1$, indicating no amplification of input perturbation

     ∴	f(x)=sqrt(x); δf(x,δx)=f(x+δx)-f(x);

     ∴	δx=10.0 .^(-3:-1); δf.(0.25,δx)./δx

     ${\displaystyle \left[\begin{array}{c}1.0000000000000009\\ 1.0000000000000009\\ 0.9999999999999998\end{array}\right]}$

     (4)

     ∴

     $\circ$as $x\to \infty$, $\hat{\kappa }\to 0$, indicating input perturbations have negligible effect upon output. This indicates incorrect identification of the variables in a problem.

     ∴	f(x)=sqrt(x); δf(x,δx)=f(x+δx)-f(x); X=floatmax(Float64);

     ∴	δx=10.0 .^(-3:-1); δf.(X,δx)./δx

     ${\displaystyle \left[\begin{array}{c}0.0\\ 0.0\\ 0.0\end{array}\right]}$

     (5)

     ∴	δx=10.0 .^(-3:-1); δf.(X/10.0,δx)./δx

     ${\displaystyle \left[\begin{array}{c}0.0\\ 0.0\\ 0.0\end{array}\right]}$

     (6)

     ∴

     $\circ$as $x\to 0_{+}$, $\hat{\kappa }\to \infty$, indicating small input perturbations have arbitrarily large effects upon output. This indicates an ill-posed problem

     ∴	f(x)=sqrt(x); δf(x,δx)=f(x+δx)-f(x); x=floatmin(Float64);

     ∴	δx=10.0 .^(-3:-1); δf.(x,δx)./δx

     ${\displaystyle \left[\begin{array}{c}31.62277660168379\\ 10.0\\ 3.162277660168379\end{array}\right]}$

     (7)

     ∴	δx=10.0 .^(-15:-12); δf.(x,δx)./δx

     ${\displaystyle \left[\begin{array}{c}3.1622776601683788e7\\ 1.0e7\\ 3.162277660168379e6\\ 1.0e6\end{array}\right]}$

     (8)

     ∴	δx=10.0 .^(-21:-18); δf.(x,δx)./δx

     ${\displaystyle \left[\begin{array}{c}3.1622776601683792e10\\ 1.0e10\\ 3.16227766016838e9\\ 1.0e9\end{array}\right]}$

     (9)

     ∴

   In general, the conditioning of the $n^{\mathrm{th}}$ root operation

   ${\displaystyle \hat{\kappa }=\frac{1}{n}\frac{1}{x^{(n-1)/n}},\mathrm{for}x>0}$

   indicates ill-conditioning for $x\to 0$, well conditioning for $x\sim 1$, incorrect model for $x\to \infty$.

   Consider now what happens when $n$ is also allowed to vary. The definition of the condition number cannot be directly applied since the limit process is not defined for $n\in \mathbb{N}$. One can however consider the problem $g:{ℝ}_{+}\times {ℝ}_{+}\to ℝ$

   ${\displaystyle g(x,y)=x^{1/y},}$

   and notice that $h:{ℝ}_{+}\times \mathbb{N}\to ℝ$, $h(x,n)=x^{1/n}$ is a restriction of $g$. The condition number of $h$ can be inferred from that of $g$

   ${\displaystyle \hat{\kappa }={lim}_{\epsilon \to 0}{sup}_{||\delta z||⩽\epsilon }\frac{|g(x+\delta x,y+\delta y)-g(x,y)|}{||\delta z||},z={\left[\begin{array}{ll}\delta x& \delta y\end{array}\right]}^T.}$

   Use the $\infty$-norm for $\delta z\in {ℝ}_{+}^2$, ${||\delta z||}_{\infty }=max(\delta x,\delta y)$. When approaching zero perturbation, $||\delta z||\to 0$ above the first bisector, the inequality

   ${\displaystyle \delta x<\delta y⩽\epsilon }$

   implies

   ${\displaystyle \hat{\kappa }=\left|\frac{\partial g}{\partial y}\right|.}$

   Conversely, for $||\delta z||\to 0$ below the first bisector

   ${\displaystyle \hat{\kappa }=\left|\frac{\partial g}{\partial x}\right|.}$

   In general, for some arbitrary path to approach zero,

   ${\displaystyle \hat{\kappa }=max(\left|\frac{\partial g}{\partial x}\right|,\left|\frac{\partial g}{\partial y}\right|)=max(\frac{1}{y}\frac{1}{x^{(y-1)/y}},\frac{x^{1/y}\ln x}{y^2}).}$

   When restricted to $y=n\in \mathbb{N}$,

   ${\displaystyle \hat{\kappa }=max(\frac{1}{n}\frac{1}{x^{(n-1)/n}},\frac{x^{1/n}\ln x}{n^2}).}$

   $\circ$

   Figure 4. Conditioning of $n^{\mathrm{th}}$-root $x^{1/n}$ with perturbations allowed in $n$.

   ∴	κ(x,n)=max(1/n/x^((n-1)/n),x^(1/n)*log(x)/n^2); clf();

   ∴	x=0.001:0.2:11; κ2=log.(κ.(x,2)); κ3=log.(κ.(x,3));

   ∴	κ4=log.(κ.(x,4)); plot(x,κ2,x,κ3,x,κ4); grid("on");

   ∴	xlabel("x"); ylabel("ln κ"); title("Condition of nth root");

   ∴	legend(["2","3","4"]);

   ∴	savefig(homedir()*"/courses/MATH661/images/H01Fig04.eps")

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   The PostScript backend does not support transparency; partially transparent artists will be rendered opaque.

   ∴

4. Completely state the mathematical problem of solving the initial value problem for an ordinary differential equation of first order. Find the absolute and relative condition numbers.

   Solution. The IVP

   ${\displaystyle y^{\text{'}}=f\left(y\right),y\left(0\right)=y_0}$

   is formulated as the mathematical problem

   ${\displaystyle F:C^{0,1}\left(ℝ\right)\times ℝ\to C\left(ℝ\right),}$

   that when evaluated for some slope function $f$ and initial condition $y_0$ gives the integral curve $y:[0,a)\to ℝ$, $a>0$. $F(f,y_0)=y$. In the above $C\left(ℝ\right)$ is the space of continuous functions and $C^{0,1}$ is the space of Lipschitz continuous functions.

   As in the $n^{\mathrm{th}}$-root problem there are two inputs to $F$, and it is useful to start with the case in which the slope function $f$ is fixed and only $y_0$ varies. The Lyapunov exponent $L$ is defined as

   ${\displaystyle \delta y\left(t\right)\cong e^{Lt}\delta y_0,}$

   to characterize this case and the condition number is simply

   ${\displaystyle \hat{\kappa }\left(t\right)=e^{Lt},}$

   i.e., the condition number and the Lyapunov exponent express the same concept.

   The condition number for the case of variation of the slope function requires a concept of taking a derivative with respect to a function $f$, a generalization of the calculus concept of taking a derivative with respect to a variable. This is known as a functional derivative and can be defined in the Fréchet sense for normed spaces and in the Gateaux sense for Banach spaces.

5. Completely state the mathematical problem of finding the roots of a cubic polynomial. Find the absolute and relative condition numbers.

   Solution. Consider the cubic $x^3+a_2{x}^2+a_{1}x+a_0=0$ with roots $x_1,x_2,x_3$ related to the polynomial coefficients by the Vieta relations

   ${\displaystyle x_1+x_2+x_3=-a_2,x_1{x}_2+x_1{x}_3+x_2{x}_3=a_1,x_1{x}_2{x}_3=-a_0.}$

   The mathematical problem of finding the roots of the cubic is $f:{ℝ}^3\to {ℂ}^3$

   ${\displaystyle f\left(𝒂\right)=𝒙,𝒂=\left[\begin{array}{l}{a}_0\\ a_1\\ a_2\end{array}\right],𝒙=\left[\begin{array}{l}{x}_1\\ x_2\\ x_3\end{array}\right]}$

   Consider the effect of small changes $\delta 𝒂$ upon the roots by taking differentials

   ${\displaystyle \begin{array}{rcl}\delta x_1+\delta x_2+\delta x_3& =& -\delta a_2\\ (x_2+x_3)\delta x_1+(x_3+x_1)\delta x_2+(x_1+x_2)\delta x_3& =& \delta a_1\\ x_2{x}_3\delta x_1+x_3{x}_1\delta x_2+x_1{x}_2\delta x_3& =& -\delta a_0\end{array}}$

   This is a linear system for $\delta 𝒙$ with matrix

   ${\displaystyle 𝑩=\left[\begin{array}{lll}1& 1& 1\\ x_2+x_3& x_3+x_1& x_1+x_2\\ x_2{x}_3& x_3{x}_1& x_1{x}_2\end{array}\right].}$

   The condition number for $f$ is the maximal amplification by the matrix $𝑩$ or

   ${\displaystyle \hat{\kappa }=||𝑩||.}$

   Consider some specific cases:

   • For $x_1=x_2=x_3=\xi$,

     ${\displaystyle 𝑩=\left[\begin{array}{lll}𝒃& 𝒃& 𝒃\end{array}\right],𝒃=\left[\begin{array}{l}1\\ 2\xi \\ {\xi }^2\end{array}\right]}$

6. Numerically compare the approximation of $b\left(t\right)=t(\pi -t)(2\pi -t)$ by linear combination of $𝒯=\{\sin t,\sin 2t,\sin 3t,\dots \}$ with that of linear combinations of $ℰ=\{1,e^t,e^{2t},e^{3t},\dots \}$. Present a study of the aproximation error as the number of terms in the linear combination increases. Estimate the order of convergence in both cases.

   Solution. The approximation is stated as

   ${\displaystyle b\left(t\right)\cong {\hat{b}}_n\left(t\right)=\sum _{k=1}^n{c}_k{a}_k\left(t\right)}$

   with the basis functions chosen either as $a_k\left(t\right)=\sin (kt)$ or $a_k\left(t\right)=e^{(k-1)t}$. The approximation error $e\left(t\right)=b\left(t\right)-{\hat{b}}_n\left(t\right)$ can be measured in various norms, e.g., the 2-norm

   ${\displaystyle {\epsilon }^2={||e\left(t\right)||}_2^2=\underset{0}{\overset{2\pi }{\int }}{(b\left(t\right)-{\hat{b}}_n\left(t\right))}^{2}dt\cong \frac{2\pi }{m}\sum _{i=1}^m{(b\left(t_i\right)-{\hat{b}}_n\left(t_i\right))}^2.}$

   $•$Based upon the code from Fig.1 in L04, define a function that returns the approximation error for given basis set $a_k\left(t\right)$, number of terms $n$, and evaluation points $m$.

   ∴	function err(m,n,a,dbg=false) h=2.0*π/m; j=1:m; t=h*j; A = a.(1,t) for k=2:n A = [A a.(k,t)] end if dbg return A end bt=t.*(π.-t).*(2*π.-t) x=A\bt; b=A*x; return norm(b-bt)*(2*pi)/m end

   err

   ∴

   $•$Define the two basis sets of interest

   ∴	s(k,t) = sin(k*t); e(k,t)=exp((k-1)*t);

   ∴

   $\circ$Verify the err function constructs the expected matrix. Note that the basis function is passed as an argument.

   ∴	err(8,1,s,true)

   ${\displaystyle \left[\begin{array}{c}0.7071067811865475\\ 1.0\\ 0.7071067811865476\\ 1.2246467991473532e-16\\ -0.7071067811865475\\ -1.0\\ -0.7071067811865477\\ -2.4492935982947064e-16\end{array}\right]}$

   (10)

   ∴	err(8,1,e,true)

   ${\displaystyle \left[\begin{array}{c}1.0\\ 1.0\\ 1.0\\ 1.0\\ 1.0\\ 1.0\\ 1.0\\ 1.0\end{array}\right]}$

   (11)

   ∴	err(8,2,s,true)

   ${\displaystyle \left[\begin{array}{cc}0.7071067811865475& 1.0\\ 1.0& 1.2246467991473532e-16\\ 0.7071067811865476& -1.0\\ 1.2246467991473532e-16& -2.4492935982947064e-16\\ -0.7071067811865475& 1.0\\ -1.0& 3.6739403974420594e-16\\ -0.7071067811865477& -1.0\\ -2.4492935982947064e-16& -4.898587196589413e-16\end{array}\right]}$

   (12)

   ∴	err(8,2,e,true)

   ${\displaystyle \left[\begin{array}{cc}1.0& 2.1932800507380152\\ 1.0& 4.810477380965351\\ 1.0& 10.550724074197761\\ 1.0& 23.140692632779267\\ 1.0& 50.75401951173493\\ 1.0& 111.31777848985621\\ 1.0& 244.15106285427498\\ 1.0& 535.4916555247646\end{array}\right]}$

   (13)

   ∴

   The numerical values within the matrix $𝑨$ are hard to interpret for large $m$, hence plot the columns in Fig. 5. The basis fubnctino plots already indicate that the exponential family is likely to lead to bad approximations since with respect to $e^{(n-1)t}$, all $e^{kt}$ for $k<n-1$ are negligibly small, hence $𝑨$ is likely to have only one independent column vector.

   $\circ$

   Figure 5. Left: Sine basis functions. Right: Exponential basis functions.

   ∴	m=256; n=5; As = err(m,n,s,true); Ae = err(m,n,e,true);

   ∴	clf(); h=2.0*π/m; j=1:m; t=h*j;

   ∴	for k=1:n plot(t,As[:,k]) end

   ∴	grid("on"); xlabel("t"); ylabel("sin(k*t)");

   ∴	title("Sine basis functions");

   ∴	cd(homedir()*"/courses/MATH661/images");

   ∴	savefig("H01Fig05a.eps");

   ∴	clf(); grid("on"); xlabel("t"); ylabel("exp(k*t)");

   ∴	title("Exp basis functions");

   ∴	for k=1:n plot(t,Ae[:,k]) end

   ∴	savefig("H01Fig05b.eps");

   ∴

   $\circ$Test the err function for larger $m$ values (more samples).

   ∴	err(1000,10,s)

   $0.0020986874338193287$

   ∴	err(1000,10,e)

   $1.4507633675453773$

   ∴	err(4,1,s,true)

   ${\displaystyle \left[\begin{array}{c}0.0\\ 1.0\\ 1.2246467991473532e-16\\ -1.0\end{array}\right]}$

   (14)

   ∴

   The convergence behavior is shown in Fig. 6. As expected as the number of terms in the linear combination increases the sine approximation converges, while that for the exponential basis has a constant error (numerically, the rank of the matrix remains 1).

   $\circ$

   Figure 6. Convergence of sine, exp basis approximation of $b\left(t\right)$ with increasing number of terms.

   ∴ m=1000; n=10:10:100; xn = log.(1.0*n); ∴ errs = err.(m,n,s); lgerrs = log.(errs); ∴ erre = err.(m,n,e); lgerre = log.(erre); ∴ clf(); plot(xn,lgerrs,"ok",xn,lgerre,"rx"); ∴ xlabel("log(n)"); ylabel("log(err)"); grid("on") ∴ title("Convergence behavior for sine, exp basis"); ∴ legend(["sine","exp"]); ∴ cd(homedir()*"/courses/MATH661/images"); ∴ savefig("H01Fig06.eps"); The PostScript backend does not support transparency; partially transparent artists will be rendered opaque. The PostScript backend does not support transparency; partially transparent artists will be rendered opaque. ∴ ∴ ∴ ∴ ∴ ∴