计算物理shooting for boundary value porblem打靶法解决边界条件微分方程问题

3 6 Boundary value and eigenvalue problems Another class of problems in physics requires knowledge of solutions to differential equations with values of physical quantities or their derivatives given at the boundaries of a specified region This applies to the solution of the Poisson equation with a given charge distribution and known boundary values of the electrostatic potential or of the stationary Schrodinger equation with a given potential and boundary conditions A typical boundary value problem in physics us usually given in a second order differential equation 3 43 xuufu Where is a function of x and are the first order and second order u u u derivatives of with respect to x and is a function of and x u xuufu u Either u or is given at each boundary Note that one can always choose a u coordinate system so that the boundaries of the system are at x 0 and x 1 without losing any generality For example if the actual boundaries are at x and for a given 1 x 2 xx problem one call always bring them back to with a transformation1 0 andxx 3 44 121 xxxxx For problems in one dimension one can have a total of four possible types of boundary conditions 1 and 0 0 uu 1 1 uu 2 and 0 0 uu 1 1 vu 3 and 0 0 vu 1 1 uu 4 and 0 0 vu 1 1 vu The boundary value problem is more difficult to solve than the similar initial value problem with the same differential equation For example if we want to solve 3 45 xuufu with the initial conditions u 0 and we can first transform the 0 u 0 0 vu differential equation into a set of first order differential equations with a redefinition of the first order derivative into a new variable The solution will follow if we adopt one of the algorithms discussed earlier in this chapter However for the boundary value problem we know only or 0 u which is not sufficient to start any algorithms for the initial value 0 u problems without ally further work Typical eigenvalue problems are even more complicated because at least one more parameter that is the eigenvalue is involved 1n the equation For example xuufu with a set of given boundary conditions defines an eigenvalue problem Here the eigenvalue can have only some selected values in order to yield acceptable solutions of the equation under the given boundary conditions Let us take the longitudinal vibrations along an elastic rod as an example here The equation describing the stationary solution of elastic waves is given by uku 2 where u x is the displacement from the equilibrium at z and the allowed values of the wavevector k are the eigenvalues of the problem The wavevector in the equation is related to the sound speed c along the rod and the allowed angular frequency by the dispersion relation ck If both ends x 0 x 1 of the rod are fixed the boundary conditions are then If one end x O is fixed and the other end x 1 is free the 0 1 0 uu boundary conditions are then u 0 0 and For this problem one can 0 1 u obtain an analytical solution For example if both ends off the rod are fixed one will have the eigenfunctions 3 49 xkxu nn sin2 as possible solutions of the differential equation Here the eigenvalues are given by 3 50 22 nkn where n The complete solution of the longitudinal waves along 2 1 the elastic rod is given by a linear combination of all the eigen functions with their associated initial solutions 3 51 n ctin nn exuctxu Where are coefficients to be determined by the initial conditions We will n c come back to this problem in Chapter 6whenwe discuss the solutions of a partial differential equation 3 7 The shooting method A simple method for solving the boundary value problem of Eq 3 43 and the eigenvalue Problem of Eq 3 46 with a set of given boundary conditions is the so called shooting method We will first discuss how it works for the boundary value Problem and then generalize it to the eigenvalue problem We first convert the second order differential equation into two first order differential equations by defining and that is uy 1 2 uy 3 52 2 1 y dx dy 3 53 21 2 xyyf dx dy with a given set of boundary conditions To illustrate the method let us assume that the boundary conditions are and Any other types of 0 0 uu 1 1 uu Boundary conditions can be solved in a similar manner The key here is to make the problem look like an initial value problem by introducing an adjustable parameter the solution is then obtained by varying the parameter Since u 0 is given already we can make a guess for the first 0rder derivative at x 0 for example Here is the parameter to be adjusted For a specific 0 u we can integrate the equation to x 1 with any of the algorithms discussed earlier for the initial value problem Since the initial choice of can hardly be the actual derivative at x 0 the value of the function 1 resulting from the u integration with to x 1 would not be the same as The idea of the 0 u 1 u shooting method is to use one of the root search algorithms to find the appropriate that ensures within a given tolerance 0 1 1 uuF Let us take an actual numerical example to illustrate the scheme Assume that we want to solve the following differential equation 3 54 1 4 2 uu With the given boundary conditions and u 1 1 we can define the new 0 0 u variables and then we have uy 1 2 uy 3 55 2 1 y dx dy r 3 56 1 4 1 2 2 y dx dy Now assume that this equation set has the initial values and 0 0 1 y 0 2 y Here is a parameter to be adjusted in order to have with 1 1 uF being the tolerance of the solution We can combine the secant method for root search and the fourth order Runge Kutta method for the initial value problem to solve the above equation The following program is an implementation of such a combined approach to the boundary value problem defined in Eq 3 54 or Eqs 3 55 and 3 56 with the given boundary conditions C The boundary value problem solved from the above program can also be solved exactly with all analytical solution 3 57 1 2 sin2 2 cos xx xu One can easily check that the above expression does satisfy the equation and the boundary conditions Here we plot both the numerical result obtained from the shooting method and the analytical solution in Fig 3 4 As one can see the shooting method provides a very accurate solution of the boundary value problem It is also a very general method for the boundary value problem Fig 3 4 Tire numerical solution of the boundary value problem of Eq 3 54 by the Shooting method compared with the analytical solution solid line of the same Problem Boundary value problems with other types of boundary conditions can be solved in a similar manner For example if and u 1 are given we can 0 0 vu 1 u make a guess on and then integrate the equation set of and to x 1 0 u 1 y 1 y the root to be sought is from Here is the numerical 0 1 1 uuF 1 u result of the equation with If given the equation 0 u 1 1 vu is solved Instead 0 1 1 vuF One can apply the shooting method to solve the eigenvalue problem too The parameter to be adjusted in the eigenvalue problem is no longer a parameter introduced but the eigenvalue of the problem For example if and are given we can integrate the equation with 0 0 uu 1 0 uu 0 u with as a small quantity Then we search for the root of 0 1 1 uuF When we obtain an approximate eigenvalue and the F corresponding eigenvector from the normalized solution of The xu introduced parameter d is not relevant here because it will be automatically modified to be the first order derivative here solution is normalized In other words one can choose the first order derivative at the boundary arbitrarily and it will not affect the results as long as the solutions are made orthonormal orthogonal 3 8 Linear equations and the Sturm Liouville problem Many eigenvalue or boundary value problems are in the form of linear equations 3 59 1 0 xbuxauxu Where d x q x and s x are functions of x Assume that the boundary conditions are and If all d x q x and s x are 0 0 uu 1 1 uu smooth one can solve the equation with the shooting method developed in the preceding section In fact one can show that an extensive search for the parameter from 0 is not necessary in this case because 0 1 1 uuF of the superposition principle of the linear equation Any linear combination of the solutions is still a solution of the equation one needs only two trial solutions and with and being two different parameters The correct 0 xu 1 xu 0 1 solution of the equation is given by 01 uxauxbux with a and b determined from u 0 and u 1 0 u 1 u Note that So we have 0 0 0 0 10 uuuu a b 1 36 0 3 61 1 1 1 10 ubuau which can easily be solved to give 3 62 1 1 1 01 11 uu uu a 3 63 1 1 1 01 01 uu uu b With a and b given from the above equation we have the solution of the differential equation from Eq 3 59 An important class of linear equations in physics is referred to as the Sturm Liouville problem defined by 3 64 xsxuxqxuxp which has the first order derivative term combined with the second order derivative term P x q x and s x are the coefficient functions of x For most actual problems s x 0 and q x with being the xxr eigenvalue of the equation r x and w x are the redefined coefficient functions The Legendre equation the Bessel equation and the related equations in physics are examples of the Sturm Liouville problem Our goal here is to construct an accurate algorithm that can integrate Sturm Liouville equation that is Eq 3 64 In Chapter 2 we obtained the three point 3 65 6 2 4 3 2 11 1 hO uh u h uu nnn And 3 36 12 2 4 4 2 2 11 1 hO uh u h uuu n n nnn Now if we multiply Eq 3 65 by and Eq 3 66 by and add them together n p n p we have 3 67 2 12 4 3 4 2 201 hOupup h uppp nnnnnnnn If we replace the first term on the right hand side with and drop the nnn uqs second term we obtain the simplest numerical algorithm for the Sturm LiouvlIle equation 3 68 24 2 2 2 1 1 nnnnnnnnnnn uqshupuhppuhpp which is accurate up to O One has to be very careful with this apparent local 4 h accuracy in the algorithm In reality because of the repeated use of the three point formula the local accuracy delivered is usually lower than O Before 4 h we discuss how to improve the accuracy of this algorithm let us illustrate it with all example The Legendre equation is given by 3 69 0 1 1 2 ull dx du x dx d With and x The solutions of the Legendre equation are the 1 0 l 1 1 Legendre polynomials Let us assume that we do not know the value of l xpl but know the first two points of x then we can treat the problem as an xpl eigenvalue problem The following program is an implementation of the simplest algorithm for the Sturm Liouville equation in combination with the bisection method for the root search to solve for the eigenvalue 1 of the Legendre l equation The eigenvalue obtained from the above program is 1 000091 which contains an error of 9 in comparison with the exact result of The error comes 5 10 1 l from both the rounding error and the inaccuracy of the algorithm If we want to have higher accuracy in the algorithm for the Sturm Liouvile Equation we can differentiate the equation twice Then we have 3 70 2 3 2 3 3 3 4 quuququpupupsuppu Where on the right hand side can be replaced with 3 u 2 1 3 quuqupups p u Which is the result of differentiating the Sturm Liouville equation once If we combine Eqs 3 67 3 70 and 3 71 we obtain a better algorithm 3 72 6 1111 hOducucuc nnnnnnn where and are given by 1 n c n c 1 n c n d 3 73 nnnnnnn nnnnnnnn ppphpqphqh phpphphqhhppc 2 4621224 3 3 3 3 32 2 22 1 3 74 4 4 2 22 2 2 12 82048 nnnn nnnnn qhpqph phpphqhpc 3 75 nnnnnnn nnnnnnnn ppphpqphqh phpphphqhhppc 2 4621224 3 3 3 3 32 2 22 1 3 76 nnnnn psphshhd 2224 4 42 which can be evaluated easily if p x q x and s x are explicitly given In the case where some the derivatives that are needed are not easily obtained analytically one can evaluate them numerically In order to maintain the high accuracy of the algorithm one needs to use comparable numerical formulas For the special case with p x 1 the above coefficients reduce to much simpler forms Without sacrificing the high accuracy of the algorithm we can apply the three point formulas to the first order and second order derivatives of q x and s x Then we have 3 77 2 1 2 1 12 1 11 2 1 nnnn qqq h c 3 78 10 1 5 4 10 1 6 5 2 11 2 nnnn qqq h c 3 79 2 1 2 1 12 1 11 2 1 nnnn qqq h c 3 80 10 12 11 2 nnnn sss h d which are slightly different from the Numerov algorithm which is an extremely accurate scheme for linear differential equations without the first order derivative term that is Eq 3 58 with d x 0 or the Sturm Liouville equation with p x 1 Many equations in physics have this form for example the Poisson equation with Spherical symmetry or the one dimensional Schrodinger equation The Numerov al