UPSC Maths Optional · Paper II · Full Proofs + PYQs 2000–2026
If \(f(r)=0\), \(f'(r)\ne0\), and \(f''\) is bounded near \(r\), then NR iteration \(x_{n+1}=x_n-\dfrac{f(x_n)}{f'(x_n)}\) converges quadratically: \(e_{n+1}\approx\dfrac{f''(r)}{2f'(r)}e_n^2\) where \(e_n=x_n-r\).
Step 1 Taylor expand \(f\) around \(r\): \(f(x_n)=f(r)+e_nf'(r)+\dfrac{e_n^2}{2}f''(\xi)\) for some \(\xi\) between \(r\) and \(x_n\). Since \(f(r)=0\): \(f(x_n)=e_nf'(r)+\dfrac{e_n^2}{2}f''(\xi)\).
Step 2 Similarly: \(f'(x_n)=f'(r)+e_nf''(\eta)\) for some \(\eta\).
Step 3 \(e_{n+1}=x_{n+1}-r=x_n-\dfrac{f(x_n)}{f'(x_n)}-r=e_n-\dfrac{e_nf'(r)+e_n^2f''(\xi)/2}{f'(r)+e_nf''(\eta)}=\dfrac{e_n^2f''(\xi)/2}{f'(r)+e_nf''(\eta)}\approx\dfrac{f''(r)}{2f'(r)}e_n^2\) (for small \(e_n\)). \(\blacksquare\)
If \(f\) is continuous on \([a,b]\) and \(f(a)f(b)<0\), bisection converges with \(|e_n|\le\dfrac{b-a}{2^n}\).
After \(n\) steps, the interval containing the root has length \((b-a)/2^n\). The midpoint \(x_n\) satisfies \(|x_n-r|\le(b-a)/2^{n+1}\). As \(n\to\infty\), \(|e_n|\to0\) (linear convergence with rate 1/2). \(\blacksquare\)
\(y(x)=y_0+s\Delta y_0+\dfrac{s(s-1)}{2!}\Delta^2y_0+\dfrac{s(s-1)(s-2)}{3!}\Delta^3y_0+\cdots\) where \(s=(x-x_0)/h\).
Use the shift operator \(E=1+\Delta\). At \(x=x_0+sh\): \(y(x)=E^sy_0=(1+\Delta)^sy_0\). Expand binomially: \(y=\left(1+s\Delta+\binom{s}{2}\Delta^2+\binom{s}{3}\Delta^3+\cdots\right)y_0=y_0+s\Delta y_0+\dfrac{s(s-1)}{2!}\Delta^2y_0+\cdots\). \(\blacksquare\)
Error in Lagrange interpolation through \((n+1)\) points: \(E(x)=\dfrac{f^{(n+1)}(\xi)}{(n+1)!}\prod_{i=0}^n(x-x_i)\) for some \(\xi\) in the interval.
Let \(\omega(x)=\prod_{i=0}^n(x-x_i)\). Define \(g(t)=f(t)-P_n(t)-\lambda\omega(t)\) where \(\lambda\) is chosen so that \(g(x)=0\). Then \(g\) has \((n+2)\) zeros: \(x_0,x_1,\ldots,x_n,x\). By Rolle's theorem applied \((n+1)\) times: \(g^{(n+1)}(\xi)=0\) for some \(\xi\). So \(f^{(n+1)}(\xi)-(n+1)!\lambda=0\Rightarrow\lambda=\dfrac{f^{(n+1)}(\xi)}{(n+1)!}\). Since \(g(x)=f(x)-P_n(x)-\lambda\omega(x)=0\): \(f(x)-P_n(x)=\dfrac{f^{(n+1)}(\xi)}{(n+1)!}\omega(x)\). \(\blacksquare\)
\(\int_a^b f(x)\,dx\approx\dfrac{h}{3}[f(x_0)+4f(x_1)+f(x_2)]\) where \(h=(b-a)/2\). Error: \(-\dfrac{h^5}{90}f^{(4)}(\xi)\).
Step 1 Fit quadratic through \((x_0,f_0),(x_1,f_1),(x_2,f_2)\) with \(x_0=a,x_1=(a+b)/2,x_2=b\): \(P_2(x)=\dfrac{(x-x_1)(x-x_2)}{(x_0-x_1)(x_0-x_2)}f_0+\dfrac{(x-x_0)(x-x_2)}{(x_1-x_0)(x_1-x_2)}f_1+\dfrac{(x-x_0)(x-x_1)}{(x_2-x_0)(x_2-x_1)}f_2\).
Step 2 Integrate: with \(s=(x-x_1)/h\), \(\int_{-1}^1(s^2-s^2)ds... \) Using \(\int_{x_0}^{x_2}P_2\,dx=h\!\int_{-1}^1\!\left[\dfrac{s(s-1)}{2}f_0-\dfrac{s^2-1}{1}f_1+\dfrac{s(s+1)}{2}f_2\right]ds\). \(\int_{-1}^1s(s-1)ds=2/3-0\)... Computing each coefficient: \(f_0:\int_{-1}^1\frac{s(s-1)}{2}ds=\frac{1}{2}\cdot\frac{2}{3}=\frac{1}{3}\). \(f_1:\int_{-1}^1(1-s^2)ds=\frac{4}{3}\). \(f_2:\int_{-1}^1\frac{s(s+1)}{2}ds=\frac{1}{3}\). Total: \(h(f_0/3+4f_1/3+f_2/3)=\frac{h}{3}(f_0+4f_1+f_2)\). \(\blacksquare\)
For \(y'=f(x,y),y(x_0)=y_0\). RK4 step: \(y_{n+1}=y_n+\dfrac{h}{6}(k_1+2k_2+2k_3+k_4)\) where \(k_1=f(x_n,y_n),k_2=f(x_n+h/2,y_n+hk_1/2),k_3=f(x_n+h/2,y_n+hk_2/2),k_4=f(x_n+h,y_n+hk_3)\).
The exact solution satisfies \(y(x_{n+1})=y(x_n)+hy'(x_n)+\dfrac{h^2}{2}y''(x_n)+\dfrac{h^3}{6}y'''(x_n)+\dfrac{h^4}{24}y^{(4)}(x_n)+O(h^5)\). The RK4 coefficients are chosen so that the weighted combination \(\frac{h}{6}(k_1+2k_2+2k_3+k_4)\) matches these Taylor terms up to \(O(h^4)\). The weights \(1,2,2,1\) (from Simpson's 1/3 rule applied to slope evaluation) ensure 4th-order accuracy. Local truncation error: \(O(h^5)\), global error: \(O(h^4)\). \(\blacksquare\)
Find \(\sqrt{3}\) using Newton-Raphson starting from \(x_0=2\). Perform 2 iterations.
Set \(f(x)=x^2-3\). \(f'(x)=2x\). NR: \(x_{n+1}=x_n-\dfrac{x_n^2-3}{2x_n}=\dfrac{x_n^2+3}{2x_n}\).
Iter 1: \(x_1=\dfrac{4+3}{4}=\dfrac{7}{4}=1.75\).
Iter 2: \(x_2=\dfrac{(7/4)^2+3}{2\cdot7/4}=\dfrac{49/16+48/16}{7/2}=\dfrac{97/16}{7/2}=\dfrac{97}{56}\approx1.73214\). True value: \(\sqrt{3}\approx1.73205\). Error \(\approx0.00009\). \(\blacksquare\)
Find the root of \(xe^x=1\) near \(x_0=0.5\) using 2 iterations of Newton-Raphson.
\(f(x)=xe^x-1,f'(x)=e^x+xe^x=(1+x)e^x\). \(x_1=0.5-\dfrac{0.5e^{0.5}-1}{1.5e^{0.5}}=0.5-\dfrac{0.5\times1.6487-1}{1.5\times1.6487}=0.5-\dfrac{0.8244-1}{2.4731}=0.5-\dfrac{-0.1756}{2.4731}=0.5+0.0710=0.5710\). \(x_2\approx0.5671\) (converges to \(0.5671\)). \(\blacksquare\)
Find a root of \(f(x)=x^3-4x-9=0\) in \([2,3]\) using 3 iterations of bisection.
\(f(2)=8-8-9=-9<0\), \(f(3)=27-12-9=6>0\). Root in \([2,3]\). Iter 1: \(x_1=2.5,f(2.5)=15.625-10-9=-3.375<0\). Root in \([2.5,3]\). Iter 2: \(x_2=2.75,f(2.75)=20.797-11-9=0.797>0\). Root in \([2.5,2.75]\). Iter 3: \(x_3=2.625,f(2.625)\approx18.07-10.5-9=-1.43<0\). Root in \([2.625,2.75]\). Approximate root \(\approx2.6875\). \(\blacksquare\)
Given \((0,1),(1,3),(2,7),(3,13)\), find \(f(1.5)\) using Lagrange interpolation.
These points lie on \(f(x)=x^2+x+1\). Let's verify: \(0+0+1=1\checkmark,1+1+1=3\checkmark,4+2+1=7\checkmark,9+3+1=13\checkmark\). Lagrange with 4 points (or just use formula): \(f(1.5)=(1.5)^2+1.5+1=2.25+1.5+1=4.75\). Lagrange formula: \(L_0=\dfrac{(1.5-1)(1.5-2)(1.5-3)}{(0-1)(0-2)(0-3)}=\dfrac{(0.5)(-0.5)(-1.5)}{(-6)}=\dfrac{0.375}{6}... \). Direct use of polynomial gives \(4.75\). \(\blacksquare\)
Using the following table, find \(f(1.5)\): \(x: 1,2,3,4\); \(y: 1,8,27,64\) (i.e., \(y=x^3\)).
Forward differences: \(\Delta y: 7,19,37\). \(\Delta^2y: 12,18\). \(\Delta^3y: 6\). \(\Delta^4y=0\). Here \(x_0=1,h=1,s=(1.5-1)/1=0.5\). Newton's forward: \(f(1.5)=1+0.5\times7+\dfrac{0.5(-0.5)}{2}\times12+\dfrac{0.5(-0.5)(-1.5)}{6}\times6=1+3.5-1.5+0.375=3.375\). And \(1.5^3=3.375\). \(\checkmark\) \(\blacksquare\)
Evaluate \(\int_0^1\dfrac{1}{1+x}\,dx\) using Simpson's rule with \(h=0.25\). Compare with exact.
\(x: 0, 0.25, 0.5, 0.75, 1\). \(y=1/(1+x): 1, 0.8, 0.6\overline{6}, 0.5714, 0.5\). Simpson's: \(\dfrac{0.25}{3}[y_0+4y_1+2y_2+4y_3+y_4]=\dfrac{0.25}{3}[1+4(0.8)+2(0.6\overline6)+4(0.5714)+0.5]=\dfrac{0.25}{3}[1+3.2+1.333+2.286+0.5]=\dfrac{0.25}{3}\times8.319=\dfrac{0.25\times8.319}{3}=0.6933\). Exact: \(\ln2=0.6931\). Error: \(0.0002\). \(\blacksquare\)
Evaluate \(\int_0^6f(x)\,dx\) using trapezoidal rule with \(h=1\) for: \(x: 0,1,2,3,4,5,6\) and \(y: 1,1.5,2,3,4,5,6\).
\(\int\approx\dfrac{1}{2}[y_0+2(y_1+y_2+y_3+y_4+y_5)+y_6]=\dfrac{1}{2}[1+2(1.5+2+3+4+5)+6]=\dfrac{1}{2}[1+2(15.5)+6]=\dfrac{1}{2}[38]=19\). \(\blacksquare\)
Solve \(y'=x+y,y(0)=1\) using RK4 with \(h=0.1\). Find \(y(0.1)\).
\(f(x,y)=x+y\). At \(x_0=0,y_0=1\): \(k_1=f(0,1)=1\). \(k_2=f(0.05,1+0.05)=f(0.05,1.05)=0.05+1.05=1.1\). \(k_3=f(0.05,1+0.05\times1.1)=f(0.05,1.055)=1.105\). \(k_4=f(0.1,1+0.1\times1.105)=f(0.1,1.1105)=1.2105\). \(y_1=1+\dfrac{0.1}{6}(1+2(1.1)+2(1.105)+1.2105)=1+\dfrac{0.1}{6}(6.6205)=1+0.11034=1.11034\). Exact: \(y=2e^x-x-1\Rightarrow y(0.1)=2e^{0.1}-1.1=2(1.10517)-1.1=1.11034\). \(\checkmark\) \(\blacksquare\)
Use Euler's method to solve \(y'=x+y,y(0)=1\) with \(h=0.2\). Find \(y(0.4)\).
\(y_{n+1}=y_n+hf(x_n,y_n)\). \(y(0.2)=1+0.2(0+1)=1.2\). \(y(0.4)=1.2+0.2(0.2+1.2)=1.2+0.28=1.48\). (Exact: \(2e^{0.4}-1.4=2(1.4918)-1.4=1.5836\); Euler has significant error due to large \(h\).) \(\blacksquare\)
Solve \(y'=x+y^2,y(0)=1\) using Heun's method with \(h=0.1\). Compute \(y(0.1)\).
Predictor: \(y^*=y_0+hf(x_0,y_0)=1+0.1(0+1)=1.1\). Corrector: \(y_1=y_0+\dfrac{h}{2}[f(x_0,y_0)+f(x_1,y^*)]=1+\dfrac{0.1}{2}[f(0,1)+f(0.1,1.1)]=1+0.05[(0+1)+(0.1+1.21)]=1+0.05[1+1.31]=1+0.05(2.31)=1+0.1155=1.1155\). \(\blacksquare\)
Find a root of \(f(x)=x^3-2x-5=0\) in \([2,3]\) using Regula Falsi. Do 2 iterations.
\(f(2)=-1,f(3)=16\). Iter 1: \(x_1=\frac{2\cdot16-3\cdot(-1)}{16-(-1)}=\frac{32+3}{17}=\frac{35}{17}\approx2.0588\). \(f(2.0588)=(2.0588)^3-2(2.0588)-5\approx8.728-4.118-5=-0.390\). Root in \([2.0588,3]\). Iter 2: \(x_2=\frac{2.0588\times16-3\times(-0.390)}{16+0.390}=\frac{32.94+1.17}{16.39}\approx2.080\). Converging to \(\approx2.0946\). \(\blacksquare\)
Solve \(x=e^{-x}\) by fixed point iteration starting with \(x_0=0.5\). Show convergence condition.
Iteration: \(x_{n+1}=e^{-x_n}\). Convergence: \(|g'(x)|=|-e^{-x}|=e^{-x}<1\) for \(x>0\). At root \(x^*\approx0.5671\): \(|g'(x^*)|=e^{-0.5671}\approx0.567<1\). \(\checkmark\) Iter: \(x_1=e^{-0.5}=0.6065,x_2=e^{-0.6065}=0.5454,x_3=0.5796,x_4=0.5604,\ldots\to0.5671\). \(\blacksquare\)
Use the secant method to find a root of \(f(x)=x^2-5\) starting with \(x_0=2,x_1=3\).
\(x_{n+1}=x_n-f(x_n)\dfrac{x_n-x_{n-1}}{f(x_n)-f(x_{n-1})}\). \(f(2)=-1,f(3)=4\). \(x_2=3-4\cdot\dfrac{3-2}{4-(-1)}=3-\dfrac{4}{5}=2.2\). \(f(2.2)=4.84-5=-0.16\). \(x_3=2.2-(-0.16)\dfrac{2.2-3}{-0.16-4}=2.2-(-0.16)\dfrac{-0.8}{-4.16}=2.2-0.031\approx2.236\). Converges to \(\sqrt5\approx2.2361\). \(\blacksquare\)
Compare the order of convergence of bisection, Newton-Raphson, and secant methods.
Bisection: linear convergence, order 1, rate 1/2. Error halves each iteration. Newton-Raphson: quadratic convergence, order 2. \(e_{n+1}\approx Ce_n^2\). Secant: superlinear, order \(\phi=(1+\sqrt5)/2\approx1.618\) (golden ratio). \(e_{n+1}\approx Ce_n^\phi\). In practice (per iteration cost): NR requires \(f'\); secant uses two function evaluations. NR most efficient per iteration; bisection always converges but slowly. \(\blacksquare\)
Describe Muller's method for finding complex roots of a polynomial.
Muller's method fits a quadratic through three points \((x_{n-2},f_{n-2}),(x_{n-1},f_{n-1}),(x_n,f_n)\) and finds its roots using the quadratic formula. Formula: \(x_{n+1}=x_n-\dfrac{2f_n}{w\pm\sqrt{w^2-4f_n\cdot\delta}}\) where \(w=\delta+\gamma,\delta=f[x_{n-1},x_n],\gamma=f[x_{n-2},x_{n-1},x_n](x_n-x_{n-1})\). Order of convergence \(\approx1.84\). Advantage: finds complex and real roots, including repeated roots. \(\blacksquare\)
Given \((0,1),(1,4),(3,40),(4,85)\), find Newton's divided difference interpolating polynomial.
Divided diff table: \(f[0]=1,f[1]=4,f[3]=40,f[4]=85\). \(f[0,1]=(4-1)/1=3\). \(f[1,3]=(40-4)/2=18\). \(f[3,4]=(85-40)/1=45\). \(f[0,1,3]=(18-3)/3=5\). \(f[1,3,4]=(45-18)/3=9\). \(f[0,1,3,4]=(9-5)/4=1\). Polynomial: \(P(x)=1+3x+5x(x-1)+1\cdot x(x-1)(x-3)=1+3x+5x^2-5x+x^3-4x^2+3x=x^3+x^2+x+1\). \(\blacksquare\)
Define Hermite interpolation and state its error formula.
Hermite interpolation: fits a polynomial matching both function values \(f(x_i)\) and derivative values \(f'(x_i)\) at nodes. For \(n+1\) nodes: degree \(\le2n+1\) polynomial. Using divided differences with repeated nodes: if \(x_0=x_1=a\): \(f[a,a]=f'(a)\). Error: \(E(x)=\dfrac{f^{(2n+2)}(\xi)}{(2n+2)!}\prod_{i=0}^n(x-x_i)^2\) for some \(\xi\). More accurate than Lagrange since derivative info reduces oscillation. \(\blacksquare\)
Explain cubic spline interpolation and state its advantages over polynomial interpolation.
Cubic spline: piecewise cubic polynomials \(S_i(x)\) on each \([x_i,x_{i+1}]\) satisfying: (1) \(S_i(x_i)=y_i,S_i(x_{i+1})=y_{i+1}\); (2) \(S_{i-1}(x_i)=S_i(x_i)\) (C0 continuity); (3) \(S'_{i-1}=S'_i\) at interior nodes (C1); (4) \(S''_{i-1}=S''_i\) (C2). Natural spline: \(S''(x_0)=S''(x_n)=0\). Advantages over high-degree polynomial: no Runge phenomenon, smooth, local changes only affect nearby intervals. \(\blacksquare\)
Given \(y(25)=0.4226,y(30)=0.5,y(35)=0.5736,y(40)=0.6428\), find \(y(38)\) using Newton's backward difference formula.
\(h=5,x_n=40,s=(38-40)/5=-0.4\). Backward differences: \(\nabla y_{40}=0.6428-0.5736=0.0692\). \(\nabla^2=0.5736-0.5+0.5-0.4226... \) \(\nabla y_{35}=0.5736-0.5=0.0736,\nabla y_{30}=0.5-0.4226=0.0774\). \(\nabla^2y_{40}=0.0692-0.0736=-0.0044\). \(\nabla^2y_{35}=0.0736-0.0774=-0.0038\). \(\nabla^3y_{40}=-0.0044-(-0.0038)=-0.0006\). \(y(38)=0.6428+(-0.4)(0.0692)+\frac{(-0.4)(-0.4+1)}{2!}(-0.0044)+\cdots=0.6428-0.02768+0.000528\approx0.6157\). \(\blacksquare\)
Given a table of \((x,y)\) values, find \(x\) for \(y=0.5\). Explain the method.
Inverse interpolation: swap roles of \(x\) and \(y\). Apply Newton's or Lagrange's formula treating \(y\) as independent variable and \(x\) as function. Alternatively: if interpolating polynomial \(P(x)=0.5\), solve numerically (Newton-Raphson on \(P(x)-0.5=0\)). For equally spaced \(y\) values: use Newton's forward difference formula with \(y_0=0.5\) to find corresponding \(x\) directly. \(\blacksquare\)
Find the maximum error in linear interpolation of \(f(x)=\sin x\) on \([0,\pi/2]\) with equally-spaced nodes of spacing \(h\).
Error in linear interpolation: \(|E(x)|\le\dfrac{h^2}{8}\max|f''(\xi)|\) for \(\xi\in[a,b]\). \(f''(x)=-\sin x\Rightarrow|f''|\le1\). Max error \(\le h^2/8\). For \(n\) equally spaced points on \([0,\pi/2]\): \(h=\pi/(2n)\), max error \(\le\dfrac{\pi^2}{32n^2}\). For error \(<0.001\): \(n>\pi/(4\sqrt{0.008})\approx27.8\Rightarrow n=28\). \(\blacksquare\)
Describe Romberg integration and explain how it achieves higher accuracy.
Romberg: systematic use of Richardson extrapolation on trapezoidal rule. Let \(T(h)\) = trapezoidal with step \(h\). Error: \(T(h)=I+c_1h^2+c_2h^4+\cdots\). Richardson: \(R_{n,1}=T(h/2^n)\), \(R_{n,m}=\dfrac{4^{m-1}R_{n,m-1}-R_{n-1,m-1}}{4^{m-1}-1}\). Each column eliminates one error term: \(R_{\cdot,2}\) eliminates \(h^2\) error (like Simpson's); \(R_{\cdot,3}\) eliminates \(h^4\); etc. Converges to \(I\) with very high order. \(\blacksquare\)
Evaluate \(\int_{-1}^1(x^4+1)\,dx\) using 2-point Gauss-Legendre quadrature. Verify with exact.
2-point GL: \(\int_{-1}^1f(x)\,dx\approx f(-1/\sqrt3)+f(1/\sqrt3)\). Nodes: \(\pm1/\sqrt3\). \(f(1/\sqrt3)=(1/3)^2+1=1/9+1=10/9\). \(f(-1/\sqrt3)=10/9\). Approximation: \(10/9+10/9=20/9\approx2.222\). Exact: \([x^5/5+x]_{-1}^1=2/5+2=12/5=2.4\). Error: \(0.178\). (2-point GL exact for polynomials of degree \(\le3\); here degree 4 — need 3-point GL for exactness.) \(\blacksquare\)
Show that the error in Simpson's 1/3 rule over \([a,b]\) with 2 panels is \(-\dfrac{h^5}{90}f^{(4)}(\xi)\).
Simpson approximates \(\int_{x_0}^{x_2}f\,dx\). Expand \(f\) in Taylor around \(x_1\): \(f(x_1+t)=f_1+tf_1'+\frac{t^2}{2}f_1''+\frac{t^3}{6}f_1'''+\frac{t^4}{24}f_1^{(4)}+\cdots\). Exact integral: \(\int_{-h}^hf(x_1+t)\,dt=2hf_1+\frac{h^3}{3}f_1''+\frac{h^5}{60}f_1^{(4)}+\cdots\). Simpson formula: \(\frac{h}{3}(f_0+4f_1+f_2)=2hf_1+\frac{h^3}{3}f_1''+\frac{h^5}{36}f_1^{(4)}+\cdots\). Error: \((\frac{1}{60}-\frac{1}{36})h^5f_1^{(4)}=-\frac{h^5}{90}f_1^{(4)}\). \(\blacksquare\)
Evaluate \(\int_0^1\int_0^1 e^{x+y}\,dx\,dy\) numerically using Simpson's rule with \(h=k=0.5\).
Grid: \(x=0,0.5,1\); \(y=0,0.5,1\). Values: \(e^{x+y}\) at 9 points. \(e^0=1,e^{0.5}=1.6487,e^1=2.7183,e^{1.5}=4.4817,e^2=7.3891\). Apply Simpson's 1/3 in both directions: weights \((1,4,1)/3\times(1,4,1)/3\). Result: \(\frac{(0.5)^2}{9}[1+4(1.6487)+2.7183+4(1.6487+4(2.7183)+4.4817)+2.7183+4(4.4817)+7.3891]\). Exact: \((e-1)^2\approx(1.7183)^2\approx2.9526\). Simpson gives \(\approx2.9527\). \(\blacksquare\)
State Weddle's rule for numerical integration and give its error.
Weddle's rule (6 panels): \(\int_{x_0}^{x_6}f\,dx\approx\dfrac{3h}{10}[f_0+5f_1+f_2+6f_3+f_4+5f_5+f_6]\). Error: \(O(h^8)\) — more accurate than Simpson's \(O(h^5)\). Requires \(n\) divisible by 6. The coefficients \((1,5,1,6,1,5,1)\times3h/10\) come from fitting a degree-6 polynomial over 7 equally-spaced points. \(\blacksquare\)
Describe the Adams-Bashforth predictor and Adams-Moulton corrector formulas.
Adams-Bashforth (explicit, predictor): \(y_{n+1}^P=y_n+\frac{h}{24}(55f_n-59f_{n-1}+37f_{n-2}-9f_{n-3})\). Order 4. Adams-Moulton (implicit, corrector): \(y_{n+1}=y_n+\frac{h}{24}(9f_{n+1}^P+19f_n-5f_{n-1}+f_{n-2})\). Order 4. Self-starting: needs 4 starting values (from RK4). Advantage: only 2 function evaluations per step (vs 4 for RK4). \(\blacksquare\)
Given \(y'=x+y,y(0)=1\), use Milne's method to find \(y(0.4)\) using starting values \(y(0)=1,y(0.1)=1.1103,y(0.2)=1.2428,y(0.3)=1.3997\).
Milne predictor: \(y_4^P=y_0+\frac{4h}{3}(2f_1-f_2+2f_3)\). \(f_1=0.1+1.1103=1.2103,f_2=0.2+1.2428=1.4428,f_3=0.3+1.3997=1.6997\). \(y_4^P=1+\frac{0.4}{3}(2(1.2103)-1.4428+2(1.6997))=1+\frac{0.4}{3}(2.4206-1.4428+3.3994)=1+\frac{0.4\times4.3772}{3}=1+0.5836=1.5836\). Milne corrector: \(y_4=y_2+\frac{h}{3}(f_2+4f_3+f_4^P)\). \(f_4^P=0.4+1.5836=1.9836\). \(y_4=1.2428+\frac{0.1}{3}(1.4428+6.7988+1.9836)=1.2428+0.3408=1.5836\). \(\blacksquare\)
What is absolute stability of a numerical method? Find the stability region of Euler's method.
Apply Euler to \(y'=\lambda y\): \(y_{n+1}=(1+h\lambda)y_n\). Method is stable iff \(|1+h\lambda|\le1\). Stability region: \(\{h\lambda\in\mathbb{C}:|1+h\lambda|\le1\}\) — a disc of radius 1 centered at \(-1\). For real negative \(\lambda\): stable iff \(h|\lambda|\le2\) (step size restriction). Implicit methods (trapezoidal, backward Euler) have much larger (or unconditional) stability regions. \(\blacksquare\)
Describe the shooting method for two-point BVP \(y''=f(x,y,y')\), \(y(a)=\alpha\), \(y(b)=\beta\).
Convert BVP to IVP: guess \(y'(a)=s\) (slope). Solve IVP \(y''=f,y(a)=\alpha,y'(a)=s\) using RK4 to get \(y(b;s)\). If \(y(b;s)=\beta\): done. If not: use root-finding (secant/NR) on \(\phi(s)=y(b;s)-\beta=0\) to find correct \(s\). For linear BVP: \(\phi(s)\) is linear in \(s\), so one interpolation step suffices. Non-linear: iterate. \(\blacksquare\)
Solve \(y''-y=0\), \(y(0)=1\), \(y(1)=e\) using finite differences with \(h=0.5\).
Grid: \(x_0=0,x_1=0.5,x_2=1\). Replace \(y''\approx(y_{i+1}-2y_i+y_{i-1})/h^2\). At \(x_1\): \(\frac{y_0-2y_1+y_2}{0.25}-y_1=0\Rightarrow4y_0-8y_1+4y_2-y_1=0\Rightarrow-9y_1=-(4y_0+4y_2)=-(4+4e)\Rightarrow y_1=\frac{4(1+e)}{9}\approx\frac{4\times3.7183}{9}=1.6526\). Exact: \(y(0.5)=e^{0.5}\approx1.6487\). Error \(\approx0.004\). \(\blacksquare\)
Solve \(10x+y+z=12\), \(x+10y+z=12\), \(x+y+10z=12\) using Gauss-Seidel. Start from \((0,0,0)\).
Dominant diagonal \((10\gg1+1)\): Gauss-Seidel converges. Iteration formulas: \(x=(12-y-z)/10,y=(12-x-z)/10,z=(12-x-y)/10\). Iter 1: \(x=12/10=1.2\), \(y=(12-1.2)/10=1.08\), \(z=(12-1.2-1.08)/10=0.972\). Iter 2: \(x=(12-1.08-0.972)/10=0.9948\), \(y=(12-0.9948-0.972)/10=1.0033\), \(z=(12-0.9948-1.0033)/10=1.0002\). Converging to \((1,1,1)\). \(\blacksquare\)
Factorize \(A=\begin{pmatrix}2&1&1\\4&3&3\\8&7&9\end{pmatrix}\) as \(LU\) and solve \(A\mathbf{x}=(4,10,24)^T\).
Gaussian elimination: \(R_2\to R_2-2R_1,R_3\to R_3-4R_1\): \(\begin{pmatrix}2&1&1\\0&1&1\\0&3&5\end{pmatrix}\). \(R_3\to R_3-3R_2\): \(\begin{pmatrix}2&1&1\\0&1&1\\0&0&2\end{pmatrix}=U\). \(L=\begin{pmatrix}1&0&0\\2&1&0\\4&3&1\end{pmatrix}\). Solve \(L\mathbf{y}=\mathbf{b}\): \(y_1=4,y_2=10-8=2,y_3=24-16-6=2\). Solve \(U\mathbf{x}=\mathbf{y}\): \(x_3=1,x_2=2-1=1,x_1=(4-1-1)/2=1\). \(\mathbf{x}=(1,1,1)^T\). \(\blacksquare\)
Find the largest eigenvalue of \(A=\begin{pmatrix}2&1\\1&2\end{pmatrix}\) by power iteration starting from \((1,1)^T\).
Iter 1: \(A\begin{pmatrix}1\\1\end{pmatrix}=\begin{pmatrix}3\\3\end{pmatrix}=3\begin{pmatrix}1\\1\end{pmatrix}\). Already eigenvector! \(\lambda_{\max}=3\), eigenvector \((1,1)^T\). Verify: eigenvalues of \(\begin{pmatrix}2&1\\1&2\end{pmatrix}\) are \(2\pm1=3,1\). Correct. \(\blacksquare\)
Given \(f(0)=1,f(0.5)=1.6487,f(1)=2.7183\), find \(f'(0)\) and \(f''(0.5)\).
\(h=0.5\). \(f'(0)\approx\frac{-3f_0+4f_1-f_2}{2h}=\frac{-3(1)+4(1.6487)-2.7183}{1}=\frac{-3+6.5948-2.7183}{1}=0.8765\). (Exact: \(e^0=1\); error due to coarse step.) \(f''(0.5)\approx\frac{f_0-2f_1+f_2}{h^2}=\frac{1-3.2974+2.7183}{0.25}=\frac{0.4209}{0.25}=1.6836\). Exact: \(e^{0.5}=1.6487\). \(\blacksquare\)
Define condition number of a matrix. When is a linear system ill-conditioned?
Condition number: \(\kappa(A)=\|A\|\cdot\|A^{-1}\|\). For SPD: \(\kappa=\lambda_{\max}/\lambda_{\min}\). Ill-conditioned: \(\kappa\gg1\). Significance: relative error in solution \(\|\delta\mathbf{x}\|/\|\mathbf{x}\|\le\kappa(A)\cdot\|\delta\mathbf{b}\|/\|\mathbf{b}\|\). Large \(\kappa\): small perturbations in \(\mathbf{b}\) cause large errors in \(\mathbf{x}\). Example: Hilbert matrix \(H_{ij}=1/(i+j-1)\) is highly ill-conditioned. Remedy: scaling, pivoting, iterative refinement. \(\blacksquare\)
Distinguish between truncation error and round-off error in numerical computation.
Truncation error: error from approximating mathematical operations (e.g., replacing infinite Taylor series by finite terms, replacing derivatives by finite differences). Reduces with smaller step \(h\). Round-off error: error from representing real numbers in finite precision (floating point). Grows as \(h\to0\) (catastrophic cancellation). Trade-off: optimal \(h\) balances both. For numerical differentiation: truncation \(\sim h\), round-off \(\sim\epsilon_M/h\); optimal \(h\sim\sqrt{\epsilon_M}\approx10^{-8}\) for double precision. \(\blacksquare\)
What are Chebyshev nodes? Why are they preferred over equally-spaced nodes?
Chebyshev nodes on \([-1,1]\): \(x_k=\cos\!\left(\frac{2k+1}{2n+2}\pi\right),k=0,\ldots,n\). Clustered at endpoints. Minimize the maximum of \(|\prod(x-x_k)|\) over \([-1,1]\): \(\max|\omega(x)|=\frac{1}{2^n}\) (minimal). Equally-spaced: \(\max|\omega|\) grows exponentially with \(n\) → Runge phenomenon. Chebyshev nodes: near-minimax interpolation error, no Runge phenomenon. Use: FFT, spectral methods, near-best polynomial approximation. \(\blacksquare\)
Why is partial pivoting important in Gaussian elimination? Illustrate with an example.
Without pivoting, small pivot causes division by small number → large round-off errors. Partial pivoting: at each step, swap rows to bring largest element in current column to pivot position. Example: \(\begin{pmatrix}0.001&1\\1&1\end{pmatrix}\begin{pmatrix}x\\y\end{pmatrix}=\begin{pmatrix}1\\2\end{pmatrix}\). Without pivoting: \(x=(1000-2)/(-999)\approx-0.998\), large error. With pivoting: swap rows, eliminate cleanly. Result: \(x=y=1\) (exact for this system). Growth factor bounded to \(2^{n-1}\) with partial pivoting. \(\blacksquare\)
Apply the Gauss-Seidel method to solve the system \(10x_1 - x_2 + 2x_3 = 6\), \(-x_1 + 11x_2 - x_3 + 3x_4 = 25\), \(2x_1 - x_2 + 10x_3 - x_4 = -11\) (ignoring \(x_4\) for a \(3\times3\) subsystem with the first three equations and \(x_4=0\)) starting from \(\mathbf{x}^{(0)}=(0,0,0)\). Perform 3 iterations.
System (3×3): \(10x_1 - x_2 + 2x_3 = 6\), \(-x_1 + 11x_2 - x_3 = 25\), \(2x_1 - x_2 + 10x_3 = -11\). Diagonal dominance holds: \(|10|>|-1|+|2|\), etc. Gauss-Seidel formulae (use latest values immediately): \[x_1 = \tfrac{1}{10}(6 + x_2 - 2x_3),\quad x_2 = \tfrac{1}{11}(25 + x_1 + x_3),\quad x_3 = \tfrac{1}{10}(-11 - 2x_1 + x_2).\] Iter 1: \(x_1=6/10=0.6\), \(x_2=(25+0.6)/11=2.3273\), \(x_3=(-11-1.2+2.3273)/10=-0.9873\). Iter 2: \(x_1=(6+2.3273+1.9746)/10=1.0302\), \(x_2=(25+1.0302-0.9873)/11=2.2766\), \(x_3=(-11-2.0604+2.2766)/10=-1.0784\). Iter 3: \(x_1=(6+2.2766+2.1568)/10=1.0433\), \(x_2=(25+1.0433-1.0784)/11=2.2695\), \(x_3=(-11-2.0866+2.2695)/10=-1.0817\). Converging to exact \((1, 2, -1)\). \(\blacksquare\)
Compare the convergence rates of Jacobi and Gauss-Seidel iteration for solving \(Ax=b\). Under what conditions does each converge?
Jacobi: \(x_i^{(k+1)} = \frac{1}{a_{ii}}\!\left(b_i - \sum_{j\neq i}a_{ij}x_j^{(k)}\right)\). Uses only old values. Iteration matrix \(B_J = -D^{-1}(L+U)\). Convergence iff \(\rho(B_J)<1\). Gauss-Seidel: uses updated values immediately. Iteration matrix \(B_{GS}=-(D+L)^{-1}U\). For many classes of matrices \(\rho(B_{GS})=[\rho(B_J)]^2\) (Stein-Rosenberg theorem), so GS converges about twice as fast per iteration. Both converge for strictly diagonally dominant or symmetric positive definite \(A\). GS converges whenever Jacobi does (for SPD and irreducibly diagonally dominant matrices). Counter-example: some matrices where Jacobi converges but GS diverges exist, but they are non-SPD. In practice, GS is preferred for speed. \(\blacksquare\)
Find the LU decomposition (without pivoting) of \(A=\begin{pmatrix}2&1&1\\4&3&3\\8&7&9\end{pmatrix}\). Express \(A=LU\) where \(L\) is unit lower triangular.
Step 1: Eliminate column 1. \(m_{21}=4/2=2\), \(m_{31}=8/2=4\). \(R_2\leftarrow R_2-2R_1=(0,1,1)\), \(R_3\leftarrow R_3-4R_1=(0,3,5)\). Step 2: Eliminate column 2. \(m_{32}=3/1=3\). \(R_3\leftarrow R_3-3R_2=(0,0,2)\). Result: \[L=\begin{pmatrix}1&0&0\\2&1&0\\4&3&1\end{pmatrix},\quad U=\begin{pmatrix}2&1&1\\0&1&1\\0&0&2\end{pmatrix}.\] Verify: \(LU=\begin{pmatrix}2&1&1\\4&3&3\\8&7&9\end{pmatrix}=A\). \(\checkmark\) \(\blacksquare\)
Compute the Cholesky decomposition \(A=LL^T\) for the symmetric positive definite matrix \(A=\begin{pmatrix}4&2&-2\\2&2&0\\-2&0&3\end{pmatrix}\).
Cholesky: \(l_{ij}=\frac{1}{l_{jj}}\!\left(a_{ij}-\sum_{k=1}^{j-1}l_{ik}l_{jk}\right)\), \(l_{jj}=\sqrt{a_{jj}-\sum_{k=1}^{j-1}l_{jk}^2}\). \(l_{11}=\sqrt{4}=2\). \(l_{21}=2/2=1\), \(l_{31}=-2/2=-1\). \(l_{22}=\sqrt{2-1^2}=1\). \(l_{32}=(0-(-1)(1))/1=1\). \(l_{33}=\sqrt{3-1-1}=1\). \[L=\begin{pmatrix}2&0&0\\1&1&0\\-1&1&1\end{pmatrix}.\] Verify \(LL^T\): \((1,1)\!: 4\). \((2,1)\!: 2\). \((3,1)\!: -2\). \((2,2)\!: 2\). \((3,2)\!: 0\). \((3,3)\!: 3\). \(\checkmark\) \(\blacksquare\)
Compute the QR factorization of \(A=\begin{pmatrix}1&1\\1&0\\0&1\end{pmatrix}\) using the classical Gram-Schmidt process.
Columns: \(\mathbf{a}_1=(1,1,0)^T\), \(\mathbf{a}_2=(1,0,1)^T\). \(\mathbf{q}_1=\mathbf{a}_1/\|\mathbf{a}_1\|=(1,1,0)^T/\sqrt{2}\). \(r_{12}=\mathbf{q}_1^T\mathbf{a}_2=(1+0)/\sqrt{2}=1/\sqrt{2}\). \(\mathbf{v}_2=\mathbf{a}_2-r_{12}\mathbf{q}_1=(1,0,1)^T-(1/2,1/2,0)^T=(1/2,-1/2,1)^T\). \(r_{22}=\|\mathbf{v}_2\|=\sqrt{1/4+1/4+1}=\sqrt{3/2}\). \(\mathbf{q}_2=(1/2,-1/2,1)^T/\sqrt{3/2}=(1/\sqrt{6},-1/\sqrt{6},2/\sqrt{6})^T\). \[Q=\begin{pmatrix}1/\sqrt2&1/\sqrt6\\1/\sqrt2&-1/\sqrt6\\0&2/\sqrt6\end{pmatrix},\quad R=\begin{pmatrix}\sqrt2&1/\sqrt2\\0&\sqrt{3/2}\end{pmatrix}.\] \(\blacksquare\)
Apply the power method to find the largest eigenvalue of \(A=\begin{pmatrix}4&1\\2&3\end{pmatrix}\), starting from \(\mathbf{x}^{(0)}=(1,1)^T\). Perform 3 iterations.
Power method: \(\mathbf{y}^{(k)}=A\mathbf{x}^{(k-1)}\), \(\mu^{(k)}=\max|\mathbf{y}^{(k)}|\), \(\mathbf{x}^{(k)}=\mathbf{y}^{(k)}/\mu^{(k)}\). Iter 1: \(\mathbf{y}=(5,5)^T\), \(\mu=5\), \(\mathbf{x}=(1,1)^T\). Iter 2: \(\mathbf{y}=(5,5)^T\), \(\mu=5\), \(\mathbf{x}=(1,1)^T\). The vector is already the eigenvector! (Since \((1,1)^T\) is the eigenvector for \(\lambda=5\): \(A(1,1)^T=(5,5)^T=5(1,1)^T\).) Exact eigenvalues: \(\det(A-\lambda I)=(\lambda-4)(\lambda-3)-2=(λ-5)(λ-2)=0\), so \(\lambda_1=5, \lambda_2=2\). Power method converges in 1 step here. \(\blacksquare\)
Describe the inverse power method for finding the smallest eigenvalue of a matrix \(A\). How is it implemented?
Smallest eigenvalue of \(A\) = largest eigenvalue of \(A^{-1}\). Algorithm: given \(\mathbf{x}^{(0)}\), at each step solve \(A\mathbf{y}^{(k)}=\mathbf{x}^{(k-1)}\) (don't actually invert; use LU decomposition), set \(\mu^{(k)}=\max|\mathbf{y}^{(k)}|\), \(\mathbf{x}^{(k)}=\mathbf{y}^{(k)}/\mu^{(k)}\). Converges to eigenvector of smallest \(|\lambda|\), with eigenvalue \(\approx1/\mu^{(k)}\). Shifted inverse: to find eigenvalue nearest to \(\sigma\), apply inverse power to \((A-\sigma I)\). Use: deflation after finding dominant eigenvalue, computing interior eigenvalues. \(\blacksquare\)
Explain the QR algorithm for computing all eigenvalues of a real symmetric matrix. State its convergence properties.
QR algorithm: set \(A_0=A\). At each step, factor \(A_k=Q_kR_k\) (QR decomposition), then set \(A_{k+1}=R_kQ_k\). The sequence \(\{A_k\}\) is orthogonally similar to \(A\) (same eigenvalues). For symmetric \(A\), \(A_k\to\) diagonal matrix with eigenvalues on diagonal. Convergence: off-diagonal entries decay; rate depends on ratio of consecutive eigenvalues \(|\lambda_{i+1}/\lambda_i|\). Practical improvements: (1) reduce to Hessenberg form first (reduces QR step cost from \(O(n^3)\) to \(O(n^2)\)); (2) use shifts \((A_k-\sigma_k I)=Q_kR_k\) for faster convergence; (3) Wilkinson shift for global convergence. \(\blacksquare\)
State and prove the Gershgorin circle theorem. Use it to bound the eigenvalues of \(A=\begin{pmatrix}5&1&0\\1&4&1\\0&1&3\end{pmatrix}\).
Theorem: Every eigenvalue of \(A\in\mathbb{C}^{n\times n}\) lies in at least one Gershgorin disk \(D_i=\{z:|z-a_{ii}|\le R_i\}\) where \(R_i=\sum_{j\neq i}|a_{ij}|\). Proof: Let \(\lambda,\mathbf{x}\) be an eigenpair; let \(x_i=\max_j|x_j|>0\). Row \(i\) of \(A\mathbf{x}=\lambda\mathbf{x}\): \(\sum_j a_{ij}x_j=\lambda x_i\), so \((\lambda-a_{ii})x_i=\sum_{j\neq i}a_{ij}x_j\). Taking absolute values: \(|\lambda-a_{ii}|\le\sum_{j\neq i}|a_{ij}||x_j|/|x_i|\le R_i\). \(\square\) Application: \(R_1=1,R_2=2,R_3=1\). Disks: \(|z-5|\le1\), \(|z-4|\le2\), \(|z-3|\le1\). All eigenvalues in \([2,6]\). Since \(A\) is SPD, all positive. \(\blacksquare\)
Define the condition number \(\kappa(A)=\|A\|\cdot\|A^{-1}\|\). Give an example of an ill-conditioned matrix and explain its practical significance.
Condition number \(\kappa(A)=\|A\|_2\cdot\|A^{-1}\|_2=\sigma_{\max}/\sigma_{\min}\) (ratio of largest to smallest singular values). For SPD: \(\kappa=\lambda_{\max}/\lambda_{\min}\). Perturbation bound: \(\frac{\|\delta\mathbf{x}\|}{\|\mathbf{x}\|}\le\kappa(A)\frac{\|\delta\mathbf{b}\|}{\|\mathbf{b}\|}\). Thus relative error in \(\mathbf{x}\) can be \(\kappa\) times larger than relative error in \(\mathbf{b}\). Ill-conditioned example: Hilbert matrix \(H_{ij}=1/(i+j-1)\). For \(n=5\), \(\kappa(H_5)\approx4.77\times10^5\). A residual of \(10^{-10}\) may give a solution error of \(10^{-5}\). Remedy: preconditioning, scaling, higher precision arithmetic. \(\blacksquare\)
Distinguish between partial and complete pivoting in Gaussian elimination. Analyse their effect on numerical stability.
Partial pivoting: at step \(k\), swap rows to bring the largest element in column \(k\) (below diagonal) to pivot position. Cost: \(O(n)\) comparisons per step. Complete pivoting: swap both rows and columns to bring the largest element in the remaining submatrix to pivot. Cost: \(O(n^2)\) comparisons per step. Stability: partial pivoting: growth factor \(g\le 2^{n-1}\) (theoretically), but in practice rarely large; standard and sufficient. Complete pivoting: growth factor \(g\le(n\prod_{k=1}^n k^{1/(k-1)})^{1/2}\), empirically much smaller. Complete pivoting is more stable but expensive; partial pivoting is the standard practice. LU with partial pivoting is backward stable: computed solution \(\hat{\mathbf{x}}\) satisfies \((A+\delta A)\hat{\mathbf{x}}=\mathbf{b}\) with \(\|\delta A\|\le C u\|A\|\) where \(u\) = machine epsilon. \(\blacksquare\)
Describe the Thomas algorithm (tridiagonal matrix algorithm) for solving a tridiagonal system \(a_i x_{i-1}+b_i x_i+c_i x_{i+1}=d_i\). State its complexity.
Tridiagonal system: \(\begin{pmatrix}b_1&c_1\\a_2&b_2&c_2\\&\ddots&\ddots&\ddots\\&&a_n&b_n\end{pmatrix}\mathbf{x}=\mathbf{d}\). Forward sweep (LU): \(w_1=b_1\); for \(i=2,\ldots,n\): \(m_i=a_i/w_{i-1}\), \(w_i=b_i-m_i c_{i-1}\), \(d'_i=d_i-m_i d'_{i-1}\). Back substitution: \(x_n=d'_n/w_n\); for \(i=n-1,\ldots,1\): \(x_i=(d'_i-c_i x_{i+1})/w_i\). Complexity: \(O(n)\) operations (only \(5n-4\) multiplications/additions). Numerically stable without pivoting if the matrix is diagonally dominant or SPD. Applications: cubic spline interpolation, Crank-Nicolson for PDEs. \(\blacksquare\)
Describe the SOR (Successive Over-Relaxation) method. Derive the optimal relaxation parameter \(\omega_{opt}\) for the model problem.
SOR: \(x_i^{(k+1)} = (1-\omega)x_i^{(k)} + \frac{\omega}{a_{ii}}\!\left(b_i - \sum_{j
Describe the conjugate gradient (CG) method for solving \(Ax=b\) with \(A\) symmetric positive definite. State its convergence theorem.
CG minimizes \(\phi(\mathbf{x})=\frac12\mathbf{x}^TA\mathbf{x}-\mathbf{b}^T\mathbf{x}\) over Krylov subspace \(\mathcal{K}_k=\text{span}\{\mathbf{r}_0,A\mathbf{r}_0,\ldots,A^{k-1}\mathbf{r}_0\}\). Algorithm: \(\mathbf{r}_0=\mathbf{b}-A\mathbf{x}_0\), \(\mathbf{p}_0=\mathbf{r}_0\). At each step: \(\alpha_k=\mathbf{r}_k^T\mathbf{r}_k/\mathbf{p}_k^TA\mathbf{p}_k\), \(\mathbf{x}_{k+1}=\mathbf{x}_k+\alpha_k\mathbf{p}_k\), \(\mathbf{r}_{k+1}=\mathbf{r}_k-\alpha_k A\mathbf{p}_k\), \(\beta_k=\mathbf{r}_{k+1}^T\mathbf{r}_{k+1}/\mathbf{r}_k^T\mathbf{r}_k\), \(\mathbf{p}_{k+1}=\mathbf{r}_{k+1}+\beta_k\mathbf{p}_k\). Convergence: \(\|\mathbf{e}_k\|_A\le2\left(\frac{\sqrt\kappa-1}{\sqrt\kappa+1}\right)^k\|\mathbf{e}_0\|_A\). Terminates in at most \(n\) steps exactly (finite precision). Effective with preconditioning. \(\blacksquare\)
Describe the CSR (Compressed Sparse Row) and CSC (Compressed Sparse Column) formats for storing sparse matrices. Give an example.
For \(A=\begin{pmatrix}5&0&1\\0&3&0\\2&0&4\end{pmatrix}\), non-zeros: 5,1,3,2,4. CSR: three arrays: \texttt{val}=[5,1,3,2,4], \texttt{col\_ind}=[0,2,1,0,2], \texttt{row\_ptr}=[0,2,3,5]. Row \(i\) has non-zeros at \texttt{val[row\_ptr[i]:row\_ptr[i+1]]\). Space: \(O(nnz+n)\). Efficient for row operations, matrix-vector products. CSC: transpose view: \texttt{val}=[5,2,3,1,4], \texttt{row\_ind}=[0,2,1,0,2], \texttt{col\_ptr}=[0,2,3,5]. Efficient for column operations. Applications: finite element, graph algorithms, iterative solvers. \(\blacksquare\)
State the backward error result for Gaussian elimination with partial pivoting. What does it imply about the computed solution?
Theorem (Wilkinson): Gaussian elimination with partial pivoting applied to \(A\mathbf{x}=\mathbf{b}\) produces a computed solution \(\hat{\mathbf{x}}\) that satisfies \((A+\delta A)\hat{\mathbf{x}}=\mathbf{b}\) where \(\|\delta A\|_\infty\le g\cdot n^3 u\|A\|_\infty\) and \(g\) is the growth factor, \(u=\epsilon_M/2\) (machine unit). With partial pivoting, \(g\le 2^{n-1}\) in theory but rarely large in practice. Implication: the computed solution is the exact solution to a nearby problem; the error \(\|\hat{\mathbf{x}}-\mathbf{x}\|/\|\mathbf{x}\|\lesssim\kappa(A)\cdot n^3 u\). The method is backward stable. \(\blacksquare\)
Prove the perturbation bound \(\|\delta\mathbf{x}\|/\|\mathbf{x}\|\le\kappa(A)\|\delta\mathbf{b}\|/\|\mathbf{b}\|\) for \(A\mathbf{x}=\mathbf{b}\), \(A(\mathbf{x}+\delta\mathbf{x})=\mathbf{b}+\delta\mathbf{b}\).
From \(A\delta\mathbf{x}=\delta\mathbf{b}\): \(\delta\mathbf{x}=A^{-1}\delta\mathbf{b}\), so \(\|\delta\mathbf{x}\|\le\|A^{-1}\|\|\delta\mathbf{b}\|\). From \(A\mathbf{x}=\mathbf{b}\): \(\|\mathbf{b}\|\le\|A\|\|\mathbf{x}\|\Rightarrow 1/\|\mathbf{x}\|\le\|A\|/\|\mathbf{b}\|\). Combining: \[\frac{\|\delta\mathbf{x}\|}{\|\mathbf{x}\|}\le\|A^{-1}\|\|\delta\mathbf{b}\|\cdot\frac{\|A\|}{\|\mathbf{b}\|}=\kappa(A)\frac{\|\delta\mathbf{b}\|}{\|\mathbf{b}\|}.\] Similarly for perturbation of \(A\): \(\|\delta\mathbf{x}\|/\|\mathbf{x}\|\le\kappa(A)\|\delta A\|/\|A\|\). These bounds are sharp (achievable). \(\blacksquare\)
Define the 1-norm, ∞-norm, 2-norm (spectral norm), and Frobenius norm of a matrix. Compute each for \(A=\begin{pmatrix}1&2\\-1&3\end{pmatrix}\).
\(\|A\|_1=\max_j\sum_i|a_{ij}|\) (max col sum) \(=\max(2,5)=5\). \(\|A\|_\infty=\max_i\sum_j|a_{ij}|\) (max row sum) \(=\max(3,4)=4\). \(\|A\|_2=\sigma_{\max}(A)=\sqrt{\lambda_{\max}(A^TA)}\). \(A^TA=\begin{pmatrix}2&-1\\-1&13\end{pmatrix}\). Char poly: \(\lambda^2-15\lambda+25=0\Rightarrow\lambda=(15\pm\sqrt{125})/2\). \(\|A\|_2=\sqrt{(15+5\sqrt5)/2}\approx3.86\). \(\|A\|_F=\sqrt{\sum|a_{ij}|^2}=\sqrt{1+4+1+9}=\sqrt{15}\approx3.87\). Note \(\|A\|_2\le\|A\|_F\le\sqrt{n}\|A\|_2\). \(\blacksquare\)
Describe the iterative refinement procedure for improving the accuracy of the solution to \(Ax=b\) computed by Gaussian elimination. When is it effective?
Given approximate solution \(\hat{\mathbf{x}}_0\) (from LU decomposition of \(A\)): 1. Compute residual \(\mathbf{r}_k=\mathbf{b}-A\hat{\mathbf{x}}_k\) in higher precision arithmetic. 2. Solve \(A\mathbf{d}_k=\mathbf{r}_k\) using the stored LU factorization: \(O(n^2)\) work. 3. Update: \(\hat{\mathbf{x}}_{k+1}=\hat{\mathbf{x}}_k+\mathbf{d}_k\). Each step multiplies the error by \(\sim u\kappa(A)\). Effective when \(u\kappa(A)<1\) (i.e., matrix not too ill-conditioned relative to machine precision); converges to full machine precision in a few steps. Requires residual computed in extra precision (e.g., 80-bit extended for 64-bit working). \(\blacksquare\)
Derive the normal equations \(A^TA\mathbf{x}=A^T\mathbf{b}\) for the least squares problem \(\min_{\mathbf{x}}\|A\mathbf{x}-\mathbf{b}\|_2\). When is the solution unique?
Minimize \(\phi(\mathbf{x})=\|A\mathbf{x}-\mathbf{b}\|^2=\mathbf{x}^TA^TA\mathbf{x}-2\mathbf{b}^TA\mathbf{x}+\mathbf{b}^T\mathbf{b}\). Setting gradient to zero: \(\nabla\phi=2A^TA\mathbf{x}-2A^T\mathbf{b}=0\Rightarrow A^TA\mathbf{x}=A^T\mathbf{b}\). Geometric interpretation: \(A\mathbf{x}-\mathbf{b}\perp\text{col}(A)\), i.e., residual is orthogonal to column space. Uniqueness: solution unique iff \(A^TA\) is nonsingular iff \(A\) has full column rank. Numerical caution: \(\kappa(A^TA)=[\kappa(A)]^2\); normal equations double the condition number. Prefer QR decomposition: \(A=QR\), solve \(R\mathbf{x}=Q^T\mathbf{b}\). \(\blacksquare\)
Explain Richardson extrapolation as applied to numerical differentiation. Derive the improved formula for \(f'(x)\) using step sizes \(h\) and \(h/2\).
Central difference: \(D(h)=\frac{f(x+h)-f(x-h)}{2h}=f'(x)+\frac{h^2}{6}f'''(x)+O(h^4)\). Error is \(O(h^2)\). With \(h/2\): \(D(h/2)=f'(x)+\frac{h^2}{24}f'''(x)+O(h^4)\). Eliminate leading error term: \(4D(h/2)-D(h)=4f'(x)+\frac{h^2}{6}f'''(x)/4\cdot4-[f'(x)+\frac{h^2}{6}f'''(x)]+O(h^4)\). \[\Rightarrow\quad f'(x)\approx\frac{4D(h/2)-D(h)}{3}+O(h^4).\] Error is now \(O(h^4)\) — two orders better than naive \(O(h^2)\). Richardson extrapolation can be iterated (Romberg-like) to achieve any order. \(\blacksquare\)
Write down the forward, backward, and central difference approximations for \(f'(x)\). State the order of error for each and explain the trade-off in choosing step size \(h\).
Forward: \(f'(x)\approx\frac{f(x+h)-f(x)}{h}\), error \(O(h)\). Backward: \(f'(x)\approx\frac{f(x)-f(x-h)}{h}\), error \(O(h)\). Central: \(f'(x)\approx\frac{f(x+h)-f(x-h)}{2h}\), error \(O(h^2)\). Derivation (Taylor): \(f(x\pm h)=f(x)\pm hf'(x)+\frac{h^2}{2}f''(x)\pm\frac{h^3}{6}f'''(x)+\cdots\). Central difference cancels \(O(h)\) term. Trade-off: reducing \(h\) reduces truncation error (\(\sim h^2\)) but increases round-off error (\(\sim\epsilon_M/h\)). Optimal \(h_{opt}\approx(3\epsilon_M/|f'''|)^{1/3}\) for central difference, giving minimum total error \(\sim\epsilon_M^{2/3}\). \(\blacksquare\)
Derive the central difference formula \(f''(x)\approx\frac{f(x+h)-2f(x)+f(x-h)}{h^2}\) and state its error. Apply it to estimate \(f''(1)\) for \(f(x)=\ln x\) with \(h=0.1\).
Derivation: \(f(x+h)=f+hf'+\frac{h^2}{2}f''+\frac{h^3}{6}f'''+\frac{h^4}{24}f^{(4)}+\cdots\) \(f(x-h)=f-hf'+\frac{h^2}{2}f''-\frac{h^3}{6}f'''+\frac{h^4}{24}f^{(4)}+\cdots\) Adding: \(f(x+h)+f(x-h)=2f+h^2f''+\frac{h^4}{12}f^{(4)}+\cdots\) \[\Rightarrow f''(x)=\frac{f(x+h)-2f(x)+f(x-h)}{h^2}+O(h^2).\] Application: \(f(x)=\ln x\), \(f''(x)=-1/x^2\). \(f(1.1)=0.09531\), \(f(1)=0\), \(f(0.9)=-0.10536\). Formula: \((0.09531-0+(-0.10536))/0.01=-0.01005/0.01=-1.005/0.01\). Wait: \(\frac{0.09531-0-0.10536}{0.01}=\frac{-0.01005}{0.01}=-1.005\). Exact: \(-1/1=-1\). Error \(\approx0.005=O(h^2)\). \(\blacksquare\)
Given \(f(0)=1,f(0.1)=1.1052,f(0.2)=1.2214,f(0.3)=1.3499\), find \(f'(0.5)\) using Newton's forward difference formula (note: use these as forward differences for \(f'(x)\) at \(x=x_0\) shifted).
Build difference table (\(h=0.1\)): \(\Delta f_0=0.1052,\Delta f_1=0.1162,\Delta f_2=0.1285\). \(\Delta^2 f_0=0.0110,\Delta^2 f_1=0.0123\). \(\Delta^3 f_0=0.0013\). Newton forward differentiation formula at \(x_0\) (\(s=0\)): \[f'(x_0)=\frac{1}{h}\left[\Delta f_0-\frac{1}{2}\Delta^2 f_0+\frac{1}{3}\Delta^3 f_0-\cdots\right]\] \[=\frac{1}{0.1}\left[0.1052-0.0055+0.000433-\cdots\right]=\frac{0.09986}{0.1}\approx0.9986.\] Exact \(f'(x)=e^x\), \(f'(0)=1\). The question asks \(f'(0.5)\): extrapolation would require more data; using the given table at \(x_0=0\), \(f'(0)\approx0.9986\) (good approximation). \(\blacksquare\)
Define absolute stability for numerical ODE methods. Find the region of absolute stability for Euler's method applied to \(y'=\lambda y\).
Test equation: \(y'=\lambda y\), \(\lambda\in\mathbb{C}\), \(\text{Re}(\lambda)<0\). Euler: \(y_{n+1}=y_n+h\lambda y_n=(1+h\lambda)y_n\). Stability requires \(|1+h\lambda|<1\) (so errors don't grow). Set \(z=h\lambda\): region is the disk \(|1+z|<1\), i.e., the disk of radius 1 centered at \(-1\) in the \(z\)-plane. For real \(\lambda<0\): need \(|1+h\lambda|<1\Rightarrow h|\lambda|<2\Rightarrow h<2/|\lambda|\). Stiff equations (large \(|\lambda|\)) require tiny \(h\) for stability — very inefficient. Implicit Euler has stability region \(|z+1|>1\) complemented, i.e., all \(\text{Re}(z)<0\) (A-stable). \(\blacksquare\)
Define a stiff ODE system. Why do explicit methods fail? Show that the backward Euler method is A-stable.
Stiff system: \(\mathbf{y}'=A\mathbf{y}\) where eigenvalues of \(A\) have widely different real parts, some very negative. Stiffness ratio: \(S=\max|\text{Re}(\lambda_i)|/\min|\text{Re}(\lambda_i)|\gg1\). Explicit failure: stability requires \(h<2/\max|\text{Re}(\lambda_i)|\), which may be astronomically small even though the "interesting" dynamics occur on a much slower time scale. Backward (implicit) Euler: \(y_{n+1}=y_n+hf(t_{n+1},y_{n+1})\). For test equation: \(y_{n+1}=y_n/(1-h\lambda)\). Stability: \(|y_{n+1}/y_n|=1/|1-h\lambda|<1\) whenever \(\text{Re}(h\lambda)<0\), i.e., for all \(z=h\lambda\) with \(\text{Re}(z)<0\). This is exactly A-stability. \(\blacksquare\)
Derive the Adams-Bashforth two-step formula \(y_{n+1}=y_n+\frac{h}{2}(3f_n-f_{n-1})\). What is its local truncation error?
Adams methods integrate the interpolating polynomial of \(f\). For two-step, fit polynomial \(P(t)\) through \((t_{n-1},f_{n-1}),(t_n,f_n)\) — linear: \(P(t)=f_n+\frac{f_n-f_{n-1}}{h}(t-t_n)\). Integrate from \(t_n\) to \(t_{n+1}\): \[\int_{t_n}^{t_{n+1}}P(t)\,dt = hf_n + \frac{h}{2}(f_n-f_{n-1}) = \frac{h}{2}(3f_n-f_{n-1}).\] Hence \(y_{n+1}=y_n+\frac{h}{2}(3f_n-f_{n-1})\). Local truncation error: expand using Taylor series: LTE \(=\frac{5}{12}h^3y'''(\xi)\) for some \(\xi\in(t_{n-1},t_{n+1})\). Order 2, explicit. Requires a starting value \(y_1\) from a one-step method. \(\blacksquare\)
Describe the Adams-Bashforth–Adams-Moulton predictor-corrector pair. What are the advantages over using only the Adams-Moulton corrector?
Predictor (AB4): \(y_{n+1}^P=y_n+\frac{h}{24}(55f_n-59f_{n-1}+37f_{n-2}-9f_{n-3})\), LTE \(=\frac{251}{720}h^5y^{(5)}\). Corrector (AM4): \(y_{n+1}^C=y_n+\frac{h}{24}(9f(t_{n+1},y_{n+1}^P)+19f_n-5f_{n-1}+f_{n-2})\), LTE \(=-\frac{19}{720}h^5y^{(5)}\). Advantage: Adams-Moulton corrector is implicit (solve nonlinear equation for \(y_{n+1}\)), but using the predictor as starting value for the implicit solve makes one corrector step sufficient (PECE mode: Predict-Evaluate-Correct-Evaluate). Cost: 2 function evaluations per step vs 4 for RK4, but higher accuracy (order 4 with error estimate). Also: error estimate \(\approx-\frac{19}{270}(y_{n+1}^C-y_{n+1}^P)\) for step-size control. \(\blacksquare\)
Use Milne's predictor-corrector method to solve \(y'=x+y\), \(y(0)=1\) with \(h=0.1\). Starting values: \(y(0)=1\), \(y(0.1)=1.1103\), \(y(0.2)=1.2428\), \(y(0.3)=1.3997\). Compute \(y(0.4)\).
\(f(x,y)=x+y\). Compute: \(f_0=0+1=1\), \(f_1=0.1+1.1103=1.2103\), \(f_2=0.2+1.2428=1.4428\), \(f_3=0.3+1.3997=1.6997\). Milne predictor: \(y_4^P=y_0+\frac{4h}{3}(2f_1-f_2+2f_3)=1+\frac{0.4}{3}(2(1.2103)-1.4428+2(1.6997))\) \(=1+\frac{0.4}{3}(2.4206-1.4428+3.3994)=1+\frac{0.4}{3}(4.3772)=1+0.5836=1.5836\). \(f_4^P=0.4+1.5836=1.9836\). Simpson corrector: \(y_4^C=y_2+\frac{h}{3}(f_2+4f_3+f_4^P)=1.2428+\frac{0.1}{3}(1.4428+4(1.6997)+1.9836)\) \(=1.2428+\frac{0.1}{3}(1.4428+6.7988+1.9836)=1.2428+\frac{0.1}{3}(10.2252)=1.2428+0.3408=1.5836\). \(\blacksquare\)
Describe Heun's method for solving \(y'=f(x,y)\). Apply it to solve \(y'=-y\), \(y(0)=1\) for one step with \(h=0.5\). Compare with exact.
Heun's method (second-order Runge-Kutta): \(k_1=f(x_n,y_n)\), \(k_2=f(x_n+h,y_n+hk_1)\), \(y_{n+1}=y_n+\frac{h}{2}(k_1+k_2)\). Application: \(f=-y\), \(y_0=1\), \(h=0.5\). \(k_1=f(0,1)=-1\). \(\tilde{y}_1=1+0.5(-1)=0.5\). \(k_2=f(0.5,0.5)=-0.5\). \(y_1=1+\frac{0.5}{2}(-1+(-0.5))=1+0.25(-1.5)=1-0.375=0.625\). Exact: \(e^{-0.5}=0.6065\). Error: \(0.625-0.6065=0.0185\). Order 2 local error (\(\sim h^3\)). vs Euler: \(y_1=0.5\), error \(0.1065\). Heun's clearly better. \(\blacksquare\)
State the classical Runge-Kutta 4th-order method for the system \(y'=f(t,y,z)\), \(z'=g(t,y,z)\). Describe how it extends to systems.
For system \(\mathbf{u}'=\mathbf{F}(t,\mathbf{u})\) with \(\mathbf{u}=(y,z)^T\), \(\mathbf{F}=(f,g)^T\): \[\mathbf{k}_1=h\mathbf{F}(t_n,\mathbf{u}_n),\quad\mathbf{k}_2=h\mathbf{F}(t_n+h/2,\mathbf{u}_n+\mathbf{k}_1/2),\] \[\mathbf{k}_3=h\mathbf{F}(t_n+h/2,\mathbf{u}_n+\mathbf{k}_2/2),\quad\mathbf{k}_4=h\mathbf{F}(t_n+h,\mathbf{u}_n+\mathbf{k}_3).\] \[\mathbf{u}_{n+1}=\mathbf{u}_n+\tfrac{1}{6}(\mathbf{k}_1+2\mathbf{k}_2+2\mathbf{k}_3+\mathbf{k}_4).\] In component form: compute \(k_1^y=hf(t_n,y_n,z_n)\), \(k_1^z=hg(t_n,y_n,z_n)\), and so on simultaneously at each stage. Local truncation error: \(O(h^5)\); global error: \(O(h^4)\). Requires 4 evaluations of \(\mathbf{F}\) per step. Extends to \(m\)-dimensional systems directly. \(\blacksquare\)
Describe the shooting method for the BVP \(y''-3y'+2y=0\), \(y(0)=1\), \(y(1)=e^2\). Find the exact solution and verify the shooting approach.
Exact solution: characteristic equation \(r^2-3r+2=(r-1)(r-2)=0\). General solution: \(y=C_1e^x+C_2e^{2x}\). BCs: \(y(0)=C_1+C_2=1\), \(y(1)=C_1e+C_2e^2=e^2\). From first: \(C_1=1-C_2\). Sub: \((1-C_2)e+C_2e^2=e^2\Rightarrow C_2(e^2-e)=e^2-e\Rightarrow C_2=1,C_1=0\). So \(y=e^{2x}\). Shooting method: guess \(y'(0)=s\), solve IVP numerically, get \(y(1;s)\). Seek \(s\) such that \(y(1;s)=e^2\). For this linear BVP: use two guesses \(s_1,s_2\), compute \(y(1;s_1),y(1;s_2)\), then by linearity: \(s^*=s_1+\frac{(e^2-y(1;s_1))(s_2-s_1)}{y(1;s_2)-y(1;s_1)}\) (one interpolation step). \(\blacksquare\)
Set up the finite difference equations for the BVP \(-y''=f(x)\), \(y(0)=y(1)=0\), using \(n\) equal subintervals of width \(h=1/n\).
Discretize \(-y''(x_i)\approx\frac{-y_{i-1}+2y_i-y_{i+1}}{h^2}=f_i\) at interior points \(x_i=ih\), \(i=1,\ldots,n-1\). With \(y_0=y_n=0\) (boundary conditions), the system is: \[\frac{1}{h^2}\begin{pmatrix}2&-1\\-1&2&-1\\&\ddots&\ddots&\ddots\\&&-1&2\end{pmatrix}\begin{pmatrix}y_1\\\vdots\\y_{n-1}\end{pmatrix}=\begin{pmatrix}f_1\\\vdots\\f_{n-1}\end{pmatrix}.\] Tridiagonal, diagonally dominant SPD matrix. Solve by Thomas algorithm in \(O(n)\). Consistency error: \(O(h^2)\) (second-order scheme). Global error in \(y_i\): also \(O(h^2)\) by maximum principle. \(\blacksquare\)
State the Lax equivalence theorem. Define consistency, stability, and convergence for finite difference methods for PDEs.
Consistency: local truncation error \(\tau_h\to0\) as \(h\to0\). Method approximates the PDE. Stability: numerical solution is bounded uniformly: \(\|u^n\|\le C\|u^0\|\) for all \(n,h\) with \(nk\le T\). Equivalently (von Neumann): all Fourier modes do not grow. Convergence: numerical solution \(u^n_j\to u(x_j,t_n)\) as \(h\to0\) (with appropriate \(k/h\) ratio). Lax equivalence theorem: For a well-posed linear PDE and a consistent finite difference scheme, stability is equivalent to convergence. Consequence: to prove convergence, it suffices to prove consistency (easy) and stability (often via energy methods or von Neumann analysis). \(\blacksquare\)
Write down the Lax-Wendroff scheme for the advection equation \(u_t+au_x=0\). Show it is second-order accurate and state the CFL stability condition.
Lax-Wendroff: \(u_j^{n+1}=u_j^n-\frac{ak}{2h}(u_{j+1}^n-u_{j-1}^n)+\frac{a^2k^2}{2h^2}(u_{j+1}^n-2u_j^n+u_{j-1}^n)\). Derivation: Taylor in time \(u_j^{n+1}\approx u_j^n+ku_t+\frac{k^2}{2}u_{tt}\); use PDE: \(u_t=-au_x\), \(u_{tt}=a^2u_{xx}\); approximate spatial derivatives by centered differences. Accuracy: LTE \(=O(k^2,h^2)\) — second order in both space and time. Stability (von Neumann): amplification factor \(|g(\xi)|^2=1-\nu^2(1-\cos\xi)\cdot[\text{terms}]\le1\) when \(|\nu|=|ak/h|\le1\). CFL condition: \(|a|k/h\le1\). \(\blacksquare\)
State the Crank-Nicolson scheme for the heat equation \(u_t=\alpha u_{xx}\). Prove it is unconditionally stable by von Neumann analysis.
Crank-Nicolson: \(\frac{u_j^{n+1}-u_j^n}{k}=\frac{\alpha}{2h^2}\left[(u_{j+1}^{n+1}-2u_j^{n+1}+u_{j-1}^{n+1})+(u_{j+1}^n-2u_j^n+u_{j-1}^n)\right]\). Von Neumann: substitute \(u_j^n=g^n e^{ij\xi}\): \[\frac{g-1}{k}=\frac{\alpha}{h^2}(e^{i\xi}-2+e^{-i\xi})\frac{g+1}{2}=\frac{-4\alpha\sin^2(\xi/2)}{h^2}\cdot\frac{g+1}{2}.\] Let \(\lambda=2\alpha k\sin^2(\xi/2)/h^2\ge0\): \(g-1=-\lambda(g+1)\Rightarrow g=\frac{1-\lambda}{1+\lambda}\). Since \(\lambda\ge0\): \(|g|=|(1-\lambda)/(1+\lambda)|\le1\) for all \(\lambda\ge0\), i.e., for all \(k,h,\xi\). Unconditionally stable. \(\blacksquare\)
Describe von Neumann stability analysis. Apply it to the forward-difference (FTCS) scheme for the heat equation to derive the stability condition \(r=\alpha k/h^2\le1/2\).
FTCS: \(u_j^{n+1}=u_j^n+r(u_{j+1}^n-2u_j^n+u_{j-1}^n)\), \(r=\alpha k/h^2\). Von Neumann method: assume \(u_j^n=g^n e^{ij\xi}\) (Fourier mode). Substitute: \[g^{n+1}e^{ij\xi}=g^n e^{ij\xi}+rg^n(e^{i(j+1)\xi}-2e^{ij\xi}+e^{i(j-1)\xi})\] \[g=1+r(e^{i\xi}-2+e^{-i\xi})=1-4r\sin^2(\xi/2).\] Stability: \(|g|\le1\) for all \(\xi\). \(\max_\xi|g|=|1-4r|\) (at \(\xi=\pi\)). Need \(|1-4r|\le1\Rightarrow 0\le r\le1/2\). Hence stability condition: \(r=\alpha k/h^2\le\tfrac{1}{2}\). \(\blacksquare\)
Describe the finite element method using linear (hat) basis functions for the BVP \(-u''=f\), \(u(0)=u(1)=0\). Derive the stiffness matrix.
Weak formulation: find \(u\in H^1_0(0,1)\) such that \(\int_0^1u'v'\,dx=\int_0^1fv\,dx\) for all \(v\in H^1_0\). Partition \([0,1]\) into \(n\) elements of size \(h=1/n\). Basis: hat functions \(\phi_i(x)\): piecewise linear, \(\phi_i(x_j)=\delta_{ij}\). Stiffness matrix \(K_{ij}=\int_0^1\phi_i'\phi_j'\,dx\). For adjacent nodes \(|i-j|=1\): \(K_{ij}=-1/h\). For \(i=j\): \(K_{ii}=2/h\) (interior). Load vector \(b_i=\int_0^1f\phi_i\,dx\approx hf(x_i)\). System: \(\frac{1}{h}\begin{pmatrix}2&-1\\-1&2&-1\\&\ddots\end{pmatrix}\mathbf{u}=\mathbf{b}\). Same as finite difference second-order scheme! \(O(h^2)\) convergence in \(L^2\). \(\blacksquare\)
Describe Lanczos's tau method for solving the BVP \(y''+y=x\), \(y(0)=y(1)=0\). How does it differ from the Galerkin method?
Tau method: seek approximate solution \(y_n(x)=\sum_{k=0}^n a_k T_k(x)\) (Chebyshev expansion) satisfying the ODE approximately. Instead of requiring the residual to be orthogonal to all basis functions (Galerkin), the tau method: (1) expands all functions in Chebyshev series; (2) the differential equation gives a relation between expansion coefficients; (3) boundary conditions replace the last 2 equations (for second-order ODE). The residual is of the form \(R_n(x)=\tau_1 T_{n+1}(x)+\tau_2 T_{n+2}(x)\) (tail terms). Key advantage over Galerkin: exact boundary condition enforcement while maintaining spectral accuracy. Produces banded algebraic systems for polynomial bases. Exponential convergence for smooth \(f\). \(\blacksquare\)
Describe the collocation method for ODEs. For the BVP \(y''=f(x,y)\), \(y(a)=y(b)=0\), how are the collocation points chosen to achieve high accuracy?
Collocation: seek \(y_n\in P_n\) (polynomial of degree \(n\)) such that \(y_n''(x_i)=f(x_i,y_n(x_i))\) at \(n-1\) interior collocation points \(x_1,\ldots,x_{n-1}\), plus two boundary conditions. This gives \(n+1\) equations for \(n+1\) coefficients. Choice of collocation points: (1) Equally spaced: only \(O(h^n)\) accuracy. (2) Gauss-Legendre points (roots of \(P_n\)): superconvergence — order \(2n\). (3) Gauss-Lobatto points (include endpoints): order \(2n-2\), useful for BVPs. Chebyshev collocation: near-minimax approximation, avoids Runge phenomenon. Spectral collocation on \(N\) Chebyshev nodes: exponential convergence for smooth solutions. \(\blacksquare\)
Write a complete algorithm for the bisection method for finding a root of \(f(x)=0\) in \([a,b]\) with stopping criterion \(|f(c)|<\varepsilon\).
Algorithm Bisection(f, a, b, ε):
1. If \(f(a)\cdot f(b)>0\), STOP: no sign change; root not guaranteed.
2. Set \(c=(a+b)/2\).
3. If \(|f(c)|<\varepsilon\) OR \((b-a)/2<\varepsilon\), RETURN \(c\) as root. (Two criteria: function value small OR interval small.)
4. If \(f(a)\cdot f(c)<0\): set \(b=c\). Else: set \(a=c\).
5. GOTO step 2.
Properties: always converges (linear convergence); after \(n\) iterations, error \(\le(b-a)/2^n\). Number of iterations needed: \(n\ge\lceil\log_2((b-a)/\varepsilon)\rceil\). Order of convergence: 1 (linear). Advantage: simple, robust. Disadvantage: slow (one bit per iteration). \(\blacksquare\)
Draw a flowchart for the Newton-Raphson iteration \(x_{n+1}=x_n-f(x_n)/f'(x_n)\) with stopping criterion \(|x_{n+1}-x_n|<\varepsilon\). State when it fails.
Flowchart steps:
[START] → Input \(x_0,\varepsilon,\text{MaxIter}\) → Set \(n=0\)
→ [Compute \(f(x_n),f'(x_n)\)]
→ [Is \(f'(x_n)=0\)?] → YES: STOP (division by zero failure)
→ NO: \(x_{n+1}=x_n-f(x_n)/f'(x_n)\)
→ [Is \(|x_{n+1}-x_n|<\varepsilon\)?] → YES: Output \(x_{n+1}\), STOP
→ NO: [Is \(n\ge\text{MaxIter}\)?] → YES: STOP (no convergence)
→ NO: \(n=n+1\), \(x_n=x_{n+1}\), go to Compute step.
Failure cases: \(f'(x_n)=0\); oscillation (e.g., \(f(x)=x^{1/3}\)); divergence if poor initial guess or flat region. \(\blacksquare\)
Write the bubble sort algorithm for sorting \(n\) elements. Analyze its time complexity.
Algorithm BubbleSort(A, n):
For \(i=1\) to \(n-1\):
For \(j=1\) to \(n-i\):
If \(A[j]>A[j+1]\): swap \(A[j]\) and \(A[j+1]\)
(After pass \(i\), the \(i\)-th largest element is in its final position.)
Complexity:
Worst case (reverse sorted): \(\sum_{i=1}^{n-1}(n-i)=n(n-1)/2\) comparisons → \(O(n^2)\).
Best case (already sorted, with flag): \(O(n)\) — one pass with no swaps, flag set.
Average: \(O(n^2)\). Space: \(O(1)\) (in-place). Stable sort. Not suitable for large \(n\). \(\blacksquare\)
Write the binary search algorithm for finding target \(T\) in a sorted array \(A[1..n]\). Analyze its time complexity.
Algorithm BinarySearch(A, n, T):
Set \(\ell=1\), \(r=n\).
While \(\ell\le r\):
\(m=\lfloor(\ell+r)/2\rfloor\)
If \(A[m]=T\): RETURN \(m\)
If \(A[m]
RETURN "not found"
Complexity: Each iteration halves the search space. After \(k\) steps, search space \(\le n/2^k\). Terminates when \(n/2^k<1\Rightarrow k>\log_2 n\). Time: \(O(\log n)\). Space: \(O(1)\) (iterative). Pre-requisite: array must be sorted. Comparison: linear search \(O(n)\); binary search \(O(\log n)\). \(\blacksquare\)
Define stack and queue. Write algorithms for push/pop (stack) and enqueue/dequeue (queue). What are their applications in numerical computing?
Stack (LIFO): Push(\(x\)): if top \(
Describe a singly linked list. Write algorithms for inserting a node after a given node and deleting a given node.
Node structure: \(\{data, \text{next}\}\). Head points to first node; last node's next is NULL.
Insert after node \(p\):
1. Create new node \(q\) with \(q.data=x\).
2. \(q.\text{next}=p.\text{next}\).
3. \(p.\text{next}=q\). Time: \(O(1)\).
Delete node after \(p\):
1. Let \(q=p.\text{next}\).
2. \(p.\text{next}=q.\text{next}\).
3. Free \(q\). Time: \(O(1)\).
Delete given node \(v\) (unknown predecessor): traverse from head to find predecessor \(p\): \(O(n)\). Then delete as above. Advantage over arrays: \(O(1)\) insert/delete if predecessor known; no shifting. Disadvantage: \(O(n)\) access by index (no random access). \(\blacksquare\)
Define Big-O notation. Give examples of algorithms with complexities \(O(1)\), \(O(\log n)\), \(O(n)\), \(O(n\log n)\), and \(O(n^2)\).
\(f(n)=O(g(n))\) if there exist \(C>0,n_0\) such that \(f(n)\le Cg(n)\) for all \(n\ge n_0\). Describes asymptotic upper bound on growth rate.
\(O(1)\) — constant: array element access \(A[i]\); stack push/pop.
\(O(\log n)\) — logarithmic: binary search; balanced BST lookup.
\(O(n)\) — linear: linear search; computing sum/max of array.
\(O(n\log n)\) — linearithmic: merge sort; heap sort; FFT.
\(O(n^2)\) — quadratic: bubble/insertion/selection sort; matrix-vector product for dense matrix; Gaussian elimination \(O(n^3)\) — but that's cubic.
Also: \(O(n^3)\) — matrix multiplication (naive); \(O(2^n)\) — exhaustive search.
Practical significance: for numerical problems (\(n=10^6\)), \(O(n\log n)\) is feasible (\(\sim2\times10^7\) ops) but \(O(n^2)\) is not (\(\sim10^{12}\)). \(\blacksquare\)
Write recursive algorithms for (i) \(n!\) and (ii) the \(n\)-th Fibonacci number. Analyse their time and space complexity.
(i) Factorial: \(n!=1\) if \(n=0\), else \(n\cdot(n-1)!\). Recursion: calls itself \(n\) times. Time: \(O(n)\). Space: \(O(n)\) (call stack depth \(n\)). Can be done iteratively in \(O(1)\) space. (ii) Fibonacci: \(F_n=F_{n-1}+F_{n-2}\), \(F_0=0,F_1=1\). Naive recursive: exponential \(O(\phi^n)\) (\(\phi=(1+\sqrt5)/2\)) — recomputes subproblems. With memoization (dynamic programming): \(O(n)\) time, \(O(n)\) space. Iterative: \(O(n)\) time, \(O(1)\) space (store only last two values). Matrix exponentiation: \(O(\log n)\) time. Tail recursion: factorial can be made tail-recursive; Fibonacci cannot easily (without accumulator). \(\blacksquare\)
Describe merge sort. Prove that its worst-case time complexity is \(O(n\log n)\). Compare with bubble sort.
Merge sort: divide array into two halves, recursively sort each, then merge. Merge(A, B): compare front elements, take smaller, repeat; \(O(n)\). MergeSort(A[1..n]): if \(n\le1\) return; split; MergeSort(left); MergeSort(right); Merge. Complexity (Recurrence): \(T(n)=2T(n/2)+O(n)\). By Master Theorem (case 2): \(T(n)=O(n\log n)\). Works for all inputs (no worst case variation). vs Bubble sort: Bubble sort \(O(n^2)\) worst/average; merge sort \(O(n\log n)\) always. For \(n=10^6\): merge sort \(\approx2\times10^7\) ops vs bubble sort \(\approx10^{12}\). Merge sort needs \(O(n)\) extra space; bubble sort \(O(1)\). Both stable. \(\blacksquare\)
Describe the quicksort algorithm with Lomuto's partition scheme. Analyse its average \(O(n\log n)\) and worst-case \(O(n^2)\) complexity.
QuickSort: choose pivot \(p=A[r]\) (last element); partition array around \(p\); recursively sort sub-arrays. Partition: maintain \(i\) (boundary of smaller elements). For \(j=\ell\) to \(r-1\): if \(A[j]\le p\): swap \(A[i+1]\) and \(A[j]\), \(i{+}{+}\). Swap \(A[i+1]\) with \(A[r]\). Return \(i+1\). Average case: pivot splits array roughly in half on average. Recurrence \(T(n)=2T(n/2)+O(n)\to O(n\log n)\). Worst case: always smallest or largest element chosen as pivot (sorted array): \(T(n)=T(n-1)+O(n)\to O(n^2)\). Remedy: random pivot or median-of-three. In-place: \(O(\log n)\) stack space on average. Fastest in practice due to cache efficiency. \(\blacksquare\)
Describe a binary search tree (BST). Write algorithms for insertion and deletion. State the complexity for balanced vs unbalanced trees.
BST property: for each node, left subtree keys \(<\) node key \(<\) right subtree keys. Insert(k): start at root; go left if \(k<\)node.key, right if \(k>\)node.key; insert when NULL reached. Delete(k): find node with key \(k\). (1) Leaf: simply remove. (2) One child: replace node with its child. (3) Two children: find in-order successor (leftmost node of right subtree); copy its key to node; delete the successor (which has at most one child). Complexity: All operations \(O(h)\) where \(h\) = height. Balanced (AVL, red-black): \(h=O(\log n)\), ops \(O(\log n)\). Degenerate (sorted input): \(h=O(n)\), ops \(O(n)\). \(\blacksquare\)
Describe Breadth-First Search (BFS) and Depth-First Search (DFS) for a graph. Compare their time complexity and applications.
BFS: Start from source \(s\); use queue. Mark \(s\) visited, enqueue. While queue nonempty: dequeue \(u\); for each unvisited neighbor \(v\): mark visited, enqueue \(v\). Explores all nodes at distance \(d\) before distance \(d+1\). DFS: Start from source; use stack (or recursion). Mark visited; explore each unvisited neighbor recursively. Complexity: Both \(O(V+E)\) for adjacency list representation. Applications: BFS: shortest path in unweighted graph; connected components; bipartiteness testing. DFS: topological sort; detecting cycles; strongly connected components (Tarjan/Kosaraju); tree/graph traversal; backtracking. Key difference: BFS = FIFO queue (levels); DFS = LIFO stack (depth-first). \(\blacksquare\)
Describe Dijkstra's algorithm for finding shortest paths from a source vertex in a weighted graph with non-negative weights. State its complexity.
Algorithm: Initialize \(d[s]=0\), \(d[v]=\infty\) for \(v\neq s\). Set \(Q=\)all vertices.
While \(Q\neq\emptyset\): extract \(u=\arg\min_{v\in Q}d[v]\) (min-priority queue); remove from \(Q\). For each neighbor \(v\) of \(u\): if \(d[u]+w(u,v)
Define the Longest Common Subsequence (LCS) problem. Derive a dynamic programming solution and state its complexity. Find LCS of "ABCBDAB" and "BDCAB".
LCS of sequences \(X=x_1\ldots x_m\) and \(Y=y_1\ldots y_n\): longest sequence \(Z\) that is a subsequence of both. Recurrence: \(L[i][j]=\begin{cases}0&i=0\text{ or }j=0\\L[i-1][j-1]+1&x_i=y_j\\\max(L[i-1][j],L[i][j-1])&x_i\neq y_j\end{cases}\). Fill \((m+1)\times(n+1)\) table; LCS length = \(L[m][n]\). Example: \(X=\)ABCBDAB (\(m=7\)), \(Y=\)BDCAB (\(n=5\)). LCS = "BCAB" or "BDAB", length 4. Complexity: \(O(mn)\) time, \(O(mn)\) space (reducible to \(O(\min(m,n))\) for length only). Reconstruction: trace back through table. Applications: diff tools, bioinformatics (DNA alignment). \(\blacksquare\)
Describe Kruskal's algorithm for finding a Minimum Spanning Tree (MST) of a connected weighted graph. Prove its correctness using the greedy exchange argument.
Algorithm: Sort all edges by weight. Initialize forest (each vertex its own component). For each edge \((u,v)\) in sorted order: if \(u\) and \(v\) in different components: add edge to MST, merge components (Union-Find). Stop when \(V-1\) edges added. Correctness (Cut property): For any cut \((S,V\setminus S)\), the minimum weight crossing edge \(e^*\) is in some MST. If we select \(e^*\) (greedy choice), it belongs to an MST. If any heavier edge \(e\) is in the current MST crosses the same cut, swapping \(e\) for \(e^*\) reduces total weight or maintains it — contradiction with minimality of chosen edge. Complexity: Sort \(O(E\log E)\); Union-Find with path compression and union by rank: \(O(E\alpha(V))\approx O(E)\). Total: \(O(E\log E)\). \(\blacksquare\)
Describe a hash table with separate chaining for collision resolution. Analyze its average-case time complexity for search, insert, and delete.
Hash table: array of \(m\) slots; hash function \(h(k)\in\{0,\ldots,m-1\}\). Each slot heads a linked list of all keys that hash to it. Insert: compute \(h(k)\); insert at front of list at slot \(h(k)\). \(O(1)\) worst case. Search: compute \(h(k)\); traverse list at \(h(k)\). Average: \(O(1+\alpha)\) where \(\alpha=n/m\) (load factor). Delete: search then remove. \(O(1+\alpha)\) average. Analysis: With simple uniform hashing and \(\alpha=O(1)\) (i.e., \(m=\Theta(n)\)): all operations \(O(1)\) average. Worst case (all keys same hash): \(O(n)\). Hash function: division method \(h(k)=k\bmod m\) (choose \(m\) prime); multiplication method. Alternative: open addressing (linear probing, quadratic probing). \(\blacksquare\)
Describe the IEEE 754 double-precision floating-point format. How are numbers stored? What is the range and precision?
IEEE 754 double (64-bit): 1 sign bit, 11 exponent bits (biased by 1023), 52 mantissa (fraction) bits. Stored value: \((-1)^s\times1.f\times2^{e-1023}\) (for normal numbers, where \(1\le e\le2046\)). Special values: \(e=0\): subnormal (no implicit 1) or zero. \(e=2047\): \(\pm\infty\) or NaN. Range: \(\approx2.2\times10^{-308}\) to \(1.8\times10^{308}\). Precision: 52-bit mantissa → \(\approx15-16\) significant decimal digits. Machine epsilon \(\epsilon_M=2^{-52}\approx2.22\times10^{-16}\). Arithmetic: operations rounded to nearest representable number (round-to-nearest-even by default). Guarantees \(|\text{fl}(a\circ b)-(a\circ b)|/|a\circ b|\le\epsilon_M/2=u\). \(\blacksquare\)
Define machine epsilon \(\varepsilon_M\). Describe an algorithm to compute it experimentally. What is its value for IEEE 754 double precision?
Definition: Machine epsilon \(\varepsilon_M\) is the smallest positive floating-point number such that \(\text{fl}(1+\varepsilon_M)\neq1\) in computer arithmetic. Equivalently, the distance from 1 to the next larger floating-point number.
Algorithm to compute \(\varepsilon_M\):
Set \(\varepsilon=1\).
While \(\text{fl}(1+\varepsilon/2)\neq1\): \(\varepsilon\leftarrow\varepsilon/2\).
Return \(\varepsilon\). (Halve until adding it to 1 no longer changes the result.)
IEEE 754 double: \(\varepsilon_M=2^{-52}\approx2.22\times10^{-16}\). Unit roundoff \(u=\varepsilon_M/2=2^{-53}\).
Single precision: \(\varepsilon_M=2^{-23}\approx1.19\times10^{-7}\).
Significance: relative rounding error in any floating-point operation \(\le u\). Governs the ultimate accuracy achievable in arithmetic. \(\blacksquare\)
Explain catastrophic cancellation. Give an example where it causes loss of significant digits and describe how to avoid it.
Catastrophic cancellation: when two nearly equal floating-point numbers are subtracted, significant digits are lost. Result has few correct digits even if operands are accurately represented. Example: \(x=1.23456789012345\), \(y=1.23456789012344\). True difference: \(10^{-14}\). In 15-digit arithmetic, both round to same 15 digits, difference may be 0 or have only 1 correct significant digit. Classic: solving quadratic \(ax^2+bx+c=0\): root \(x=(-b+\sqrt{b^2-4ac})/(2a)\). If \(b>0\) and \(b^2\gg4ac\): \(-b+\sqrt{b^2-4ac}\approx0\) — catastrophic cancellation. Avoidance: compute \(x_1=-b-\text{sgn}(b)\sqrt{b^2-4ac}/(2a)\), then \(x_2=c/(ax_1)\) (Vieta's). Other remedies: use extended precision; rearrange formulas; compute \(\ln(1+x)\) instead of \(\log(1+x)\) for small \(x\) using \texttt{log1p()}. \(\blacksquare\)
Define forward error, backward error, and backward stability for a numerical algorithm. Illustrate with the example of computing the inner product \(\mathbf{x}^T\mathbf{y}\).
Forward error: \(\|\hat{f}(x)-f(x)\|\) — how far the computed answer \(\hat{f}(x)\) is from the true answer \(f(x)\). Backward error: smallest \(\|\delta x\|\) such that \(\hat{f}(x)=f(x+\delta x)\) — how much we must perturb the input to make our computed answer exactly right. Backward stability: algorithm is backward stable if \(\|\delta x\|/\|x\|=O(u)\) (machine epsilon) for all \(x\). Backward stable \(\Rightarrow\) forward error \(\le\kappa(f)\cdot O(u)\). Inner product \(\hat{s}=\sum_{i=1}^n\hat{x}_iy_i\): floating-point evaluation satisfies \(\hat{s}=\mathbf{x}^T\mathbf{y}+O(nu\|\mathbf{x}\|\|\mathbf{y}\|)\). Relative error \(\le nu\), so backward stable for small \(n\). For large \(n\): pairwise summation or compensated summation (Kahan) reduces error to \(O(u)\) regardless of \(n\). \(\blacksquare\)
Describe the Gauss-Seidel iteration method for solving the linear system \(Ax=b\). Apply it to solve \(4x+y=9\), \(x+3y=10\) starting from \((0,0)\). Perform 3 iterations.
Gauss-Seidel uses the latest updated values immediately: \(x^{(k+1)}=\dfrac{9-y^{(k)}}{4}\), \(y^{(k+1)}=\dfrac{10-x^{(k+1)}}{3}\).
Iter 1: \(x=\frac{9-0}{4}=2.25\), \(y=\frac{10-2.25}{3}=2.583\).
Iter 2: \(x=\frac{9-2.583}{4}=1.604\), \(y=\frac{10-1.604}{3}=2.799\).
Iter 3: \(x=\frac{9-2.799}{4}=1.550\), \(y=\frac{10-1.550}{3}=2.817\).
Converges to exact \(x=\frac{17}{11}\approx1.545\), \(y=\frac{31}{11}\approx2.818\). Diagonal dominance (\(4>1\), \(3>1\)) guarantees convergence. \(\blacksquare\)
Convert \((125.75)_{10}\) to binary and hexadecimal. Convert \((1101.011)_2\) to decimal.
\(125=(1111101)_2\) (repeated division by 2). \(0.75=(0.11)_2\) (since \(0.75\times2=1.5\to1\), \(0.5\times2=1.0\to1\)).
So \((125.75)_{10}=(1111101.11)_2\). Grouping in nibbles from point: \(0111\ 1101.\ 1100=(7D.C)_{16}\).
For \((1101.011)_2\): \(8+4+0+1+0+\frac{1}{4}+\frac{1}{8}=13+0.25+0.125=13.375\). \(\blacksquare\)
Simplify the Boolean expression \(F=AB'+A'B+AB\) and draw the corresponding logic circuit.
\(F=AB'+A'B+AB=(AB'+AB)+A'B=A(B'+B)+A'B=A\cdot1+A'B=A+A'B\).
By absorption: \(A+A'B=(A+A')(A+B)=1\cdot(A+B)=A+B\).
The circuit is simply an OR gate with inputs \(A\) and \(B\). \(\blacksquare\)
Minimize the Boolean function \(F(A,B,C,D)=\sum m(0,1,3,5,7,8,9,11,13,15)\) using a Karnaugh map.
Plot 10 minterms on a 4-variable K-map (rows: \(AB\), columns: \(CD\) in Gray code).
Groups identified:
• \(\{0,1,8,9\}\): \(A\) varies, \(B=0,C=0\) → \(\bar{B}\bar{C}\)
• \(\{1,3,5,7\}\): \(A=0,D=1\) → \(\bar{A}D\)
• \(\{5,7,13,15\}\): \(B=1,D=1\) → \(BD\)
• \(\{9,11,13,15\}\): \(A=1,D=1\) → \(AD\)
Minimal SOP (using essential primes \(\bar{B}\bar{C}\) and \(AD\), plus \(BD\) to cover \(\{5,7\}\)): \[F=\bar{B}\bar{C}+AD+BD=\bar{B}\bar{C}+D(A+B).\] \(\blacksquare\)
Use Lagrange's interpolation formula to find \(f(3)\) given: \(f(0)=2,\,f(1)=3,\,f(4)=5,\,f(5)=7\).
Nodes: \(x_0=0,x_1=1,x_2=4,x_3=5\). At \(x=3\):
\(L_0=\frac{(3-1)(3-4)(3-5)}{(0-1)(0-4)(0-5)}=\frac{(2)(-1)(-2)}{(-1)(-4)(-5)}=\frac{4}{-20}=-\frac{1}{5}\).
\(L_1=\frac{(3)(3-4)(3-5)}{(1)(1-4)(1-5)}=\frac{3(-1)(-2)}{1(-3)(-4)}=\frac{6}{12}=\frac{1}{2}\).
\(L_2=\frac{(3)(3-1)(3-5)}{(4)(4-1)(4-5)}=\frac{3(2)(-2)}{4(3)(-1)}=\frac{-12}{-12}=1\).
\(L_3=\frac{(3)(3-1)(3-4)}{(5)(5-1)(5-4)}=\frac{3(2)(-1)}{5(4)(1)}=\frac{-6}{20}=-\frac{3}{10}\).
\(f(3)=2(-\tfrac{1}{5})+3(\tfrac{1}{2})+5(1)+7(-\tfrac{3}{10})=-0.4+1.5+5-2.1=4.0\). \(\blacksquare\)
State the trapezoidal rule. Use it with \(h=0.25\) to evaluate \(\int_0^1 e^{-x^2}dx\). Compare with Simpson's rule. Find the error bound.
Trapezoidal rule: \(\int_a^b f\,dx\approx\dfrac{h}{2}[f(x_0)+2f(x_1)+\cdots+2f(x_{n-1})+f(x_n)]\). Error bound: \(|E_T|\le\dfrac{(b-a)h^2}{12}\max|f''|\).
Points \(x=0,0.25,0.5,0.75,1\). Values: \(f(0)=1,\;f(0.25)\approx0.9394,\;f(0.5)\approx0.7788,\;f(0.75)\approx0.5698,\;f(1)\approx0.3679\).
Trapezoidal: \(\frac{0.25}{2}[1+2(0.9394+0.7788+0.5698)+0.3679]=\frac{0.25}{2}[1+4.5760+0.3679]=\frac{0.25}{2}(5.9439)\approx0.7430\).
For Simpson's 1/3 (same 5 pts): \(\frac{0.25}{3}[1+4(0.9394)+2(0.7788)+4(0.5698)+0.3679]=\frac{0.25}{3}(8.9613)\approx0.7468\).
Exact \(\approx0.7468\). Trapezoidal error \(\approx0.0038\); Simpson error \(\approx0\) (4th-order). Error bound for trapezoidal: \(|f''(x)|=|2e^{-x^2}(2x^2-1)|\le2\), so \(|E_T|\le\frac{(1)(0.25)^2}{12}(2)\approx0.0104\). \(\blacksquare\)
Describe Newton-Raphson method. Prove that it has quadratic convergence near a simple root. Apply it to find \(\sqrt[3]{28}\) correct to 4 decimal places.
NR formula: \(x_{n+1}=x_n-\dfrac{f(x_n)}{f'(x_n)}\). Quadratic convergence: error \(e_{n+1}\approx\dfrac{f''(\xi)}{2f'(\xi)}e_n^2\) (see proof in §Proofs above).
For \(\sqrt[3]{28}\): \(f(x)=x^3-28\), \(f'(x)=3x^2\). Formula: \(x_{n+1}=x_n-\dfrac{x_n^3-28}{3x_n^2}=\dfrac{2x_n^3+28}{3x_n^2}\).
Start \(x_0=3\): \(x_1=\dfrac{2(27)+28}{3(9)}=\dfrac{82}{27}\approx3.0370\).
\(x_2=\dfrac{2(3.0370)^3+28}{3(3.0370)^2}=\dfrac{2(28.04)+28}{27.67}=\dfrac{84.08}{27.67}\approx3.0366\).
\(x_3\approx3.0366\) (converged). \(\sqrt[3]{28}\approx3.0366\). \(\blacksquare\)
What are logic gates? Define AND, OR, NOT, NAND, NOR, XOR. Draw truth tables for NAND and XOR. Show NAND is universal.
Logic gates implement Boolean operations on binary inputs.
NAND truth table: \(00\to1,\;01\to1,\;10\to1,\;11\to0\). XOR truth table: \(00\to0,\;01\to1,\;10\to1,\;11\to0\).
NAND is universal (any Boolean function buildable from NAND alone):
• NOT: \(\overline{A}=A\text{ NAND }A\)
• AND: \(A\cdot B=(A\text{ NAND }B)\text{ NAND }(A\text{ NAND }B)\)
• OR: \(A+B=(A\text{ NAND }A)\text{ NAND }(B\text{ NAND }B)\)
Since AND, OR, NOT form a functionally complete set, and all three are realizable with NAND, NAND alone suffices for any Boolean function. \(\blacksquare\)
Represent \(-37\) in 8-bit signed magnitude, one's complement, and two's complement. Also represent \(1.5\) in IEEE 754 single-precision format.
\(37=(00100101)_2\).
• Signed magnitude \(-37\): \(\mathbf{1}0100101\) (flip sign bit).
• One's complement: \(\mathbf{11011010}\) (invert all bits of \(00100101\)).
• Two's complement: \(\mathbf{11011011}\) (one's complement + 1).
IEEE 754 single for \(1.5=1.1_2\times2^0\): sign bit \(=0\), exponent \(=127+0=127=(01111111)_2\), mantissa (after implicit 1) \(=1\underbrace{00\cdots0}_{22}\).
Full 32-bit: \(0\;01111111\;10000000000000000000000\). \(\blacksquare\)
Convert \((2A3F)_{16}\) to binary and decimal. Convert \((0.1)_{10}\) to binary (show it is a repeating fraction).
\(2=0010,\;A=1010,\;3=0011,\;F=1111\). So \((2A3F)_{16}=(0010\,1010\,0011\,1111)_2\).
Decimal: \(2\times16^3+10\times16^2+3\times16+15=8192+2560+48+15=10815\).
For \((0.1)_{10}\) to binary (multiply-by-2 method):
\(0.1\times2=0.2\to0;\;0.2\times2=0.4\to0;\;0.4\times2=0.8\to0;\;0.8\times2=1.6\to1;\;0.6\times2=1.2\to1;\;0.2\times2=0.4\to0\;\cdots\) (repeats).
\((0.1)_{10}=(0.\overline{0001100110011\ldots})_2=(0.0\overline{0011})_2\). It is a non-terminating repeating binary fraction. \(\blacksquare\)
Use Simpson's 1/3 rule with 6 intervals to evaluate \(\int_0^6\dfrac{1}{1+x^2}dx\). Compare with the exact value \(\arctan 6\).
\(h=1\), \(n=6\) (even ✓). Points \(x=0,1,2,3,4,5,6\). \(f(x)=1/(1+x^2)\):
\(f(0)=1,\;f(1)=0.5,\;f(2)=0.2,\;f(3)=0.1,\;f(4)=\tfrac{1}{17}\approx0.0588,\;f(5)=\tfrac{1}{26}\approx0.0385,\;f(6)=\tfrac{1}{37}\approx0.0270\).
Simpson's 1/3: \(\dfrac{h}{3}[f_0+4f_1+2f_2+4f_3+2f_4+4f_5+f_6]\)
\(=\dfrac{1}{3}[1+2+0.4+0.4+0.1176+0.154+0.0270]=\dfrac{1}{3}(4.0986)\approx1.3662\).
Exact: \(\arctan 6\approx1.4056\). Error \(\approx0.0394\). (Large error due to coarse \(h=1\); function has significant curvature near \(x=0\).) \(\blacksquare\)
Derive the RK4 method. Apply it with \(h=0.1\) to solve \(y'=x+y\), \(y(0)=1\), for \(y(0.1)\).
RK4 (weights from Simpson's 1/3 on slope estimates at \(x_n,\,x_n+h/2,\,x_n+h/2,\,x_n+h\)): \[k_1=hf(x_n,y_n),\quad k_2=hf\!\left(x_n+\tfrac{h}{2},y_n+\tfrac{k_1}{2}\right),\] \[k_3=hf\!\left(x_n+\tfrac{h}{2},y_n+\tfrac{k_2}{2}\right),\quad k_4=hf(x_n+h,y_n+k_3),\] \[y_{n+1}=y_n+\tfrac{1}{6}(k_1+2k_2+2k_3+k_4).\] Local truncation error \(O(h^5)\), global error \(O(h^4)\).
At \(x_0=0,y_0=1\), \(f(x,y)=x+y\):
\(k_1=0.1(0+1)=0.1\).
\(k_2=0.1(0.05+1+0.05)=0.1(1.1)=0.11\).
\(k_3=0.1(0.05+1+0.055)=0.1(1.105)=0.1105\).
\(k_4=0.1(0.1+1+0.1105)=0.1(1.2105)=0.12105\).
\(y(0.1)=1+\tfrac{1}{6}(0.1+0.22+0.221+0.12105)=1+\tfrac{0.66205}{6}=1+0.11034=1.11034\).
Exact: \(y=2e^x-x-1\Rightarrow y(0.1)=2e^{0.1}-1.1\approx1.11034\). \(\checkmark\) \(\blacksquare\)