2.6 SOLUTION OF TRIDIAGONAL LINEAR SYSTEMS
Exercises:
1 Use the tridiagonal algorithm in this section to compute the solution to the following system of equations:Solution: After the elimination step is completed, the triangular system isand the solution is
2 Write a computer code to solve the previous problemSolution: The following is a MATLAB script that produced the previous solutionfunction [delta, f, x] = trisol(l,d,u,b) n = length(d); x = zeros(n,1); for k=2:n d(k) = d(k) - u(k-1)*l(k)/d(k-1); b(k) = b(k) - b(k-1)*l(k)/d(k-1); end x(n) = b(n)/d(n); for k=(n-1):(-1):1 x(k) = (b(k) - u(k)*x(k+1))/d(k); end delta = d; f = b;
3 Use the algorithm of this section to solve the following system of equations:You should get the solution .
4 Write a computer code to solve the previous problem.
5 The diagonal dominance condition is an example of a sufficient but not necessary condition. That is, the algorithm will often work for systems that are not diagonally dominant. Show that the following system is not diagonally dominant, but then use the tridiagonal algorithm in this section to compute the solution to it:You should get the solution .Solution: The matrix fails to be diagonally dominant because of the values in the second, third, and fourth rows. However, if we apply the algorithm, we get the triangular systemand then the correct values for the solution, . (It always helps to type in the correct right‐side vector when checking these things.)
6 Use the tridiagonal algorithm in this section to compute the solution to the following system of equations:Note that this is a very small change from the previous problem, since the only difference is that has changed by only . How much has the answer changed?Solution: The solution is nowConsidering the small change in the problem, this is a substantial change in the solution.
7 Write a computer code to do the previous two problems.
8 Verify that the following system is diagonally dominant and use the algorithm of this section to find the solution.
9 Use the algorithm of this section to find the solution to this system:Note that the right side here is different from that in the previous problem by only a small amount in the component. Comment on your results here as compared to those in the previous problem.Solution: The solution is nowAgain, this is a very large change in the solution for so modest a change in the problem.
10 Write a computer code to do the previous two problems.
11 Write a code that carries out the tridiagonal solution algorithm, and test it on the following system of equations,where is withand for all . Check your results by computing the residual . What is the largest component (in absolute value) of ? (You could also check your results by using MATLAB's backslash operator to solve the system.)
12 Extend the tridiagonal algorithm to a pentadiagonal matrix, i.e., one with five non‐zero diagonals. Write a program to carry out this solution algorithm, and apply it to the systemCheck your results by again computing the residual vector, or by using MATLAB's backslash operation.Solution: A MATLAB script for doing this is given below. On the example in the exercise it returns the (exact) solution .function [delta, f, x] = pentasol(k,l,d,u,v,b) n = length(d); x = zeros(n,1); d(2) = d(2) - l(2)*u(1)/d(1); u(2) = u(2) - l(2)*v(1)/d(1); b(2) = b(2) - l(2)*b(1)/d(1); for j=3:n l(j) = l(j) - k(3)*u(j-2)/d(j-2); d(j) = d(j) - k(j)*v(j-2)/d(j-2); b(j) = b(j) - k(j)*b(j-2)/d(j-2); d(j) = d(j) - l(j)*u(j-1)/d(j-1); if j < n u(j) = u(j) - l(j)*v(j-1)/d(j-1); end b(j) = b(j) - l(j)*b(j-1)/d(j-1); end delta = d; f = b; x = f; % x(n) = f(n)/d(n); x(n-1) = (b(n-1) - u(n-1)*x(n))/d(n-1); for j=(n-2):(-1):1 x(j) = (b(j) - u(j)*x(j+1) - v(j)*x(j+2))/d(j) end
13 Consider the family of tridiagonal problems defined by the matrix , withand a randomly defined right‐hand‐side vector. (Use rand to generate the random vectors.) Solve the system over the range ; use the appropriate timing functions to estimate the CPU time required for each case, and plot the result as a function of . Comment on your results. (You will probably want to use a semilog plot.)Solution: I got the plot in Fig. 2.5. While this is a “noisy” plot, it is very much generally a logarithmic curve, which suggests the actual timing numbers are a straight line, which is what we would expect, given that the cost of doing the solution is linear in the size of the matrix.Figure 2.5 Estimated timing cost for Problem 13, semilog scale.
2.7 APPLICATION: SIMPLE TWO‐POINT BOUNDARY VALUE PROBLEMS
Exercises:
1 Solve, by hand, the two‐point BVP (2.28) and (2.29) when , using . Write out the linear system explicitly prior to solution. You should get the following system:Solution: The approximate solution is
2 Repeat the above, this time using . What is the system now?
3 Repeat it again, this time using . What is the system now? (For this problem, you probably will want to use a computer code to actually solve the system.)Solution: The approximate solution is
4 Write a program which solves the two‐point BVP (2.28) and (2.29) where is as given below:;;;.The exact solutions are as given. Using , do we get the same kind of accuracy as in Table 2.10? Explain why or why not.Solution: The approximate solutions will display the kind of accuracy suggested in the text, except for (d); in this case, the singularity in the logarithm function in the solution affects the accuracy of the approximation.
5 Try to apply the ideas of this section to approximating the solution of the two point boundary value problemCan we get a tridiagonal system that Algorithm 2.6 can be applied to? Hint: Consider some of the approximations from §2.2; use the most accurate ones that can be easily used.Solution: Use the approximationto construct the approximation. Using the approximationwill result in a loss of accuracy.
6 Solve the two‐point boundary value problem problemusing a range of mesh sizes, starting with , and going as far as . Comment on your results.Solution: For larger values of , the approximate solution is erratic and wildly oscillatory. It isn't until about or so that the solution begins to settle down.
7 Generalize the solution of the two‐point boundary value problem to the case where and . Apply this to the solution of the problemwhich has exact solutionSolve this for a range of values of the mesh. Do we get the expected accuracy?
8 Consider the problem of determining the deflection of a thin beam, supported at both ends, due to a uniform load being placed along the beam. In one simple model, the deflection as a function of position along the beam satisfies the boundary value problemHere is a constant that depends on the material properties of the beam, is the length of the beam, and depends on the material properties of the beam as well as the size of the load placed on the beam. For a six‐foot‐long beam, with and , what is the maximum deflection of the beam? Use a fine enough grid that you can be confident of the accuracy of your results. Note that this problem is slightly more general than our example (2.28)‐(2.29); you will have to adapt our method to this more general case.Solution: See Figure 2.6.
9 Repeat the above problem, except this time use a three‐foot‐long beam. How much does the maximum deflection change, and is it larger or smaller?
10 Repeat the beam problem again, but this time use a 12‐foot‐long beam.
11 Try to apply the ideas of this section to the solution of the nonlinear boundary value problem defined byFigure 2.6 Exact solution for Exercise 2.7.8.Write out the systems of equations for the specific case of . What goes wrong? Why can't we proceed with the approximate solution?Solution: The problem is that the system of equations that is produced by the discretization process is nonlinear, so