Heath - Scientific Computing (523150), страница 97
Текст из файла (страница 97)
For each type of problem in partsa, b, and c to follow, describe how you coulduse this routine to solve it. In each case, youranswer should address the following points:1. What is the function f for the ODE subproblem?2. How would you obtain the necessary initialconditions?3. What special properties, if any, would theODE subproblem have that would affectthe choice of ODE method?(a) Compute the definite integralZg(t) dt.(b) Solve the two-point boundary value problemy 00 = y 2 + t, 0 ≤ t ≤ 1,with boundary conditionsand y(1) = 1.(c) Solve the heat equationut = cuxx ,......................................................................................................................................................................................................................................................................................................................................................1•0•••0•1x11.3 Give examples to show that neither consistency nor stability alone is sufficient to guarantee convergence of a finite difference schemefor solving a partial differential equation numerically.bay(0) = 0a second-order accurate, centered finite difference scheme to compute the approximate valueof the solution at the center of the square.y0 ≤ x ≤ 1,t ≥ 0,with initial conditionu(0, x) = g(x)and boundary conditionsu(t, 0) = 0,u(t, 1) = 0.11.2 Consider a finite difference solution ofthe Poisson equation uxx + uyy = x + y onthe unit square using the boundary conditionsand mesh points shown in the drawing.
Use11.4 Draw pictures to illustrate the nonzeropattern of the matrix resulting from a finitedifference discretization of the Laplace equation on a d-dimensional grid, with k grid pointsin each dimension, for d = 1, 2, and 3, as described at the end of Section 11.3.1. Use avalue of k that is large enough to show thegeneral pattern clearly. In each case, what arethe numerical values of the nonzero entries?11.5 (a) For the matrix of Example 11.5, verify that the nonzero pattern of the Choleskyfactor, including fill, is as shown in Figs. 11.6–11.8.(b) Verify that the reorderings given inFigs. 11.7–11.8 are indeed minimum degreeand nested dissection orderings, respectively.In each case, write out the permutation matrix P corresponding to the reordering, suchthat the reordered matrix is given by P AP T ,where A is the matrix in the original naturalordering.11.6 Prove that the Jacobi iterative methodfor solving the linear system Ax = b con-362CHAPTER 11.
PARTIAL DIFFERENTIAL EQUATIONSverges if the matrix A is diagonally dominantby rows. (Hint: Use the ∞-norm.)11.7 Prove that the SOR method diverges ifω does not lie in the interval (0, 2).11.8 Prove that the successive A-orthogonalsearch directions generated by the conjugategradient method satisfy a three-term recurrence, so that each new gradient need be orthogonalized only against the previous two.11.9 Show that the subspace spanned by thefirst m search directions in the conjugate gradi-ent method is the same as the Krylov subspacegenerated by the sequence r0 , Ar0 , A2 r0 , .
. .,Am−1 r0 .11.10 On the basis of the background factsand motivation outlined in Section 11.5.5, fillin the details of the derivation of the conjugate gradient algorithm stated there. Youwill probably find it easiest first to derive astraightforward version of the algorithm, thentry to improve its efficiency by reducing thenumber of matrix-vector multiplications andother arithmetic required.Computer Problems11.1 Consider the heat equationut = uxx ,0 ≤ x ≤ 1,t ≥ 0,with boundary conditionsu(t, 0) = 0,u(t, 1) = 0,and initial condition2xif 0 ≤ x ≤ 0.5u(0, x) =.2 − 2x if 0.5 ≤ x ≤ 1(a) Using the full discretization given in Example 11.1, with ∆x = 0.05 and ∆t = 0.0012,evolve the solution from t = 0 to t = 0.06.Plot the solution at the initial time, the finaltime, and periodically in between (say, everyten time steps or so).
If you use interactivegraphics, you may prefer to plot the solutionat every time step, giving the effect of an animated movie.(b) Repeat part a, but use a time step of∆t = 0.0013. Plot the solution in the samemanner as before. Can you explain the difference in results?(c) Solve the same equation again, thistime using the implicit method given in Section 11.2.4 based on the backward Eulermethod. Again use ∆x = 0.05, but try a muchlarger time step, say, ∆t = 0.005 to advancethe solution from t = 0 to t = 0.06. Plot thesolution at each time step.
How do the resultscompare with those in parts a and b? Explain.(d ) Solve the same equation again, this timeusing the Crank-Nicolson method given in Section 11.2.4. Again use ∆x = 0.05 and ∆t =0.005 to advance the solution from t = 0 tot = 0.06.
Plot the solution at each time step.How do the results compare with your previousresults? Explain.(e) Form the semidiscrete system for this equation given in Section 11.2.1, using a finite difference spatial discretization with ∆x = 0.05,and use an ODE solver to integrate it fromt = 0 to t = 0.06. Plot the solution as before. How do the results compare with yourprevious results?11.2 (a) Use the method of lines and the ODEsolver of your choice to solve the heat equationut = uxx ,0 ≤ x ≤ 1,t ≥ 0,with Dirichlet boundary conditionsu(t, 0) = 0,u(t, 1) = 0,and initial conditionu(0, x) = sin(πx).Integrate from t = 0 to t = 0.1.
Plotthe computed solution, preferably as a threedimensional surface over the (t, x) plane. Ifyou do not have three-dimensional plotting capability, plot the solution as a function of x fora few values of t, including the initial and final times. Determine the maximum error inthe computed solution by comparing with theexact solutionu(t, x) = exp(−π 2 t) sin(πx).COMPUTER PROBLEMS363Experiment with various spatial mesh sizes h,and try to characterize the error as a functionof h.
On a log-log scale, plot the maximumerror as a function of h.(b) Repeat part a, but this time with Neumannboundary conditionsux (t, 0) = 0,ux (t, 1) = 0,and initial conditionu(0, x) = cos(πx),and compare with the exact solutionu(t, x) = exp(−π 2 t) cos(πx).11.3 Use the method of lines and the ODEsolver of your choice to solve the wave equationutt = uxx ,0 ≤ x ≤ 1,t ≥ 0,with boundary conditionsu(t, 0) = 0,u(t, 1) = 0,with boundary conditionu(t, 0) = 1,and initial conditionu(0, x) = 0.Integrate from t = 0 to t = 2, and plotthe computed solution, preferably as a threedimensional surface over the (t, x) plane. Ifyou do not have three-dimensional plotting capability, plot the solution at x = 1 as a function of t.
Try both one-sided (upwind) andcentered finite difference schemes for the spatial discretization. The exact solution is a stepfunction of height 1 moving to the right withvelocity 1. Does either difference scheme comeclose to this ideal? How would you describethe difference between the computed solutionsgiven by the two schemes? Which computedsolution is smoother? Which is more accurate?11.5 Use the method of lines and the ODEsolver of your choice to solve the nonlinearBurgers’ equationut = −uux + vuxx ,and initial conditions0 ≤ x ≤ 1,t ≥ 0,with boundary conditionsu(0, x) = sin(πx),ut (0, x) = 0.Integrate from t = 0 to t = 2.Plotthe computed solution, preferably as a threedimensional surface over the (t, x) plane.
Ifyou do not have three-dimensional plotting capability, plot the solution as a function of x fora few values of t, including the initial and final times. Determine the maximum error inthe computed solution by comparing with theexact solutionu(t, x) = cos(πt) sin(πx).Note that the solution is periodic with period2, so the solution should be the same at t = 0and t = 2. Is this true for your computed solution? Experiment with various spatial meshsizes h, and try to characterize the error asa function of h.
On a log-log scale, plot themaximum error as a function of h.11.4 Use the method of lines and the ODEsolver of your choice to solve the first-orderhyperbolic equationut = −ux ,0 ≤ x ≤ 1,t ≥ 0,u(t, x) = 1/(1 + exp(x/(2v) − t/(4v)))for x = 0 and x = 1, and initial conditionu(0, x) = 1/(1 + exp(x/(2v))),where v is a scalar parameter representing viscosity. Note that this equation has both a nonlinear convective (hyperbolic) term, −uux , anda linear diffusion (parabolic) term, vuxx , andthe balance between them is determined by therelative value of v.Integrate from t = 0 to t = 2, and plotthe computed solution, preferably as a threedimensional surface over the (t, x) plane.
Ifyou do not have three-dimensional plotting capability, plot the solution at x = 1 as a function of t. Experiment with various values forv, say, v = 1, 0.5, 0.1, and 0.01. Try both onesided (upwind) and centered finite differenceschemes for the spatial discretization. Compare your solutions with the exact solution,which is the function used previously to define the boundary conditions. Compare theaccuracy and smoothness of the solutions computed by the two schemes as v varies.364CHAPTER 11. PARTIAL DIFFERENTIAL EQUATIONS11.6 Use the standard five-point finite difference discretization to solve the Poisson equation−uxx − uyy = 2on the L-shaped region in the accompanyingdiagram, with boundary conditions u(x, y) = 0on all boundaries.......................................................1.0...
. . . . .......... . . . . ............ . . . . ........ . . . . ..........0.5 .................................................... .. .. .. .. ........... . . . . . . . . ......... . . . . . . . . ....... . . . . . . . . ........0.00.00.51.0(a) Choose a mesh size h (for example, h = 0.1in the diagram) and set up the appropriatematrix and right-hand-side vector. Solve theresulting linear system by a direct methodand plot the resulting solution, preferably as athree-dimensional surface over the (x, y) plane.Experiment with various values for the meshsize h and observe the effect on the solution.(b) Repeat part a, but this time use the SORmethod to solve the linear system. Note thatwith an iterative method you need not set upthe matrix explicitly but can work directlywith the approximate solution at each meshpoint, updating the values at each iteration.For a given mesh size h, experiment with various values for the SOR parameter ω in therange 0 < ω < 2.