Heath - Scientific Computing (523150), страница 86
Текст из файла (страница 86)
Describe the nonzero patternof the matrix of the resulting system of linearalgebraic equations.10.13 Suppose you are using the shootingmethod to solve a two-point boundary valueproblem for an ODE on an interval [a, b]. Ifthe ODE in question is unstable on some portion of the interval, then the resulting sequenceof initial value problems may be very sensitiveto initial conditions, making it difficult to hitthe required boundary condition.(a) How could you cope with such illconditioning?(b) How would this affect the nonlinear algebraic equation to be solved?10.14 In solving a two-point boundary valueproblem for a second-order ODE numerically,does the approximate solution produced by finite element collocation at a finite set of n discrete points always agree with the exact solution at those n points?EXERCISES321Exercises10.1 Consider the two-point boundary valueproblem for the second-order ODEy 00 = y 3 + t,a ≤ t ≤ b,with boundary conditionsy(a) = α,y(b) = β.To use the shooting method to solve this problem, one needs a starting guess for the initialslope y 0 (a).
One way to obtain such a startingguess for the initial slope is, in effect, to doa “preliminary shooting” in which we take asingle step of Euler’s method with h = b − a.(a) Using this approach, write out the resulting algebraic equation for the initial slope.(b) What starting value for the initial sloperesults from this approach?10.2 Suppose that the altitude of the trajectory of a projectile is described by the secondorder ordinary differential equation y 00 = −4.Suppose that the projectile is fired from position t = 0 and height y(0) = 1 and is tostrike a target at position t = 1, also of heighty(1) = 1.(a) Solve this boundary value problem by theshooting method:1. To determine the initial slope at t = 0 required to hit the desired target at t = 1,use the trapezoid rule with stepsize h = 1to derive a system of two equations for theunknown initial slope s0 = y 0 (0) and finalslope s1 = y 0 (1).2.
What are the resulting values for the initialand final slopes?3. Using the initial slope just determined anda stepsize of h = 0.5, use the trapezoidrule once again to compute the approximate height of the projectile at t = 0.5.(b) Solve the same boundary value problem again, this time using a finite differencemethod with h = 0.5. What is the resulting approximate height of the projectile at thepoint t = 0.5?(c) Solve the same boundary value problemonce again, this time using collocation at thepoint t = 0.5, together with the boundary values, to determine a quadratic polynomial u(t)approximating the solution.
What is the resulting approximate height of the projectile atthe point t = 0.5?Computer Problems10.1 Solve the two-point boundary valueproblemy 00 = 10y 3 + 3y + t2 ,0 ≤ t ≤ 1,with boundary conditionsy(0) = 0,y(1) = 1,by each of the following methods.(a) Shooting method. Use a one-dimensionalnonlinear equation solver to find an initialslope y 0 (0) such that the solution of the resulting initial value problem hits the target valuefor y(1).
Solve each required initial value problem using a library ODE solver or one of yourown design. Plot the sequence of solutions youobtain.(b) Finite difference method. Divide the giveninterval 0 ≤ t ≤ 1 into n+1 equal subintervals,0 = t0 < t1 < · · · < tn < tn+1 = 1,with each subinterval of length h = 1/(n + 1).Let yi , i = 1, . .
. , n, represent the approximatesolution values at the n interior points. Obtaina system of n algebraic equations for the yi byreplacing the second derivative in the differential equation by the finite difference approximationyi00 (t) ≈yi+1 − 2yi + yi−1,h2i = 1, . . . , n. Use a library routine, or one ofyour own design, to solve the resulting systemof nonlinear equations. A reasonable starting guess for the nonlinear solver is a straight322CHAPTER 10. BOUNDARY VALUE PROBLEMS FOR ODESline between the boundary values. Plot the sequences of solutions you obtain for n = 1, 3,7, and 15.(c) Collocation method. Divide the given interval 0 ≤ t ≤ 1 into n − 1 equal subintervals,0 = t1 < t2 < · · · < tn−1 < tn = 1,with each subinterval of length h = 1/(n − 1).Take the approximate solution u(t) to be apolynomial of degree n − 1.
Forcing u(t) tosatisfy the boundary conditions at the endpoints and to satisfy the ODE at the n − 2interior points yields a system of n equationsthat determine the n coefficients of the polynomial u(t). Use a library routine, or one of yourown design, to solve this system of nonlinearalgebraic equations.
The resulting polynomialcan then be evaluated at any point in the interval to obtain an approximate solution value atthat point. Print the polynomial coefficientsand plot the solutions you obtain for n = 3, 4,5, and 6.10.2 Solve the two-point boundary valueproblemy 00 = −(1 + ey ),0 ≤ t ≤ 1,with boundary conditions(a) Use both the shooting and finite differencemethods to determine the curve of the ropewhen the boundary conditions are 00.750 0 y(0) = , y(1) = .0011These conditions correspond to a slack rope.Plot the solution curve you obtain for eachmethod.(b) Use both the shooting and finite differencemethods to determine the curve of the ropewhen the boundary conditions are 00.850 0.50 y(0) = , y(1) = .0011These conditions correspond to a taut rope.Plot the solution curve you obtain for eachmethod.10.4 The deflection of a horizontal beam supported at both ends and subjected to axialand transverse loads can be described by thesecond-order ODEy 00 = λ(−t2 − 1)y,with boundary conditionsy(−1) = 0,y(0) = 0,y(1) = 0.y(1) = 1,using each of the methods in the previous exercise.10.3 The curve of a hanging rope is describedby the system of ODEsy10y20y30y40−1 ≤ t ≤ 1,= cos(y3 ),= sin(y3 ),= (cos(y3 ) − sin(y3 )| sin(y3 )|)/y4 ,= sin(y3 ) − cos(y3 )| cos(y3 )|,where y1 (t) and y2 (t) are the horizontal andvertical coordinates of the rope, y3 (t) is theangle between the tangent to the rope and thehorizontal axis, y4 (t) is the tension in the rope,and the variable t is the arc length along therope, with the length of the rope normalizedso that 0 ≤ t ≤ 1.The eigenvalues and eigenfunctions for thistwo-point boundary value problem determinethe frequencies and modes of vibration of thebeam.
Use a finite difference discretizationof the ODE to derive an algebraic eigenvalueproblem whose eigenvalues and eigenvectorsapproximate those of the ODE, then computethe eigenvalues and eigenvectors using a library routine (see Section 4.6). Experimentwith various mesh sizes and observe how theeigenvalues behave.10.5 The time-independentequation in one dimension,Schrödinger−ψ 00 (x) + V (x)ψ(x) = Eψ(x),where we have chosen units so that the quantities are dimensionless, describes the wave function ψ of a particle of energy E subject to apotential function V .
The square of the waveCOMPUTER PROBLEMS323function, |ψ(x)|2 , can be interpreted as theprobability of finding the particle at positionx.Assume that the particle is confined to a onedimensional box, say, the interval [0, 1], withinwhich it can move freely. Thus, the potential is zero within the unit interval and infiniteelsewhere. Since there is zero probability offinding the particle outside the box, the wavefunction must be zero at its boundaries. Thus,we have an eigenvalue problem for the secondorder ODE00−ψ (x) = Eψ(x),0 ≤ x ≤ 1,subject to the boundary conditionsψ(0) = 0and ψ(1) = 0.Note that the discrete eigenvalues E are theonly energy levels permitted; this feature givesquantum mechanics its name.Use a finite difference discretization of theODE to derive an algebraic eigenvalue problem whose eigenvalues and eigenvectors approximate those of the ODE, then computethe eigenvalues and eigenvectors using a library routine (see Section 4.6).
Experimentwith various mesh sizes and observe how theeigenvalues behave.An analytical solution to this problem is easilyobtained, which gives the eigenvaluesEk = k 2 π 2and corresponding eigenfunctionsψk (x) = sin(kπx),k = 1, 2, . . . .How do your computed eigenvalues and eigenvectors compare with these analytical values asthe mesh size of your discretization decreases?Try to characterize the error as a function ofthe mesh size.Note that a nonzero potential V would notseriously complicate the numerical solution ofthe Schrödinger equation, but would generallymake an analytical solution much more difficult to obtain.324CHAPTER 10.
BOUNDARY VALUE PROBLEMS FOR ODESChapter 11Partial Differential Equations11.1Partial Differential EquationsWe turn now to partial differential equations (PDEs), where many of the numerical techniques we saw for ODEs, both initial and boundary value problems, are also applicable.The situation is more complicated with PDEs, however, because there are additional independent variables, typically one or more space dimensions and possibly a time dimensionas well. Additional dimensions significantly increase computational complexity. Problemformulation also becomes more complex than for ODEs, as we can have a pure initial valueproblem, a pure boundary value problem, or a mixture of the two. Moreover, the equationand boundary data may be defined over an irregular domain in space.First, we establish some notation.
For simplicity, we will deal only with single PDEs(as opposed to systems of several PDEs) with only two independent variables (either twospace variables, which we denote by x and y, or one space and one time variable, which wedenote by x and t). In a more general setting, there could be any number of dimensionsand any number of equations in a coupled system of PDEs. We denote by u the unknownsolution function to be determined and its partial derivatives with respect to the independentvariables by appropriate subscripts: ux = ∂u/∂x, uxy = ∂ 2 u/∂x∂y, etc.11.1.1Classification of Partial Differential EquationsPartial differential equations are classified by the value of the discriminant, b2 − 4ac, in thegeneral linear, two-dimensional, second-order PDEauxx + buxy + cuyy + dux + euy + f u + g = 0,b2 − 4ac > 0:b2 − 4ac = 0:b2 − 4ac < 0:hyperbolic,parabolic,elliptic.In practice, this classification is not always so clean and simple.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
