Erdelyi note on orthogonal polys

Erdelyi note on orthogonal polys PhL 1.11.05 The time has now come to clean up this subject. I think Erdelyi's presentation looks good, and I can supplement it as needed. I did take some earlier notes while reading some Jim Ball papers. In these many pages I review in great detail the opening "general" pages of this chapter, and I skip the details of the specific functions described, such as Jacobi's etc. The main function of this document is to have a repository for the derivations of some quite messy formulas that I have never derived before. And in going through all this, I got much more familiar with the notation used, especially for the many constants. In a separate document I give an "executive overview" of what I learned in this chapter, and what I did while reading it. This document is mostly just fine details. Things are done in book section order. 10.1 Systems of Orthogonal Functions Given an interval (a,b) and a non-negative weight function w(x), we can define a scalar product <f|g> and thus a Hilbert space. The weight function has to be non-negative so that <f|f> = ||f||2 will be positive for any function f. The integral here is associated by Erdelyi with the name Stieltjes. Comment: Web says that the usual integrals we deal with are Riemann integrals, and that Lebesgue integrals are not part of the ordinary world. The Stieltjes integral derives from this strange sum idea: i f(xi) [ g(xi+1) - g(xi) ] generalizing i f(xi) [ xi+1 - xi ] where xi* is some arbitrary point in the interval (xi, xi+1). In the limit, this is f(x)dg , and so you can think of this as a generalization of f(x)dx. If things are reasonable, you write dg = g'(x) dx = w(x) dx. Now suppose you use g(x) = wi (x - xi) + wi+1 (x - xi+1) + wi+2(x - xi+2) + ... which is just a sum of steps of height wi. and width [ xi+1 - xi ]. Then w(x) = g'(x) = wi (x - xi) + wi+1 (x - xi+1) + wi+2 (x - xi+2) .. Then you get f(x)dg = i wi f(xi). So the claim is that the Stieltjes idea encompasses both continuous and discrete scalar products. I think in the second case, g(x) is a distribution rather than a function. I think the scalar product as defined above gives the L2 Hilbert Space of square integrable functions, and we can take the powers of x as a basis for a space of infinite dimension, or we can stop at some power and have a finite dimensional Hilbert Space. The basis functions are always linearly independent, perhaps not orthogonal. You can orthonormalize using G-S and this is mentioned bottom of page 154 where you start with functions which are not orthonormal, and you create functions which are. Lemma 1.1 Consider the following weighted and summed determinant expression Q = i fi We claim that we can move the i fi into the bottom row and let it act on the terms there, so we have Q = Proof: To see this, just write the thing out as follows Q = i fi [ m a bm(i)] = m a [ i fi bm(i)] = as shown in second expression above. Corollary 1.1 We could also state this replacing i fi = dx f(x) in which case we get: Q = dx f(x) = These things help us to understand the comments on page 155. Suppose instead of doing the G-S process as described bottom of 154 and top of 155, we construct the n(x) as shown in (10), where we are using a 2x2 matrix inside the det for 1(x) and a 3x3 for 2(x) and an (n+1) x (n+1) for n(x). WARNING: Notice that this author uses and as different things. These are both the same Greek letter "phi", so that is quite confusing. The are given by (10), whereas the are going to be normalized! At this point, each statement really requires work, so here we go: Claim 1.1: The n(x) given by (10) is orthogonal to k for k < n. Proof: If we compute (k, n) using our Corollary 1.1 above, the bottom row of the determinant becomes the following: (k, 0), (k, 1), (k, 2).... . But this duplicates row k of the determinant, as long as we restrict to k = 0, 1, 2...n-1. QED. Claim 1.2: The n(x) given by (10) is orthogonal to k for k < n. Proof: k for k < n is given by a lincom of k for k < n. This is what (10) says. But according to Claim 1, we know that n(x) is orthogonal to all these k for k < n. Thus, n(x) is orthogonal to any lincom you can make of the k for k < n, and therefore n(x) is orthogonal to k for k < n. Definition 1.1: Define ajk = (j, k) = dx w(x) j(x)k(x) Definition 1.2: The Gram Determinant is Gn = det(ajk) where we let j,k range 0 to n. Thus, Gn is the determinant of a matrix of dimension n+1. Claim 1.3: This matrix Gn is the upper left n+1 size subdeterminant of the matrix for n+1. It is therefore the cofactor of the lower right element which is n+1. Proof: Cross out the last row and column of the matrix for n+1, and look at the subdeterminant so formed. The top row has a00 through a0n . The first column has a00 through an0. This matrix is exactly Gn. Corollary 1.3: In the matrix expression (10) for n(x), the coefficient of n(x) is Gn-1. Later when we choose the specific form n(x) = xn , we will have that Gn-1is the coefficient of the leading power xn in the matrix expression for n(x). Claim 1.4: The expression shown with the coefficients is a quadratic form. Proof: Write it in the form jk jT ajk k = T a where matrix a is n+1 x n+1. Claim 1.5: Gn is a discriminant of this quadratic form. Proof: According to M&M page 323, this matrix and all submatrices are called discriminants of a quadratic form, to Gn is the largest one. Theorem 1.1: The || k ||2 = Gk-1Gk. From this, it follows that || k ||2 = 1 as shown in (11). Proof: This took me a while to see, but now it's obvious. Consider: || k ||2 = (k, k) Now look at the expansion (10) for k = 00 + 11 + ... + kk where the are the messy cofactors which arise from doing the Laplace expansion of the determinant along the bottom row. However, we know that k = Gk-1 from Claim 3 above. We also know from Claim 1 that n is orthogonal to the lower k. Therefore: || k ||2 = (k, k) = (00 + 11 + ... + Gk-1k, k) = Gk-1 (k, k). But, consider (k, k) which we can form by applying (k, to expansion (1) for k, using the idea of Corollary 1.1 above. We then find that (k, k) = Gk ! That is, we get the entire regular matrix axy which is of dimension k+1. Thus we have shown that || k ||2 = Gk-1 Gk. QED. Comment: Expansion (10) says that 0 = 0. Our equation (11) then says that || 0 ||2 = || 0 ||2 /(G0G-1) We take G-1 = 1 as they suggest, and G0 = a00 = (0,0) = (0,0) , so things work right for the 0 case. That is to say, || 0 ||2 = || 0 ||2/ (0,0) = 1. I am going to ignore equations (12) and (13) unless I find they are needed later! Theorem 1.2: The interval (a,b) and the weight w(x) determine a unique set of orthonormal functions. The dimension N has not really been mentioned. But suppose we start with the non-orthogonal powers xn, then we can just start orthonormlizing them as far as we want to go. I think I downtick Erdelyi here for being vague about Hilbert Space and the dimension of that space, and so on. Maybe yet to come. He ends this section by commenting that there are a few ways people decide to normalize functions. You might want to have them be monic, or to have a specified value at some point. 10.2 The Approximation Problem: Fitting a function with a lincom of orthog polys. The question brought up here is this: We want to "fit" a general function f(x) as a lincom of ortho polys with some coefficients. How do we select the coefficients to make a best least squares fit? This takes us back to the nice discussion in Scheid on pages 242-243, Chapter 21. If you are trying to best fit a vector f in a Hilbert space with a vector q in a subspace (eg, in a space of some finite number of ortho polys for your basis), then the best you can do is get to the point p in the subspace which is closest to f, as expressed in terms of the norm in the space, that is, minimize || f - q ||. This point of closest approach to f is reached ( in the case of an orthog basis) when you use the Fourier coefficients (f,k) as your lincom coefficients. Remember this whole business that, if you use ortho polys, then you improve your fit by adding more terms, and you never have to alter any of the earlier coefficients, and recall how this is not the case if you use non-orthog polys. Erdelyi shows explicitly how choosing the Fourier coefficient minimizes your error, which is those extra terms in the equation after (2). We then get (4) which implies (5) which is the famous Bessel's Inequality ai2 dx w(x) f2(x). He shows this with the sum going to infinity. But we know that for any "reasonable" basis functions (that is, closed basis), when you add up all the pieces you get equality, which is called Parseval's Formula (which I am familiar with in normal Fourier integral analysis), in the case f(x) is in L2. In this case, our poly approx to f(x) is "converging in the mean". I think this has to do with the space being complete and therefore really a Hilbert Space, and we associate the Riesz-Fischer theorem with this concept. We also have a Weierstrauss Theorem on this subject. When the interval (a,b) is finite, you always have things being "closed" and complete. Infinite intervals are more worrisome. Just in passing, recall from Scheid that if you attempt a least squared fit with a lincom of xk powers (which is a non orthog basis), you end up with those "normal equations", as on page 241 there. You have to solve these normal equations to find your coeffs, and as you add more terms, all the coeffs "move". 10.3 General Properties of Orthog Polys Our theorem is restated: interval and weight determine the polys uniquely apart from constant in each. Now in (1) we define the moments of a weight function, cn = dx w(x)xn. Equation (2) is trivially true for our case k(x) = xk. We can now apply our general theory of Section 10.1 to get Gn as shown, since we have that ajk = (j, k) = cj+k. So to get the subscript on c in (3), just add row and column. For now, I ignore the right equation shown in (3). Equation (4) is now just (10) of the first section for k(x) = xk , where we have added a particular normalizing factor kn / Gn-1. We are free to select kn for different purposes, BUT notice the following fact: kn is always the coefficient of xn in pn(x). This fact is used later in deriving the recursion formula. Theorem: pn(x) = (kn/Gn-1) * matrix where kn is the coefficient of xn. [ This is equation (4). ] Proof: if we look at only the matrix part, we know that coeff of xn is precisely Gn-1. This cancels the factor shown above, and you are then left with kn as the coefficient of xn. Fact: Notice that p0(x) = k0 and the entire matrix determinant equals 1. Think of k0x0 if you want. Theorem 3.1: Equation (6) for || pk ||2 is true with pn(x) as normalized in (4). Proof: In (4) the determinant is n and we shows above in Theorem 1.1 that || k ||2 = Gk-1Gk, therefore we have that || pk ||2 = (kk/Gk-1)2 || k ||2 = (kk/Gk-1)2Gk-1Gk = kk2 Gk/Gk-1 . Corollary 3.1. If we want to have || pk ||2= 1, we need to choose kk2 Gk/Gk-1 = 1 or kk2 = Gk-1/Gk. Theorem 3.2. All n zeros of pn(x) lie inside the interval (a,b), and between each consecutive pair of zeros of pn(x) lies exactly one zero of pn+1(x). The first claim is proved in a way I don't really like much, I am too lazy to write the whole thing up using my Erdelyi decoding method. The second claim is only claimed, refer to Szego he says. I think I have a simpler and clearer proof of the first claim elsewhere using Rolle's Theorem or some such. Theorem 3.3 (Recursion Formula) The recursion formula shown in (7) page 158 is true for the pn(x) as defined in (4), where the leading terms of pn(x) are knxn + kn'xn-1 , and where A,B,C are as given on top of page 159 (and these are restated below). Proof: This proof is long and tedious. You have your set of pn(x) and you construct this combination: fn(x) = pn+1(x) - An x pn(x) An = kn+1 / kn where the kn is the coefficient of the highest power in pn(x). I agree completely that fn(x) must be of degree n or less, because we have canceled out the xn+1 term by construction. fn(x) = { kn+1 xn+1 + lower terms} - ( kn+1 / kn) x { kn xn + lower terms} Now, I also agree that you can therefore write fn(x) as a lincom of the pn(x) as shown. For the moment, let's skip the claim that all i = 0 except 0 and 1 and come back to this claim very soon. Meanwhile, I do agree with the next line on page 159 and here is my proof: ( pn-1, fn ) = ( pn-1, pn+1 - An x pn) = - An ( pn-1, x pn) on the one hand, but then ( pn-1, fn ) = 1 ( pn-1, pn-1) on the other hand QED Suppose I try doing this with some lower power on p ? ( pn-2, fn ) = ( pn-2, pn+1 - An x pn) = - An ( pn-2, x pn) on the one hand, but then ( pn-2, fn ) = 2 ( pn-2, pn-2) on the other hand The trick is now to write ( pn-2, x pn) = (xpn-2, pn) , that is, we just move the x over to the other factor in the integral. But of course xpn-2 is a poly of degree n-1 or less which can be expanded in pr or order n-1 or less, and all these terms will be orthog to pn. Thus, we have 0 = (xpn-2, pn) = ( pn-2, x pn) => 2 = 0. This line or reasoning applies for n-3 and 3, and to n-4 and 4 and so on. What about 0 ? ( pn, fn ) = ( pn, pn+1 - An x pn) = - An ( pn, x pn) on the one hand, but then ( pn, fn ) = 0 ( pn, pn) on the other hand This one does not vanish! Therefore we have arrived at this point: pn+1 - An x pn = 0 pn + 1 pn-1 so right now we have the basic proof, we rewrite the above as pn+1 = pn ( 0 + Anx) + 1 pn-1 and now we just have to find the expressions for 0 and 1. The claim is that (a) Bn = 0 and (b) Cn = - 1 . Let's now prove these two facts. (b) From above we have that - An ( pn-1, x pn) = 1 ( pn-1, pn-1) = 1 hn-1. Now imitate the first expansion as follows: fn-1(x) = pn(x) - An-1 x pn-1(x) Therefore ( pn-1, x pn) = ( x pn-1, pn) = (pn - fn-1, pn) / An-1 = (pn, pn) / An-1 = hn/An-1 . This tells us that 1 hn-1 = - An * hn/An-1 or 1 = - (An/An-1) * (hn/hn-1) = - Cn. That is, 1 = - Cn. We note that there is a sign error in his statement of this fact on page 159, but he has the sign correct in (7). (a) From above we have that - An ( pn, x pn) = 0 ( pn, pn) = 0 hn. Now go back to an earlier line where we wrote out fn(x) = pn+1(x) - An x pn(x), fn(x) = { kn+1 xn+1 + k'n+1xn + lower terms} - ( kn+1 / kn) x { kn xn + k'n xn-1 + lower terms} where now we show the second term in each expansion. We already saw that the leading terms cancel, now what is the coefficient of the next term? fn(x) = [ k'n+1 - ( kn+1 / kn) k'n ] xn + lower terms. Now as our first step, write x pn = (pn+1 - fn)/An so that we have ( pn, x pn) = - (pn, fn)/An. Now in the remaining scalar product, only the leading xn term above survives (by the usual argument), so we have ( pn, x pn) = - (1/An) [ k'n+1 - ( kn+1 / kn) k'n ] (pn, xn) But if we expand xn we get xn = pn + pn-1 + ..., and comparing powers of xn gives 1 = kn so we conclude that = 1/kn so we have then ( pn, x pn) = - (1/An) [ k'n+1 - ( kn+1 / kn) k'n ] (pn, pn) /kn = - (1/An) [ k'n+1 - ( kn+1 / kn) k'n ] hn/kn But this is equal to - 0 hn/An from our first line above, so we conclude that -0/An = - (1/An) [ k'n+1 - ( kn+1 / kn) k'n ] (1/kn ) = - kn / kn+1 [ k'n+1 - ( kn+1 / kn) k'n ] (1/kn ) = - [ k'n+1 - ( kn+1 / kn) k'n ] (1/kn+1) = - [ ( k'n+1/kn+1) - ( k'n/kn) ] = - [ rn+1 - rn ] and therefore 0 = An [ rn+1 - rn ] = Bn as claimed! Therefore, our recursion relation is: pn+1 = pn ( 0 + Anx) + 1 pn-1 = pn(Bn + An x) - Cnpn-1 // agrees with (7) QED Comments: This recursion relation is a VERY general result. We have not assumed any particular interval (a,b) or weight w(x). We have the pn normalized to ||pn||2 = hn and in this normalization, we are required to know the leading two coefficients kn and kn' . We have not had to say how many polys we have in a set of polys. Only one fact is needed to show the recursion relation, and that one fact is the orthogonality of the functions pn . The factor x is in there because you need it so you can subtract out the highest order term in the proof. Only x can do this, not some other g(x). So ortho polys always have "x" in the recursion relation. This is not necessarily true for other special functions. Comment on p-1(x). If we make up this item and set it = 0, then we don't have to alter the recursion formula for application to n = 1. There is no such function in our set, or course, since we defined the set to have indices 0,1,2.... and on upwards. Theorem 3.4 (The Christoffel-Darboux formulas). The first claim is this: !Syntax Error, I(1/hi) pi(x)pi(y) = ( kn/kn+1) (1/hn) (1/(x-y)) [ pn+1(x)pn(y) - pn+1(y)pn(x) ] Once you have the above, you can write y = x+dx so that (x-y) = -dx in the denominator, and the numerator becomes [ pn+1(x)pn(x+dx) - pn+1(x+dx)pn(x) ] = dx [ pn+1(x)pn'(x) - p'n+1(x)pn(x) ], the leading terms canceled, and we get the second term first due to the -dx downstairs, so result is !Syntax Error, I(1/hi) [pi(x)] 2 = ( kn/kn+1) (1/hn) [ pn(x) p'n+1(x) - pn'(x) pn+1(x) ] Notice that the upper limit of the sum on the left matches the ortho poly labels on the right. So this is some kind of fancy "sum rule" that will no doubt get used a lot. Proof of the first claim above: Start with the recursion formula which is this: pn+1 = pn (Bn + An x) - Cnpn-1 Isolate the xpn term to get An x pn = pn+1 - Bnpn + Cnpn-1 Now write this explicitly for x and y: An x pn(x) = pn+1(x) - Bnpn(x) + Cnpn-1(x) An y pn(y) = pn+1(y) - Bnpn(y) + Cnpn-1(y) Multiply first by pn(y) and second by pn(x) An x pn(x) pn(y) = pn+1(x) pn(y) - Bnpn(x) pn(y) + Cnpn-1(x) pn(y) An y pn(x)pn(y) = pn+1(y)pn(x) - Bnpn(y)pn(x) + Cnpn-1(y)pn(x) Subtract An (x-y) pn(x) pn(y) = [ pn+1(x) pn(y) + Cnpn-1(x) pn(y)] - [ pn+1(y)pn(x) + Cnpn-1(y)pn(x)] where the two middle terms canceled. Regroup terms on the right to get An (x-y) pn(x) pn(y) = [ pn+1(x) pn(y) - pn+1(y)pn(x)] - Cn[ pn-1(x) pn(y) - pn-1(y)pn(x)] Reorder factors on the right to put n first, An (x-y) pn(x) pn(y) = [ pn+1(x) pn(y) - pn+1(y)pn(x)] - Cn[ pn(y) pn-1(x) - pn(x) pn-1(y)] We want 1/hn on the left, so multiply both sides by 1/(hnAn) (1/hn) (x-y) pn(x) pn(y) = [ pn+1(x) pn(y) - pn+1(y)pn(x)]/(hnAn) - Cn[ pn(y) pn-1(x) - pn(x) pn-1(y)]/(hnAn) Now compute: Cn/(hnAn) = 1/(hn-1An-1). Then we have (1/hn) (x-y) pn(x) pn(y) = [ pn+1(x) pn(y) - pn+1(y)pn(x)]/(hnAn) - [ pn(y) pn-1(x) - pn(x) pn-1(y)]/(hn-1An-1) But write this as (1/hn) (x-y) pn(x) pn(y) = Fn+1 - Fn where Fn+1 = [ pn+1(x) pn(y) - pn+1(y)pn(x)]/(hnAn) Notice that F0 = [ p0(x) p-1(y) - p0(y)p-1(x)]/(h0A0) = 0 since p-1= 0 as discussed above. Now change index to i to be consistent with theorem: (1/hi) (x-y) pi(x) pi(y) = Fi+1 - Fi Sum both sides i = 1 to n to get a simple telescoping result, !Syntax Error, I(1/hi) (x-y) pi(x) pi(y) = Fn+1 - F0 = Fn+1 Move the (x-y) to the right to get !Syntax Error, I(1/hi) pi(x) pi(y) = (1/(x-y)) [ pn+1(x) pn(y) - pn+1(y)pn(x)]/(hnAn) Finally, compute 1/(hnAn) = kn/ (hnkn+1) to get !Syntax Error, I(1/hi) pi(x) pi(y) = (kn/kn+1) (1/hn)(1/(x-y)) [ pn+1(x) pn(y) - pn+1(y)pn(x)] and this is our claimed result. No rocket science, just lots of fiddling. Hildebrand p 342 was my guide for this proof. QED. The rest of this section I don't think is useful to me right now, so we move on. 1.4. Mechanical Quadrature (meaning Gaussian Quadrature). In this section, we have statements of facts that I am now very familiar with from my Scheid reading. The "weights" which Scheid calls Ai are here called n where the first index is like i, and the second reminds you that you are fitting n points in your interval (a,b). All expressions from (1) through (6) on page 161 apply to an arbitrary set of x-values x (as they are called here). However, it is noted in mid-page that the exactness of the fit J to the integral I greatly improves if you select the x as the zeros of pn(x). The Gaussian quadrature weights n= A are here called the Christoffel numbers. The results (6) and (7) are only true when we select the x as the zeros of pn(), because then we can write pn(x) = constantn * n(x). By the way, we know that pn(x) has n zeros in the interval (a,b) and we thus know the zeros are all real, so pn(x) must factor completely over R. These two (6) and (7) results appear in Scheid page 129 where the integrands are really just Li and Li2 . It is interesting that both integrals are the same. (Of course in Scheid they contain instead of p). The Big Result of this section, however, is equation (8) for the weights, again only true when x are the zeros of pn(x). You very rarely see this formula stated in the general case, but here it is! Theorem 4.1: (weight formula) The Gaussian quadrature weights are given by n = - hn (kn+1/kn) / [ pn'(x) pn+1(x) ] = + hn-1(kn/kn-1) (1/[pn'(x)pn-1(x)]) Remember that hn = (pn, pn) = the normalization factor, and kn is the coefficient of xn in pn(x). The second form does not appear in Erdelyi, but I derive it below. Proof: You can see that the leading factors on the right are the same as appear in the C-D formula, so we expect that to be used in the proof! To start off, write the C-D formula setting y = x. !Syntax Error, I(1/hi) pi(x) pi(x) = (kn/kn+1) (1/hn)(1/(x-x)) [ pn+1(x) pn(x) - pn+1(x)pn(x)] = - (kn/kn+1) (1/hn)(1/(x-x))pn+1(x)pn(x) because the first term on the right vanishes since pn(x)= 0. We next apply w(x)dx to both sides. We can think of this as (1/k0) dx w(x) p0(x)dx. When we apply this to the LHS, all terms except i=0 vanish because (pi, p0) = 0. For i=0, the term is (p0, p0)/h0 = 1. So we get this result: 1 = - (kn/kn+1) (1/hn)pn+1(x) dx w(x) (1/(x-x))pn(x) so we have by this means evaluated the (certainly non-obvious) integral on the right to be: dx w(x) (1/(x-x))pn(x) = - [ hn/pn+1(x) ] (kn+1/kn) and this is true for all x which are zeros of pn(x). This is the result we now need because we have from (7) that n = dx w(x) L(x) = (1/pn'(x)) dx w(x) (1/(x-x))pn(x) and the integral here is just the one evaluated, so we have n = - (1/[pn'(x)pn+1(x)])hn (kn+1/kn) QED. Now, from our recursion formula on page 158 we have, when x = a zero, pn+1(x) = -Cn pn-1(x) since the middle term vanishes. We can then rewrite our formula as: n = + (1/[pn'(x)pn-1(x)])(hn /Cn) (kn+1/kn) but now compute (hn /Cn) (kn+1/kn) = hn-1(kn/kn-1), so we get this simpler result n = hn-1(kn/kn-1) (1/[pn'(x)pn-1(x)]) Comment: for the Classical Poly cases, we get formulas for pn'(x) which let us replace pn'(x) with pn-1(x) times other factors, and then we get a formula for the above where we don't need any derivatives of pn(x) in the weight formula. However, I don't think there is a simple general formula, but we shall see later! We now skip the rest of this section. 10.5 Continued Fractions I don't even understand the very first equation (1)'s notation with the vertical bar, nor do I care about continued fractions right now, so skip this whole section. 10.6 The "Classical" Orthogonal Polynomials. Start by assuming the R formula (1) on page 164 and see what we get. Note that X(x) = degree k. We find that the only possibilities are k = 0,1,2 as he now shows. Case 0: [ k = 0 => Hermite polynomials. ] Suppose X = X0 a constant, meaning k = 0. Assume that p1(x) = k1x + k1'. Then (2) tells us that K1 ( k1x + k1') = X0 w'(x)/w(x) => w'(x)/w(x) = (K1/X0) ( k1x + k1') Now define a change of variable y = x + so that dy = dx. Then define W(y) = w(x) => W'(y) = dW/dy = dw/dy = dw/dx * dx/dy = w'(x) (1/). Then we have W'(y)/W(y) = w'(x)/w(x) * (1/) = (K1/X0) ( k1x + k1') = -2y We want to select and so that this last expression is in fact -2y. Thus, we have y = -(K1/2X0) ( k1x + k1') from which we can read off that = -(k1K1/2X0) and = -(K1/2X0) k1'. The first equation says that 2 = -(k1K1/2X0) => = i sqrt(k1K1/2X0) and then = whatever it is. The point is that we have change variables from x to y, and we end up with W'(y)/W(y) = -2y and this then tells us dW/W = -2ydy and then lnW(y) = -y2 and then W(y) = exp(-y2). Our change of variables is linear but a bit ugly, but it does seem to lead to the Hermite weight. We then write everything in terms of the y variable, including the Rodriguez formula. You can see that for the Hermite's, the R formula is exact as Erdelyi gives it on page 164, but we write it in y, and we are done if we set X = 1 and w(y) = exp(-y2) and call the functions pn(y) = Hn(y). This shows up on page 785 of A&S, where they use n for Erdelyi's Kn. Of course A&S nail down the value of n. Case 1: [ k = 1 => Laguerre polynomials. ] Here we have to do some kind of more complex variable change (I don't think it is just linear as he says), to get to a form W'(y)/W(y) = -1 + /y . Here from (2) we get (3) and we can then write w'(x)/w(x) = (Ax+B)/(x+) = c [ - 1+ d/(x+e)] with apropo constants. As our first change of variables, let z = fx so dz = f dx and then W'(z)/W(z) = [w'(x)/w(x)] (1/f). Then we get W'(z)/W(z) = (1/f) c [ - 1+ d/((z/f)+e)] . So choose f = c and get W'(z)/W(z) = -1 + g/(z+h). Now do a second change of variables from z to y of the form y = z+h and get W'(y)/W(y) = -1 + g/y, and this is the desired form, but he writes it W'(y)/W(y) = -1 + /y. If you solve this equation dW/W = (-1 + /y)dy by integrating to lnW = -y + lny then expo to get W = exp(-y)exp(lny) = exp(-y) y and we have then arrived at the Laguerre poly world. Without loss of generality, you can pick X(y) = y. As before, Rodriguez comes out exactly as you think, see A&S page 785. Case 2: [ k = 2 => Hypergeometric (Jacobi) polynomials. ] I don't follow things here but accept the results. Looking at the table on page 164, we can see that the second line is just a special case of this case, and the first is a special case of the second line. Authors show on page 165 (but I don't follow) that if k > 2, you find that p2(x) is not of degree 2, so things are inconsistent, so the only allowed (and distinct) Rodriguez choices are X = 1, x and (1-x)2. 10.7 General Properties of the classical ortho polys. This section opens with this claim Theorem 7.1: The quantity Dk[w(x)Xn] vanishes at a and b for the three lower cases of the table on page 164 provided > -1 and > -1 and > -1 in the Laguerre case, as long as k = 0,1,2...up to (n-1). Proof: Let's first look at the Jacobi case. There we have Dk [ (1-x)+n (1+x)+n ] as our quantity of interest. There are many terms in this expansion, but each time D hits something like (1-x)+n, it lowers the power by one. So the worse case is a term like Dn-1 (1-x)+n = (+n)(+n-1).... (1-x)+1 . This vanishes at x = 1 as long as > -1, as they claim. Same idea on the side for x = -1, so I accept the Theorem as true for this case. What about the Hermite case? There we have Dk [ e-x*x ] which is just e-x*x x poly, and this vanishes at the two endpoints which are infinite. What about the Laguerre case? Here we have Dk [x+n e-x ]. As in our first case, the worse case for x=0 occurs when we grind down the first exponent, but we will have (x)+1 again, and we will be OK then for > -1. The infinite side will take care of itself since e-x will be in every term. QED. Theorem 7.2: The set of polynomials defined by the Rodriguez formula forms an orthonormal system. Proof: Using the above theorem, author shows that (f,pn) = 0 for any poly f of degree < n, where pn is the Rodriguez poly. But since all pk for k < n are candidates for this f, we know that (pk,pn) = 0 k < n, and this means that we have an orthonormal system! Comment: This is a very nice result. First we saw that we could only have k = 0,1,2 in the Rodriguez and we took standard forms for each case to get our three cases. Now we see that in each case, we get an orthogonal system on (a,b) with the w(x) and X shown in the table on page 164. This is of course only a small subset of the universe of possible orthogonal systems. One wonders if there is a more general R formula that can generate ALL orthogonal systems? Lemma 7.3. Derive the 4th equation from the bottom of page 166 using Leibniz' Formula. Proof: That formula is this: Dn(fg) = i=0n (n,i) Dif * Dn-ig where (n,i) is the binomial coefficient. Our application is going to have one higher power, so write Dn+1(fg) = i=0n+1 (n+1,i) Dif * Dn+1-ig Now we take f = X and g = D(wXn). Notice that Dn-1g = Dn(wXn) = Knwpn from R formula. Now in the sum, since K is at most degree 2, only the first three terms will survive. They are: (n+1,0) X Dn+1g + (n+1,1) X' Dng + (n+1,2) X" Dn-1g = [ 1*X*D2 + (n+1)*X' D + n(n+1)/2 * X" ] (Dn-1g) = [ X D2 + (n+1) X' D + n(n+1)/2 X" ] (Knwpn ) and this is the claimed result, QED. Lemma 7.4. Derive the 3rd equation from the bottom of page 166. Proof: This equation has two lines. The first line is proven if we can show that XD(wXn) = (wXn)[ K1p1 + (n-1)X' ] First, do the left side to get XD(wXn) = X [ wnXn-1 + w' Xn ] = wnXnX' + w' Xn+1 = wXn [ nX' + (w'/w)X ] Now use 10.6.3 which says (w'/w) = (K1p1 - X')/X to get XD(wXn) = wXn [ nX' + (K1p1 - X') ] = (wXn) [ K1p1 + (n-1)X' ] QED. Now go back to our Leibniz formula with n+1 Dn+1(fg) = i=0n+1 (n+1,i) Dif * Dn+1-ig and this time take f = [ K1p1 + (n-1)X' ] and g = (wXn). Since f is of degree 1 at most, only the first two terms survive and we get (n+1,0) f Dn+1g + (n+1,1) Df Dng = [ 1*f* D + (n+1)*Df ] Dng = { f D + (n+1) (Df) } (Knwpn ) = { [ K1p1 + (n-1)X' ] D + (n+1) [ K1p1' + (n-1)X" ] } (Knwpn ) and this is the second line QED. **************************************************************************** Theorem 7.3. Show that the polynomial pn(x) satisfies the ODE shown as (1) on bottom of page 166 where n is given by (2) in which k1 is the leading coefficient of p1(x), that is, p1(x) = k1x + k1'. Proof: It all works, but a HUGE amount of algebra is required! I don't know why it is so messy. Perhaps another day I can let Maple do this. STEP 1. Compare results of previous equations as they say to get: [ X D2 + (n+1) X' D + n(n+1)/2 X" ] (Knwpn ) = { [ K1p1 + (n-1)X' ] D + (n+1) [ K1p1' + (n-1)X" ] } (Knwpn ) Cancel the Kn and gather like terms to the left side: [ X D2 +{ (n+1) X' - [ K1p1 + (n-1)X' ] } D + { n(n+1)/2 X" - (n+1) [ K1p1' + (n-1)X" ] } ] (wpn) = 0 Simplify some factors { (n+1) X' - [ K1p1 + (n-1)X' ] } = { (n+1) X' - K1p1 - (n-1)X' ] } = { 2X' -K1p1} { n(n+1)/2 X" - (n+1) [ K1p1' + (n-1)X" ] } = (n+1){ n/2 X" - [ K1p1' + (n-1)X" ] } = (n+1){ n/2 X" - K1p1' - (n-1)X" } = (n+1)/2{ n X" - 2K1p1' - 2(n-1)X" } = (n+1)/2{ (-n+2) X" - 2K1p1' } = - (n+1){ (n-2) X" + 2K1p1' }/2 Therefore we get [ X D2 +{2X' - K1p1 } D - (n+1) { (n-2)/2 X" + K1p1' } ] (wpn) = 0 STEP 2: This is still a pretty messy result! Now do this: D(wpn) = w'pn + wDpn D2(wpn) = w"pn + w'Dpn + wD2pn + w' Dpn = w"pn + 2w'Dpn + wD2pn Jam this in and get X{ w"pn + 2w'Dpn + wD2pn} + {2X' - K1p1 } {w'pn + wDpn } - (n+1){ (n-2)/2 X" + K1p1' }wpn = 0 Now group terms: XwD2pn + [ 2Xw' + {2X' - K1p1 }w] Dpn + [ Xw" + {2X' - K1p1 }w' - (n+1){ (n-2)/2 X" + K1p1' }w]pn = 0 Now divide by w to make the first term match: XD2pn + [ 2Xw'/w + {2X' - K1p1 }] Dpn + [ Xw"/w + {2X' - K1p1 }w'/w - (n+1){ (n-2)/2 X" + K1p1' }]pn = 0 STEP 3: The first term now has the right form. To verify the second term, we need to show that: [ 2Xw'/w + {2X' - K1p1 }] = K1p1 We can use (3) to replace 2Xw'/w = 2(K1p1- X') and then we have LHS = [ 2(K1p1- X') + {2X' - K1p1 }] = K1p1 hurray! STEP 4: Now it only remains to show that the third term is right. It's factor is this [ Xw"/w + {2X' - K1p1 }w'/w - (n+1){ (n-2)/2 X" + K1p1' }] = LHS Let's differentiate (2) page 164 to get K1p1' = X" + [ -X(w'/w)2 + Xw"/w + X'w'/w ] Put K1p1' in for the last factor and we get LHS = [ Xw"/w + {2X' - K1p1 }w'/w - (n+1){ (n-2)/2 X" + (X" + [ -X(w'/w)2 + Xw"/w + X'w'/w ]) }] = [ Xw"/w + {2X' - K1p1 }w'/w - (n+1){ (n-2)/2 X" + X" -X(w'/w)2 + Xw"/w + X'w'/w }] Now replace K1p1 = X' + Xw'/w in the first term to get {2X' - K1p1 } = X' - Xw'/w, then we have LHS = [ Xw"/w + {X' - Xw'/w }w'/w - (n+1){ (n-2)/2 X" + X" -X(w'/w)2 + Xw"/w + X'w'/w }] If we now group factors like (w"/w) and so on [ I did this on scratch paper] we simplify to: LHS = -n { X(w"/w) + X'(w'/w) - X(w'/w)2 + (n+1)X"/2 } Looking now at the predicted answer: we have to show that -n { X(w"/w) + X'(w'/w) - X(w'/w)2 + (n+1)X"/2 } = -n { k1K1 + (n-1)/2X" } which is to say, we must show that X(w"/w) + X'(w'/w) - X(w'/w)2 + (n+1)X"/2 = k1K1 + (n-1)/2X" or we must show that X(w"/w) + X'(w'/w) - X(w'/w)2 + X" = k1K1 where I presume k1 is the leading coefficient in p1 (the only way we have ever used this symbol). Now suppose we write p1(x) = k1x + k1', then the R formula says k1x + k1' = (1/K1w) D[wX] Differentiate both sides then to get: K1k1 = D{ w-1 D[wX] } = D{ X' + (w'/w)X } = X" + X'(w'/w) + (w"/w)X - (w'/w)2X But this is exactly what we want, so we are finally done!!! QED. Comments: In the first part of this chapter, we talk about orthogonal polynomials in general. There is no differential equation at all. Only when we start talking about the restricted class of ortho polys given by the Rodriguez formula do we end up with an ODE. Deriving the form of the ODE is non-trivial, but we have done it in full detail. ***************************************************************************** Lemma 7.5. The expression Xpn' - (n/2)X" x pn = f where poly f is of degree n. Proof: For X = zeroth order = constant (k=0), we have pn' = degree n-1 and X" = 0, so deg(f) = n-1. For X = first order, Xpn' has a leading power of xn and again X" = 0, so deg(f) = n. For X = second order, things are a little more complicated. Assume that X = ax2 + o(x), so that X" = 2a. Our expression is then [ ax2 + o(x)] pn' - na x pn. If we write pn = knxn + lower, then we have the leading power of the first term being ax2 * nkn xn-1 = anknxn+1. But the second term has -na knxn+1 so these two terms cancel and the leading term must then be order pn. Notice that our Lemma is true for any quadratic form of X, such as X = (1-x)2 or 1 - x2. This shows that the lemma is true in all three cases k = 0,1 and 2 (where k is the power of X in Rodriguez ). Notice *************************************************************************** Lemma 7.6. Show that I = ( Xpn', pn-k ) = 0 where k = 2,3,4...n. Proof: Write this out as dx w pn-k X (Dpn) and do parts once, parts vanish hopefully with same conditions as in the parts done at the top of page 66. This gives I = dx w pn-k X (Dpn) = - dx D(wX pn-k )pn Now write D(wX pn-k ) = [ D(wX) pn-k + (wX) pn-k' ] so we how have I = - dx [ D(wX) pn-k + (wX) pn-k' ] pn The second term has the form (Xpn-k' , pn) . The left factor has degree n-k -1+2 = n -k+1 at most. Since we are dealing here only with k 2, left factor has degree n-1 at most. This scalar product thus vanishes since we can write this factor as a sum of pn of degree n-1 max. We are then left with I = - dx [ D(wX) pn-k ] pn We can now write D(wX) = K1 w p1 from Rod, so we have I = - K1 dx w p1 pn-k pn = -K1 (p1 pn-k, pn) Now the product p1 pn-k is at most of degree n-k+1 which is the same situation as before with the "left factor" and by the same argument we get 0 for this scalar product. Therefore I = 0, QED. ******************************************************************************** Theorem 7.4. Derive the differentiation formula stated on page 167 as (4). Proof: Note: This proof assumes only that X is a quadratic with leading term X = ax2. STEP 1: ( Show i = 0) Start as suggested and we then have that Xpn' - (n/2)X" x pn = f(x) = npn + npn-1 + sum(i=2,n) i pn-i We get to this point because we know f(x) is degree n or less, and we know we can expand any such poly in terms of the full set of pn as shown. The coefficients are given strange names, fine. So right off the bat we have the following version of our formula: Xpn' = [n + (n/2)X" x ] pn + npn-1 + sum(i=2,n) i pn-i So our task is to evaluate n and n and to show that the i = 0. Let's first show that 2 = 0. Close both sides of second last equation with ( ,pn-2) to get ( Xpn' - (n/2)X" x pn, pn-2) = 2 (pn-2, pn-2) Write the LHS as ( Xpn', pn-2 ) - (n/2)X" (x pn, pn-2) but we know that (x pn, pn-2) = (pn, x pn-2) = 0 because xpn-2 is of degree n-1 and can be expanded in terms of pn for n-1 and below, all of which are orthogonal to pn. So we are then left with ( Xpn', pn-2 ) = 2 (pn-2, pn-2) and we want to show that the LHS = 0. But this is what our Lemma 7.6 above shows, so we have 2 = 0. The same discussion causes all the k = 0 for k 2. The thing we get is in that case ( Xpn', pn-k ) = k (pn-k, pn-k) k 2 and our Lemma applies to all these cases at once, so k = 0. We then have this result so far: Xpn' = [n + (n/2)X" x ] pn + npn-1 which shows the essential fact that the derivative mixes pn with one adjacent poly. We now have to isolate the two parameters and evaluate them! STEP 2: (Solve for n) If we close with pn-1 we get (Xpn', pn-1) = (n/2)X" (xpn, pn-1) + n(pn-1, pn-1) and we have something we might solve for n. We can evaluate this: (xpn, pn-1) = (pn, xpn-1) = (pn, x kn-1xn-1 + o(xn-1) ) where we can drop lower terms in the expansion = kn-1 (pn, xn) = kn-1 (pn, (1/kn)pn + o(xn-1)) and again drop extra terms = (kn-1/kn) (pn, pn) OK So we then have I = (Xpn', pn-1) = (n/2)X" (kn-1/kn) (pn, pn) + n(pn-1, pn-1) and we have now to deal with the LHS. Following the start of our Lemma above, we can do parts to get - I = dx [ D(wX) pn-1 + (wX) pn-1' ] pn = I1 + I2 Use the Rod formula D(wX) = K1 w p1 in the first term to get I1 = K1 (p1, pn-1pn) = K1 (pn, pn-1p1) = K1 (pn, pn-1 [k1x + o(x0)] ) = K1k1(pn, xpn-1) = K1k1 (xpn, pn-1) but we just evaluated this bracket above so we can go on to say I1 = K1k1 (kn-1/kn) (pn, pn) Now consider I2 = (Xpn-1' ,pn) = (X {(n-1)kn-1 xn-2 + o(xn-3)}, pn) = (n-1) kn-1 (X xn-2, pn) The rightmost bracket is nonzero only if X is the k=2 case. In that case we have (X xn-2, pn) = ( [ ax2 + o(x)] xn-2, pn) = a(xn, pn) = a( pn/kn + o(xn-1) , pn) = (a/kn)(pn, pn) Therefore we have shown that I2 = a(n-1)(kn-1/kn)(pn, pn) for this one case, and it is zero for the lower cases. Thus, since X" = 0,0,2a for our 3 cases, we can add a factor (X"/2a) and say I2 = (n-1)(kn-1/kn)(pn, pn)(X"/2) and we have then covered all three cases. So we have now learned that I = (Xpn', pn-1) = -I1 - I2 = - K1k1 (kn-1/kn) (pn, pn) - (n-1)(kn-1/kn)(pn, pn) (+ X"/2) = - (kn-1/kn) [ K1k1 +(n-1)(X"/2) ] (pn, pn) but from above we also have that I = (n/2)X" (kn-1/kn) (pn, pn) + n(pn-1, pn-1). Therefore: - (kn-1/kn) [ K1k1 +(n-1)(X"/2) ] (pn, pn) = (n/2)X" (kn-1/kn) (pn, pn) + n(pn-1, pn-1) which we can write as - (kn-1/kn) [ K1k1 +(n-1)(X"/2) ] hn = n (X"/2) (kn-1/kn) hn + nhn-1 and we can group the hn terms to get nhn-1 = -[K1k1 + (2n-1)(X"/2) ] hn (kn-1/kn) and finally n = - (hn/hn-1) (kn-1/kn) [K1k1 + (2n-1)(X"/2) ] // my result We can now compare this with Erdelyi's formula for n n = -(Cn/An)[ k1K1 + (n - 1/2)X" ] but I know that (Cn/An) = hn/(An-1hn-1) = (hn/hn-1)(kn-1/kn) , so his answer then becomes n = -(hn/hn-1)(kn-1/kn)[ k1K1 + (2n-1)(X"/2)] / his which agrees with my result, so happy now about n. STEP 3: (Find n for Jacobi polys) As an application, let's look at his Jacobi differential formula on p 170 top (2n++)X pn' = first term + 2(n+)(n+)pn-1 and our general formula said, Xpn' = [n + (n/2)X" x ] pn + npn-1 Let's divide his through to get X pn' = first term' + 2(n+)(n+)/(2n++)pn-1 so this would suggest that he uses n = 2(n+)(n+)/(2n++). This result is explicitly stated in his constants section, good! Let's now look at our now-derived formula for n n = - (hn/hn-1)(kn-1/kn)[ k1K1 + (2n-1)(X"/2) ] The constants on page 169 are just a huge mess! I just computed some things, (hn/hn-1) = (1/n) (kn-1/kn) = 2n The product of these two things is then, (hn/hn-1)(kn-1/kn) = 2 And we have k1 = (1/2) (2++) and K1 = -2 so that k1K1 = -(2++). Then we get n = - 2 [ -(2++) + (2n-1)(X"/2) ] If we now assume (for the first time ever) that X = 1 - x2 so X" = -2, this becomes n = - 2 [ -(2++) - (2n-1)] = 2 [ (2++) + (2n-1)] = 2 [ (2n++ +1)] = 2 = the right answer. STEP 4: (Compute n). First here are three little Lemmas we will be using below: _____________________________________________________________________ Lemma S1 (xn, pn) = (1/kn)hn Proof: (xn, pn) = ( pn/kn + o(xn-1) , pn) = (1/kn)(pn, pn) QED. _____________________________________________________________________ Lemma S2 (pn, xn+1) = -rn+1 (pn, xn) = -rn+1 (1/kn)hn Proof: Write pn+1/kn+1 = xn+1 + (k'n+1/kn+1) xn + o(xn-1) and solve for xn+1 (pn, xn+1) = (pn, pn+1/kn+1) - (k'n+1/kn+1)(pn, xn) - terms like (pn, xn-1) and lower = - (k'n+1/kn+1)(pn, xn) = - rn+1(pn, xn) which gives the first result. Then use lemma (1)'s result to get the final result. ______________________________________________________________________ Lemma S3 (pn, xpn) = kn[ rn - rn+1](pn, xn) = ( rn - rn+1) hn Proof: Start with (pn, xpn) = (pn, x [ knxn + kn'xn-1 + o(xn-2) ] ) = kn (pn, xn+1) + kn' (pn, xn) Use lemma 2 to get = kn {-rn+1 (pn, xn)} + kn' (pn, xn) = [ kn' - kn rn+1 ] (pn, xn) = kn[ rn - rn+1](pn, xn) = kn[ rn - rn+1](1/kn)hn = [ rn - rn+1] hn where we used lemma 1 for (pn, xn). QED. _______________________________________________________________________ (Compute n) From the end of STEP 1 we quote this result: Xpn' = [n + (n/2)X" x ] pn + npn-1 and we have already found n. Let's now close with pn to get: (pn, Xpn' ) = n (pn, pn) + n(X"/2) (pn, xpn) = nhn + n(X"/2) (pn, xpn) Let's start with (pn, Xpn' ). If we were to do a parts integration at this point, we pick up constants like k1 and K1 just as we did in computing n . This method is reviewed in Appendix A. We end up with a complex result for n and we cannot even tell if it is right or wrong! So instead, to get the simple form result which appears in Erdelyi, we decide here to NOT do parts integration at this step, and instead we do this: (pn, Xpn' ) = (pn, [ X(0) + x X'(0) + (x2/2) X"] [ nknxn-1 + (n-1)k'n xn-2 + o(xn-3) ] ) = X(0) * 0 + X'(0) * nkn (pn, xn) + (X"/2)*[ nkn (pn, xn+1) + (n-1)k'n (pn, xn) ] = (X"/2)* nkn (pn, xn+1) + [ X'(0) * nkn + (n-1)k'n (X"/2)] (pn, xn) Now let's use Lemma S2 (pn, xn+1) = - (k'n+1/kn+1)(pn, xn) = -rn+1 (pn, xn) so we then have (pn, Xpn' ) = (X"/2)* nkn { -rn+1 (pn, xn) } + [ X'(0) * nkn + (n-1)k'n (X"/2)] (pn, xn) = [ X'(0) * nkn + (n-1)k'n (X"/2) - rn+1 (X"/2)* nkn ](pn, xn) = kn [ X'(0) * n + (n-1)rn (X"/2) - rn+1 (X"/2)* n ](pn, xn) = kn [ X'(0) * n + (n-1)rn (X"/2) - rn+1 (X"/2)* n ](pn, xn) From Lemma S1 we have (xn, pn) = ( pn/kn + o(xn-1) , pn) = (1/kn)(pn, pn) = (1/kn)hn so our non-parts version of (pn, Xpn' ) is (pn, Xpn' ) = [ X'(0) * n + (n-1)rn (X"/2) - rn+1 (X"/2)* n ] hn // non-parts Meanwhile, the RHS is this (pn, Xpn' ) = nhn + n(X"/2) (pn, xpn) and from Lemma S3 we have (pn, xpn) = ( rn - rn+1) hn so we have for the RHS (pn, Xpn' ) = nhn + n(X"/2) ( rn - rn+1) hn Now here then is our total equation (in the non-parts path) nhn + n(X"/2) ( rn - rn+1) hn = [ X'(0) * n + (n-1)rn (X"/2) - rn+1 (X"/2)* n ] hn and we can cancel the hn at once to get n + n(X"/2) ( rn - rn+1) = [ X'(0) * n + (n-1)rn (X"/2) - rn+1 (X"/2)* n ] so we then have n = - n(X"/2) ( rn - rn+1) + [ X'(0) * n + (n-1)rn (X"/2) - rn+1 (X"/2)* n ] = rn { - n(X"/2) + (n-1) (X"/2) } + rn+1 { n(X"/2) - (X"/2)* n } + nX'(0) = rn(X"/2) { -n +n - 1} + 0 + nX'(0) = nX'(0) - rn(X"/2) which is the right short-form answer, as shown in 10.7.5. Theorem 7.5. Derive the "other" formula for X pn' . Proof: Our starting point is the result on page 167 in equations (4) and (5) which we have derived in Theorem 7.4 previous. Now we use the recurrence formula on page 158 pn+1 = (Anx + Bn) pn - Cn pn-1 and we solve for pn-1 pn-1 = (x An/Cn + Bn/Cn) pn - pn+1/Cn We then jam this into (4) on page 167 to get Xpn' = (n + n(X"/2) x) pn + n { (x An/Cn + Bn/Cn) pn - pn+1/Cn } = ( [ n+ n(Bn/Cn) ] + [ n(X"/2) + n(An/Cn) ] x) pn - (n/Cn) pn+1 = ( a + bx) pn + d pn+1 This result does appear more complicated, but not really sure. We prefer the other formula because we are happier having pn-1 instead of pn+1 in any formula! Theorem 7.6. Derive the two formulas for the Gaussian Weights shown in page 167 equation (7). Proof: These are (luckily) very easy to show. I am going to leave out the arguments of the pn and X here but they are all x , the zeros of pn(x). We start off with n = - hnAn/[ pn' pn+1] as we already derived on page 162 (8) Our recursion relation says that pn+1 = -Cnpn-1. So we can then have this alternative form n = + hn(An/Cn) / [ pn' pn-1] = hn-1 An-1 / [ pn' pn-1] Next, (4) says that X pn' = n pn-1. This then gives us n = hn-1 An-1 ( X/n) / [ pn-1] 2 which is perhaps the best form, all eval at x = x The other alternative is to instead replace pn-1 to get n = hn-1 An-1 ( n/X) / [ pn'] 2 and this concludes proof of formulas (7) A few comments on the Differential Equation notations. As presented in Erdelyi, we have in general X pn" + K1p1 pn' + npn = 0 and for the zeros of pn we have Xpn' = n pn-1 where nn = n(Cn/An) // x = xv only Now how does this relate to A&S? From page 781, it certainly seems that they use g2 = X since we only see the three values 1, x and 1-x2 in this entire table. In the Rodriquez, we see that they have written g2 = g and en = Kn . Then they go on to define these things kn and kn' and hn same as in Erdelyi bn = kn+1/kn = An an = bn (rn+1 - rn) = An (rn+1 - rn) = Bn cn = Cn g2 = X en = Kn g1 = K1p1 g0 = n for the pure orthog poly cases in Table page 781 = n more or less in Table page 783 For sure, on page 783 we have g0/g2 = n/X. Now let's write the Gaussian Weights in A&S notation. Start from above, n = hn-1 An-1 ( X/n) / [ pn-1] 2 = hn-1 (kn/kn-1) ( X/n) / [ pn-1] 2 = hn-1 (kn/kn-1) ( g2/g0) / [ pn-1] 2 and this agrees with my separate derivation of this result in my 1978 or older notes. ******************************************************************************* Appendix A: In developing an expression for the n of equation 10.7.5, it seems possible to start the analysis by doing a parts integration, in a manner similar to how we evaluated n. That method leads to a very complex answer for n which is either wrong (because parts is not justified or I goofed), or is right but in a way that is hard to see. In this appendix, I just want to record that calculation since I did the whole thing. Then there are comments at the end. (Compute n). From the end of STEP 1 we quote this result: Xpn' = [n + (n/2)X" x ] pn + npn-1 and we have already found n. Let's now close with pn to get: (pn, Xpn' ) = n (pn, pn) + n(X"/2) (pn, xpn) = nhn + n(X"/2) (pn, xpn) Let's start with (pn, Xpn' ). Do the usual parts integration to find that - (pn, Xpn' ) = + dx pn D[ wXpn] = + dx pn { (wX)pn' + pnD(wX) } The first term is the same object that is on the LHS, and for the second D(wX) = K1wp1 as we used earlier from the Rod formula, so we have - (pn, Xpn' ) = + (pn, Xpn' ) + K1 dx w pn pn p1 = (pn, Xpn' ) + K1 (pn, pn p1) which tells us that -2(pn, Xpn' ) = K1 (pn, pn p1) so that (pn, Xpn' ) = -(K1/2)(pn, pn p1) Now examine (pn, pn p1) = (pn, pn [ k1x + k1']) = k1(pn, xpn) + k1' hn. Thus we have found: (pn, Xpn' ) = -(K1/2) [ k1(pn, xpn) + k1' hn ] = -(k1K1/2) (pn, xpn) - (k1'K1/2) hn Inserting this result as the LHS of our starting closed equation, (pn, Xpn' ) = nhn + n(X"/2) (pn, xpn) we get -(k1K1/2) (pn, xpn) - (k1'K1/2) hn = nhn + n(X"/2) (pn, xpn) which we can rewrite as - nhn = [ k1K1/2 + n(X"/2 ) ] (pn, xpn) + (k1'K1/2) hn Next, we evaluate (pn, xpn): (pn, xpn) = (pn, x [ knxn + kn'xn-1 + o(xn-2) ] ) = kn (pn, xn+1) + kn' (pn, xn) Let's start with the left term on the right using pn+1/kn+1 = xn+1 + (k'n+1/kn+1) xn + o(xn-1) and solve for xn+1 (pn, xn+1) = (pn, pn+1/kn+1) - (k'n+1/kn+1)(pn, xn) - terms like (pn, xn-1) and lower = - (k'n+1/kn+1)(pn, xn) We can now combine this with the right term on the right to get (pn, xpn) = - kn (k'n+1/kn+1)(pn, xn) + kn' (pn, xn) = kn [ (kn'/kn) - (kn+1'/kn+1) ] (pn, xn) Then repeating a piece of work done way back, we have for the rightmost factor: (xn, pn) = ( pn/kn + o(xn-1) , pn) = (1/kn)(pn, pn) = (1/kn)hn and this tells us that (pn, xpn) = [ (kn'/kn) - (kn+1'/kn+1) ] hn = ( rn - rn+1) hn and we can now complete the solution using from above - nhn = [ k1K1/2 + n(X"/2 ) ] (pn, xpn) + (k1'K1/2) hn = [ k1K1/2 + n(X"/2 ) ]( rn - rn+1) hn + (k1'K1/2) hn to find that n = - [ k1K1/2 + n(X"/2 ) ] [ rn - rn+1 ] + (k1'K1/2) which, unfortunately, is not the correct answer! The book's answer is this: n= nX'(0) - (X"/2) rn which is very much simpler. There are three possibilities for this disagreement (1) the parts integration stage is unjustified. That would be true if the following did not vanish: pn2 wX | ab 0 ?? How would we check this? Consider Jacobi at x = 1, for example. w(x) ~ (1-x) and X ~ 2(1-x) and the pn is finite, so it would seem that the parts vanish at a for > -1. At the b = -1 endpoint we have w(x) ~ (1+x) and X ~ 2(1+x) so they should vanish here for > -1. In the Hermite case we have w being a huge expo at both end points, so parts = 0 in that case. In Laguerre, our only concern would be at a = 0 and there we have w ~ x and X ~ x so again we are fine if > -1. So my conclusion is that the parts integration I did was probably OK! (2) I made one or more algebraic errors. I did check everything once. Regardless of errors, my answer is never going to have X'(0) in it, and it IS going to have things like k1 and K1 in it. (3) The results really are the same. I suspect this is where the truth lies. Once again our two forms are: n = - [ k1K1/2 + n(X"/2 ) ] [ rn - rn+1 ] + (k1'K1/2) // me n = nX'(0) - (X"/2) rn // them Regardless of any errors I may have made, my parts method brings in things like k1 and K1 which don't appear in the second form of the result. Similarly, my method does not bring in X'(0). We know that k1 and k1' and X'(0) and K1 are connected by the equation p1K1w = D[wX] , but the connection is non-trivial. So I don't think you can just state flatly that the above two forms do not agree. Who knows what strange relations my be between all these constants. Appendix B: What would n have been had I used the non-parts approach there? I = (Xpn', pn-1) = (pn-1, Xpn'), so (pn-1, Xpn' ) = (pn-1, [ X(0) + x X'(0) + (x2/2) X"] [ nknxn-1 + (n-1)k'n xn-2 + (n-2) kn"xn-3 + o(xn-4) ] ) Here you see the problem. We pick up a term of this form: (X"/2) (n-2) kn" (pn-1, xn-1) and then we are stuck forever with this kn" parameter which is not one of Erdelyi's chosen parameters. Doing the parts integration gets is mixed up with k1 and K1 but at least these are an "the list" of OK parameters to use.