Una introducción a los polinomios ortogonales reales y a la cuadratura numérica de tipo gaussiano

Traballo Fin de Grao Una introducción a los polinomios ortogonales reales y a la cuadratura numérica de tipo gaussiano Eva López Troitiño Julio, 2022 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA GRAO DE MATEMÁTICAS Traballo Fin de Grao Una introducción a los polinomios ortogonales reales y a la cuadratura numérica de tipo gaussiano Eva López Troitiño Julio, 2022 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA iii iv Trabajo propuesto Área de Coñecemento: Matemática Aplicada Título: Una introducción a los polinomios ortogonales reales y a la cuadratura numérica de tipo gaussiano Breve descrición do contido El trabajo abordaría los siguientes puntos: i . Generalidades sobre polinomios ortogonales reales: concepto de función peso y existencia de un sistema de polinomios ortogonales respecto a una función peso. Propiedad de recurrencia. Propiedades de las raíces de los polinomios ortogonales. ii . Los polinomios ortogonales de Jacobi y sus propiedades elementales. Casos particulares: los polinomios de Legendre y los polinomios de Tchebyche. iii . Concepto de fórmula de cuadratura y de orden de una fórmula de cuadratura. iv . Fórmulas de cuadratura de Gauss, Gauss-Radau y GaussLobatto. Casos particulares: pesos de Jacobi, de Legendre y de Tchebyche. Recomendacións Outras observacións Índice Resumen vii Introducción ix 1. Generalidades sobre polinomios ortogonales reales 1 2. Los polinomios ortogonales de Jacobi y sus propiedades elementales 9 2.1. PolinomiosdeLegendre ................................ 9 2.2. Polinomios de Tchebyche . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.3. PolinomiosdeJacobi.................................. 21 3. Fórmulas de cuadratura de tipo gaussiano 25 3.1. Nociones básicas sobre fórmulas de cuadratura . . . . . . . . . . . . . . . . . . . 25 3.2. Fórmulas de cuadratura de Gauss, Gauss-Radau y Gauss-Lobatto . . . . . . . . . 28 3.2.1. FórmulasdeGauss............................... 28 3.2.2. Fórmulas de Gauss-Radau . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 3.2.3. Fórmulas de Gauss-Lobatto . . . . . . . . . . . . . . . . . . . . . . . . . . 34 3.3. Ejemplos de fórmulas de cuadratura de tipo gaussiano . . . . . . . . . . . . . . . 39 I. Función Gamma 47 Bibliografía 49 v Resumen El objetivo de este trabajo es presentar una introducción a las diferentes propiedades de los polinomios ortogonales reales y a la cuadratura numérica gaussiana. Empezaremos describiendo las propiedades generales de los polinomios ortogonales reales con respecto a una función peso, como la relación de recurrencia y las propiedades de sus raíces. Nos centraremos en estudiar las propiedades de los polinomios ortogonales de Jacobi y algunos de sus casos particulares (polinomios de Legendre y de Tchebyche). Después, presentaremos nociones básicas sobre fórmulas de cuadratura y nalmente nos centraremos en las fórmulas de cuadratura de tipo gaussiano. Para comprender mejor dichas fórmulas presentaremos algúnos ejemplos de estas cuadraturas para el peso de Legendre y el de Tchebyche. Abstract The objective of this work is to present an introduction to the dierent properties of real orthogonal polynomials and Gaussian numerical quadrature. We will begin by describing the general properties of real orthogonal polynomials with respect to a weight function, such as the recurrence relation and the properties of their roots. We will focus on studying the properties of orthogonal Jacobi polynomials and some of their particular cases (Legendre and Tchebyche polynomials). Then, we will present basic notions about quadrature formulas and nally we will focus on quadrature formulas of Gaussian type. In order to better understand these formulas we will present some examples of these quadratures for the Legendre and Tchebyche weights. vii 4 1. Generalidades sobre polinomios ortogonales reales Denición 1.9. Se dice que dos funciones f, g ∈ L2 ω(a, b) son ortogonales si se cumple que (f, g) = 0. Usaremos la notación f⊥g para expresar que f y g son ortogonales. De modo más general, si S es un subconjunto de L2 ω(a, b) , la notación f⊥S signica que f⊥g para todo g∈S . Denotaremos con Pn el espacio de polinomios con coecientes reales de grado menor o igual que n , es decir, Pn:= p:x∈R7→ p(x) = a0+a1x+. . . +an−1xn−1+anxn, ai∈R,0≤i≤n. También usaremos la notación deg(p) para referirnos al grado del polinomio p . Teorema 1.10. (Existencia y unicidad salvo una constante multiplicativa de los polinomios ortogonales) Existe una sucesión {pn}∞ n=o de polinomios mutuamente ortogonales con deg(pn) = n para todo n∈N . Además, si {˜pn}∞ n=o es otra sucesión de polinomios mutuamente ortogonales con deg(˜pn) = n , entonces para cada entero n≥0 existe una constante an∈R\ {0} tal que ˜pn=anpn . Demostración. Existencia. Tomamos la base canónica del espacio de los polinomios, {1, x, x2, . . .} , y aplicamos el algoritmo de ortogonalización de Gram-Schmidt. Sea p0= 1 . Buscamos p1(x) de la forma p1=x+a10p0 y tal que cumpla p0⊥p1 ; imponemos pues 0 = (p0, p1)=(p0, x) + a10(p0, p0), lo que equivale a a10 =−(p0,x) (p0,p0) ya que (p0, p0)>0 . Por construcción, p1 es de grado 1 y ortogonal a p0 . Supongamos construída la sucesión nita de polinomios ortogonales {p0, p1, ..., pn−1} tal que deg(pj) = j para todo j= 0,1, ..., n −1 . Buscamos pn de la forma pn=xn+ n−1 X i=0 anipi y tal que pn⊥pj para todo j= 0,1, . . . , n −1 . Imponemos pues 0 = (pn, pj)=(xn, pj) + n−1 X i=0 ani(pi, pj)=(xn, pj) + anj(pj, pj), j = 0,1, . . . , n −1. 5 Dado que (pj, pj)>0 ; podemos obtener los coecientes anj . Por construción, pn es mónico, de grado n y ortogonal al conjunto {p0, p1, . . . , pn−1} . De modo que salvo, el cálculo explícito de aij con 0≤j < i ≤n se establece la existencia de una sucesión de polinomios ortogonales. Unicidad. Supongamos que {˜pn}n n=0 es otra sucesión de polinomios ortogonales, tal que deg(˜pn) = n , para todo n≥0 . Nótese que como deg(pj) = j para todo j= 0,1, . . . , n , el conjunto {p0, p1, . . . , pn} es una base de Pn . Dado que ˜pn∈Pn , ˜pn se puede escribir de la forma ˜pn(x) = n X j=0 ajpj(x). Además, ˜pn⊥pj para j= 0,1, ..., n −1 , por lo que 0 = (˜pn, pj) = n X i=0 aipi, pj!=aj(pj, pj) =⇒aj= 0, j = 0,1, ..., n −1, obteniendo así que ˜pn(x) = anpn(x) ; an6= 0 ya que deg(˜pn) = n . Con esto queda probada la unicidad (salvo una constante multiplicativa) de una sucesión de polinomios ortogonales. Observación 1.11 . Con las notaciones del teorema 1.10, si además {pn}∞ n=0 y {˜pn}∞ n=0 son ortonormales , es decir, ||pn|| =||˜pn|| = 1 , entonces ˜pn(x) = ±pn(x) . En efecto, 1 = (˜pn,˜pn)=(anpn, anpn) = a2 n(pn, pn) = a2 n=⇒an=±1, entonces concluímos que ˜pn(x) = ±pn(x),∀n∈N . Corolario 1.12. (pn, p)=0 para todo p∈Pn−1 . Demostración. Es consecuencia inmediata de que (pn, pj)=0 , para todo j= 0,1, . . . , n −1 y <{p0, p1, . . . , pn−1}>=Pn−1 . A continuación veremos dos propiedades generales de los polinomios ortogonales. Teorema 1.13. La sucesión de polinomios ortogonales mónicos {pn}∞ n=0 con deg(pn) = n para todo n≥0 verica la relación de recurrencia pn+1(x)=(x−λn+1)pn(x)−µn+1pn−1(x),∀n≥0, donde p−1≡0, (1.2) λn+1 =(xpn, pn) (pn, pn), n ≥0 y µn+1 =(pn, pn) (pn−1, pn−1), n ≥1. 6 1. Generalidades sobre polinomios ortogonales reales Demostración. Como {pn}∞ n=0 es una sucesión de polinomios ortogonales, pn(x) es un polinomio ortogonal de grado n . Tenemos que demostrar que el polinomio mónico ortogonal de grado n+1 , pn+1(x) , verica (1.2). Dado que pn(x) también es mónico, pn+1 −xpn∈Pn y por tanto pn+1 se puede expresar de la forma pn+1(x) = xpn(x) + n X i=0 aipi(x) = (x+an)pn(x) + n−1 X i=0 aipi(x). Como (pn, pj)=0 para j= 0,1, ..., n −1 , tenemos que (pn+1, pn) = (x+an)pn+ n−1 X i=0 aipi, pn!= (xpn, pn) + an(pn, pn) = 0. Puesto que (pn, pn)>0 podemos despejar an , obteniendo que an=−(xpn, pn) (pn, pn). Si n= 0 , deniendo λ1=−a0 se obtiene (1.2) ( µ1 es irrelevante ya que p−1≡0 ). Supongamos pues, que n≥1 . Para j= 0,1, ..., n −1 , 0=(pn+1, pj) = (x+an)pn+ n−1 X i=0 aipi, pj!= (xpn, pj) + an(pn, pj) + n−1 X i=0 ai(pi, pj) = = (pn, xpj) + an(pn, pj) + n−1 X i=0 ai(pi, pj)=(pn, xpj) + aj(pj, pj). Si 0≤j≤n−2 , xpj tiene grado j+ 1 ≤n−1 ; por tanto (pn, xpj)=0 , luego aj= 0 . Consideremos ahora j=n−1 (pn, xpn−1)=(pn, pn)−(pn, pn−xpn−1)=(pn, pn) ya que deg(pn−xpn−1)≤n−1 . Por tanto, (pn, pn) + an−1(pn−1, pn−1) = 0 y como (pn−1, pn−1)>0 , obtenemos an−1=−(pn, pn) (pn−1, pn−1). Finalmente, (para n≥1 ) deniendo λn+1 := −an y µn+1 := −an−1 se obtiene el resultado. Teorema 1.14. Sea q(x) un polinomio de grado n tal que q⊥Pl , donde 0≤l≤n−1 . Entonces q tiene al menos l+ 1 raíces de multiplicidad impar en (a, b) . 7 Demostración. Veamos en primer lugar que q cambia de signo en (a, b) . Si q(x) es de signo positivo en (a, b) , entonces (en virtud de la condición (c) de la denición de función peso) (q, p0) = Zb a q(x)ω(x)dx > 0, pero esto es una contradicción, ya que (q, p0)=0 . Análogamente, si q(x) es de signo negativo en (a, b) , entonces (q, p0)<0 y se tiene una contradicción. Por consiguiente q cambia de signo en (a, b) ; por tanto, q(x) tiene al menos una raíz de multiplicidad impar en (a, b) . Supongamos ahora que q(x) sólo tiene k raíces de multiplicidad impar en el intervalo (a, b) , x1, ..., xk , con k < l + 1 . Sea r(x) = Qk j=1(x−xj) , entonces q(x)r(x) = q(x)(x−x1)···(x−xk) = sn−k(x)(x−x1)2···(x−xk)2, donde sn−k(x) no cambia de signo en (a, b) . Luego (q, r) = Zb a q(x)r(x)ω(x)dx =Zb a sn−k(x) k Y j=1 (x−xj)2ω(x)dx 6= 0 (1.3) ya que q(x)r(x) es un polinomio que no cambia de signo en (a, b) y ω es una función peso en (a, b) . Por otro lado, como q(x) es un polinomio ortogonal a Pl y deg(r) = k < l + 1 tenemos que (q, r) = 0 , en contradicción con (1.3). Por lo tanto k=l+ 1 y las raíces de q(x) son de multiplicidad impar en (a, b) . Corolario 1.15. Sea pn un polinomio ortogonal de grado n con n≥1 . Las raíces x1, ..., xn de pn son reales, simples y están en el intervalo (a, b) . Demostración. Es consecuencia inmediata del teorema anterior, ya que su hipótesis se verica con l=n−1 . Como además deg(pn) = n , pn tiene exactamente n raíces simples en (a, b) . 8 1. Generalidades sobre polinomios ortogonales reales Capítulo 2 Los polinomios ortogonales de Jacobi y sus propiedades elementales Los polinomios de Jacobi son los polinomios ortogonales para el peso ω(x) = (1 −x)α(1 + x)β en el intervalo (−1,1), donde α > −1 y β > −1. Presentaremos una serie de propiedades de los mismos, pero para mayor claridad, comenzaremos con los casos particulares de los polinomios de Legendre y después los de Tchebyche. Demostraremos sus propiedades elementales y enunciaremos las de los de Jacobi. Para elaborar este capítulo, hemos consultado las referencias bibliográcas [1] y [5]. 2.1. Polinomios de Legendre Denición 2.1. Los polinomios de Legendre , que denotaremos por Pn(x) , son los polinomios ortogonales entre sí con respecto a la función peso ω(x)=1 en el intervalo (−1,1) tales que para todo entero n≥0 , Pn(x) tiene grado n y satisface Pn(1) = 1 . Denotaremos por kn el coeciente principal de Pn(x) . Estos polinomios satisfacen las siguientes propiedades: 1. SIMETRÍA Pn(−x)=(−1)nPn(x),∀x∈R (2.1) 2. ECUACIÓN DIFERENCIAL (1 −x2)P00 n(x)−2xP0 n(x) + n(n+ 1)Pn(x) = 0 (2.2) 9 10 2. Los polinomios ortogonales de Jacobi y sus propiedades elementales 3. FÓRMULA DE RODRIGUES Pn(x) = (−1)n 2nn! dn dxn(1 −x2)n (2.3) 4. COEFICIENTE PRINCIPAL kn=(2n)! 2n(n!)2 (2.4) 5. NORMA Z1 −1 (Pn(x))2dx =2 2n+ 1 (2.5) 6. FÓRMULA DE RECURRENCIA (P0(x)=1, P1(x) = x (n+ 1)Pn+1(x) = (2n+ 1)xPn(x)−nPn−1(x), para todo n≥1. (2.6) 7. EXPRESIÓN EXPLÍCITA Pn(x) = 1 2n bn 2c X m=0 (−1)mn m2n−2m nxn−2m, (2.7) donde bsc denota la parte entera de s∈R , es decir, bsc=max{k∈Z/k ≤s} . Demostraremos estas propiedades: Simetría. Sea ϕ∈Pn−1 arbitrario. Como Pn es ortogonal a cualquier polinomio de grado ≤n−1 , haciendo el cambio de variable x=−y , tenemos que 0 = Z1 −1 Pn(x)ϕ(x)dx =−Z−1 1 Pn(−y)ϕ(−y)dy =Z1 −1 Pn(−y)ϕ(−y)dy. Por lo tanto, Pn(−x) también es ortogonal a Pn−1 . Además, por la unicidad (salvo una constante multiplicativa) del polinomio ortogonal de grado n , tenemos que Pn(−x) = cnPn(x) con cn∈R, cn6= 0. (2.8) Sabemos que podemos escribir Pn(x) de la forma Pn(x) = knxn+ términos de grado ≤n−1, con lo cual Pn(−x) = (−1)nknxn+ términos de grado ≤n−1; entonces (−1)nknxn=cnknxn. 2.1. Polinomios de Legendre 11 Por lo tanto cn= (−1)n ; subtituyendo cn en (2.8) obtenemos (2.1). Ecuación diferencial Sea ϕ∈Pn−1 arbitrario. El polinomio d dx n(1 −x2)P0 no tiene grado a lo sumo n . Por otra parte, integrando por partes dos veces obtenemos Z1 −1 d dx n(1 −x2)P0 no(x)ϕ(x)dx = (1 −x2)P0 n(x)ϕ(x) 1 −1−Z1 −1 (1 −x2)P0 n(x)ϕ0(x)dx = =−Z1 −1 (1 −x2)P0 n(x)ϕ0(x)dx = =−Pn(x)(1 −x2)ϕ0(x) 1 −1+Z1 −1 Pn(x)d dx n(1 −x2)ϕ0o(x)dx = =Z1 −1 Pn(x)d dx n(1 −x2)ϕ0o(x)dx = 0 ya que deg d dx n(1 −x2)ϕ0o≤n−1 y Pn⊥Pn−1 . Por lo tanto podemos escribir d dx n(1 −x2)P0 no como d dx n(1 −x2)P0 no=anPn, an∈R. (2.9) Además, como Pn(x) = knxn + términos de grado ≤n−1 , tenemos que d dx n(1 −x2)P0 no=d dx (1 −x2)(nknxn−1+ términos de grado ≤n−1)= =d dx −nknxn+1 + términos de grado ≤n)= =−n(n+ 1)knxn+ términos de grado ≤n−1 y comparando los coecientes principales en (2.9), deducimos que −n(n+ 1)kn=ankn. Por lo tanto an=−n(n+ 1) , y substituyendo en (2.9) obtenemos la ecuación (2.2). Fórmula de Rodrigues. Consideramos dn dxn(1 −x2)n , que es un polinomio de grado n . Veamos que es ortogonal a cualquier polinomio de grado n−1 . Sea ϕ∈Pn−1 arbitrario. Integrando por partes Z1 −1 dn dxn(1 −x2)nϕ(x)dx =dn−1 dxn−1(1 −x2)nϕ(x) 1 −1−Z1 −1 dn−1 dxn−1(1 −x2)nϕ0(x)dx = =−Z1 −1 dn−1 dxn−1(1 −x2)nϕ0(x)dx. 12 2. Los polinomios ortogonales de Jacobi y sus propiedades elementales Tras haber integrado por partes n veces, obtenemos Z1 −1 dn dxn(1 −x2)nϕ(x)dx = (−1)nZ1 −11−x2ndnϕ dxndx = 0 ya que dnϕ dxn= 0 debido a que ϕ∈Pn−1 . Por lo tanto, el polinomio dn dxn(1 −x2)n es ortogonal a Pn−1 ; como además tiene grado n , existe un número real bn tal que dn dxn(1 −x2)n=bnPn(x). (2.10) Tenemos que Pn(1) = 1 , por lo que bn=dn dxn(1 −x2)n(1) . Usando la fórmula de Leibnitz 1 obtenemos que dn dxn1−x2n= n X k=0 n kdk dxk{(1 −x)n}dn−k dxn−k{(1 + x)n}. Observamos que (1−x)n tiene una raíz en x= 1 de multiplicidad n , por lo tanto hdk dxk{(1 −x)n}i(1) = 0 , para k= 0,1, ..., n −1 , y por consiguiente dn dxn1−x2n(1) = 2ndn dxn{(1 −x)n}(1). Es fácil ver que dn dxn{(1 −x)n}= (−1)nn!. En consecuencia, bn= (−1)n2nn! . Substituyendo bn en (2.10) tenemos que dn dxn1−x2n= (−1)n2nn!Pn(x) y despejando Pn obtenemos (2.3), como queriamos demostrar. Coeciente principal Partimos de (2.3). Derivando n veces 1−x2n y jándonos en el coeciente del término de mayor grado, tenemos que dn dxn(−1)nx2n+ términos de grado ≤2n−1= = (−1)n2n(2n−1) ·. . . ·(n+ 1)xn+ términos de grado ≤n−1. 1 Para dos funciones f , g n veces diferenciables, la derivada de orden n cumple que dn dxn(fg) = n X k=0 n k!dkf dxk dn−kg dxn−k. 2.1. Polinomios de Legendre 13 Además, como 2n(2n−1) ·. . . ·(n+ 1) = (2n)! n! , obtenemos (2.4) como queriamos demostrar. Norma Usando (2.3) e integrando n veces por partes tenemos que Z1 −1 (Pn(x))2=Z1 −1 Pn(x)Pn(x)dx =(−1)n 2nn!Z1 −1 dn dxn1−x2nPn(x)dx =··· = =1 2nn!Z1 −11−x2ndnPn dxndx =1 2nn!Z1 −11−x2nknn!dx = =kn 2nn!Z1 −11−x2ndx =(2n)! 22n(n!)2Z1 −11−x2ndx. (2.11) Ahora tenemos que calcular R1 −11−x2ndx . Usando que el integrando es una función par y haciendo el cambio x= cos θ tenemos que Z1 −11−x2ndx = 2 Z1 01−x2ndx = 2 Z0 π/2 sin2nθ(−sin θ)dθ = 2 Zπ/2 0 (sin θ)2n+1 dθ. (2.12) Integrando por partes Zπ/2 0 (sin θ)2n+1 dθ =−cos θsin2nθ π/2 0+ 2nZπ/2 0 sin2n−1θcos2θdθ = =2nZπ/2 0 sin2n−1θ1−sin2θdθ = 2nZπ/2 0 sin2n−1θ−2nZπ/2 0 sin2n+1 θdθ. Denotando In=Rπ/2 0sin2n+1 θdθ y despejando obtenemos In=2n 2n+ 1In−1. (2.13) Cambiando n por n−1 , tenemos que In−1=2n−2 2n−1In−2 y remplazando en (2.13) In=2n(2n−2) (2n+ 1)(2n−1)In−2. Reiterando el procedimiento resulta que In=2n(2n−2) ·. . . ·2 (2n+ 1)(2n−1) ·. . . ·3I0. Por otra parte, I0=Zπ/2 0 sin θdθ =−cos θ|π/2 0= 1 y entonces In=2nn! (2n+ 1)(2n−1) ·. . . ·3. (2.14) 20 2. Los polinomios ortogonales de Jacobi y sus propiedades elementales que es un polinomio de grado ≤n . Un cálculo parecido prueba que para 0≤j≤n−1 , dj dxj(1 −x2)n−1 2 es una combinación lineal de (1 −x)n−k−1 2(1 + x)n−1 2−(j−k) , k= 0,1, . . . , j ; por consiguiente dj dxjn1−x2n−1 2o(±1) = 0. Veamos que Qn es ortogonal a Pn−1 . Sea ϕ∈Pn−1 arbitrario. Integrando por partes Z1 −1 Qn(x)ϕ(x)1−x2−1/2dx =Z1 −11−x21/2dn dxnn(1 −x2)n−1 2oϕ(x)1−x2−1/2dx = =Z1 −1 dn dxnn(1 −x2)n−1 2oϕ(x)dx =dn−1 dxn−1n(1 −x2)n−1 2oϕ(x) 1 −1− −Z1 −1 dn−1 dxn−1n(1 −x2)n−1 2oϕ0(x)dx =−Z1 −1 dn−1 dxn−1n(1 −x2)n−1 2oϕ0(x)dx. Tras haber integrado por partes n veces, obtenemos que Z1 −1 dn dxnn1−x2n−1 2oϕ0(x)dx = (−1)nZ1 −11−x2n−1 2dnϕ dxn= 0 ya que dnϕ dxn= 0 debido a que ϕ∈Pn−1 . Por lo tanto, Qn⊥Pn−1 y como además deg(Qn)≤n y Qn6≡ 0 , necesariamente Qn tiene grado n y lo podemos escribir como Qn(x) = 1−x21/2dn dxnn1−x2n−1 2o=bnTn(x), (2.28) donde bn∈R . Evaluando Qn(1) con ayuda de (2.27) y dado que (1 −x)n−k se anula en x= 1 para todo k= 0,1, . . . , n −1 , tenemos que bn=Qn(1) = (−1)n2nΓ(n+1 2) Γ(1 2). Substituyendo bn en (2.28), resulta 1−x21/2dn dxnn1−x2n−1 2o= (−1)n2nΓ(n+1 2) Γ(1 2)Tn(x) y despejando Tn(x) obtenemos (2.22) como queriamos demostrar. Simetría Procedemos por inducción utilizando la fórmula de recurrencia. 2.3. Polinomios de Jacobi 21 El enunciado es cierto para n= 0 y n= 1 , puesto que T0(x)=1 y T1(x) = x . Supongamos cierto el enunciado para todo natural ≤n y probémoslo para n+ 1 . Tenemos que Tn+1(−x) = −2xTn(−x)−Tn−1(−x) = −2x(−1)nTn(x)−(−1)n−1Tn−1(x) = =(−1)n+1 (2xTn(x)−Tn+1(x)) = (−1)n+1Tn+1(x). Por lo tanto, la propiedad de simetría se verica. Partiendo de T0(x) = 1 y T1(x) = x , y usando la fórmula de recurrencia podemos calcular sucesivamente los polinomios de Tchebyche, obteniendo así la siguiente lista: T0(x) = 1 T1(x) = x T2(x) = 2x2−1 T3(x) = 4x3−3x T4(x) = 8x4−8x2+ 1 T5(x) = 16x5−20x3+ 5x T6(x) = 32x6−48x4+ 18x2−1 2.3. Polinomios de Jacobi Consideramos el intervalo (−1,1) y la función peso ω(x) = (1 −x)α(1 + x)β donde α > −1 y β > −1 . Para s∈R y k entero, k≥0 , denotaremos por s k el coeciente binomial dado por la fórmula s k=s(s−1) ···(s−k+ 1) k! si k≥1, y s 0= 1. Denición 2.4. Los polinomios de Jacobi , que denotamos por P(α,β) n(x) , son los polinomios mutuamente ortogonales con respecto al peso ω(x) = (1 −x)α(1 + x)β con α > −1 y β > −1 en el intervalo (−1,1) tales que para todo entero n≥0 , P(α,β) n(x) tiene grado n y cumple la condición P(α,β) n(1) = n+α n . Denotamos por k(α,β) n el coeciente principal de P(α,β) n(x) . Los polinomios P(α,β) n(x) satisfacen las siguientes propiedades: 22 2. Los polinomios ortogonales de Jacobi y sus propiedades elementales 1. FÓRMULA DE RODRIGUES P(α,β) n(x) = (−1)n 2nn!(1 −x)−α(1 + x)−βdn dxnn(1 −x)n+α(1 + x)n+βo (2.29) 2. EXPRESIÓN EXPLÍCITA P(α,β) n(x) = 1 2n n X m=0 n+α mn+β n−m(x−1)n−m(x+ 1)m (2.30) 3. SIMETRÍA P(α,β) n(−x) = (−1)nP(β,α) n(x),∀x∈R (2.31) 4. COEFICIENTE PRINCIPAL k(α,β) n=1 2n2n+α+β n (2.32) 5. NORMA Y ORTOGONALIDAD Z1 −1 P(α,β) n(x)P(α,β) m(x)(1 −x)α(1 + x)βdx = =                  2α+β+1 (2n+α+β+1)n! Γ(n+α+1)Γ(n+β+1) Γ(n+α+β+1) si n =m y n6= 0 ó α+β+ 1 6= 0 2α+β+1 Γ(α+1)Γ(β+1) Γ(α+β+2) si n =m= 0 y α+β+ 1 = 0 0si n 6=m (2.33) 6. ECUACIÓN DIFERENCIAL ((1 −x2)y00 + (β−α−(α+β+ 2)x)y0+n(n+α+β+ 1)y= 0, y=P(α,β) n(x),∀n≥0 (2.34) 7. FÓRMULA DE RECURRENCIA 2(n+ 1)(n+α+β+ 1)(2n+α+β)P(α,β) n+1 (x) = (2n+α+β+ 1)[α2−β2+ (2n+α+ +β+ 2)(2n+α+β)x]P(α,β) n(x)−2(n+α)(n+β)(2n+α+β+ 2)P(α,β) n−1(x),∀n≥0 (2.35) La demostración de estas propiedades se puede ver en [5]. Cuando α=β estos polinomios se denominan polinomios ultraesféricos . Entre ellos tenemos algunos polinomios que tienen nombres especiales: 2.3. Polinomios de Jacobi 23 Cuando α=β= 0 , se llaman polinomios de Legendre ; los hemos estudiado anteriomente denotándolos por Pn(x) . Cuando α=β=−1 2 , se llaman polinomios de Tchebyche de primera especie . Dado que P(−1 2,−1 2) n(1) = n−1/2 n6= 1 = Tn(1) , obviamente P(−1 2,−1 2) n(1) 6=Tn(1) , pero como ambos polinomios son de grado n y ortogonales a Pn−1 , en virtud del Teorema 1.10 se tiene que P(−1 2,−1 2) n(x) = n−1/2 nTn(x),∀x∈R. Cuando α=β=1 2 , se llaman polinomios de Tchebyche de segunda especie . Para este peso es habitual considerar la sucesión de polinomios ortogonales {Un(x)}∞ n=0 tal que para todo entero n≥0 , deg(Un) = n y Un(1) = n+1 . En virtud del teorema 1.10, tenemos que P(1 2,1 2) n(x) = cnUn(x), con cn=n+1/2 n n+ 1 . 24 2. Los polinomios ortogonales de Jacobi y sus propiedades elementales Capítulo 3 Fórmulas de cuadratura de tipo gaussiano El objectivo de este capítulo es el estudio de los conceptos de las fórmulas de cuadratura de tipo gaussiano. Previamente, introducimos los conceptos de fórmulas de cuadratura, fórmulas de cuadratura de tipo interpolatorio polinómico y grado de precisión de una fórmula de cuadratura. Posteriormente, estudiaremos las fórmulas de cuadratura de tipo gaussiano (fórmulas de Gauss, Gauss-Radau y Gauss-Lobatto), que son aquellas que para un cierto número de nodos intentan dar un grado de precisión máximo. También veremos algunos ejemplos de estas cuadraturas gaussianas, para los casos del peso de Legendre y el de Tchebyche. Como referencias bibliográcas hemos usado [2], [3] y [6]. 3.1. Nociones básicas sobre fórmulas de cuadratura Intentamos aproximar la integral denida I(f) = Zb a f(x)ω(x)dx (3.1) donde f(x) e una función continua en el intervalo acotado [a, b] y ω es una función peso en (a, b) . Observación 3.1 . Si f es una función continua en [a, b] , entonces fω ∈ L1(a, b) , donde L1(a, b) = nf: (a, b)→R medible /Rb a|f(x)|ω(x)dx < +∞o . En efecto, la función peso ω es medible y ω≥0 en el intervalo acotado (a, b) . La condición (b) de la denición de función peso implica que Rb aω(x)<+∞ , es decir, que ω∈ L1(a, b) . Como f es continua en el compacto [a, b] , f está acotada, con lo cual fω ∈ L1(a, b) . Sean xn,0, xn,1, . . . , xn,n ∈[a, b] con xn,i 6=xn,j para i6=j . 25 26 3. Fórmulas de cuadratura de tipo gaussiano Denición 3.2. Se llama fórmula de cuadratura a una aproximación numérica de la integral (3.1), de la forma In(f) := n X i=0 wn,if(xn,i). (3.2) Los números wn,i se llaman coecientes o pesos de cuadratura y los puntos xn,i , i= 0,1, . . . , n se llaman puntos o nodos de cuadratura . Denimos el error de cuadratura como En(f) := I(f)−In(f). Otro concepto importante es el grado de precisión o de exactitud de una fórmula de cuadratura. Denición 3.3. Se dice que una fórmula de cuadratura tiene grado de precisión m si y sólo si En(xk) = 0 para k= 0,1, . . . , m y En(xm+1)6= 0 . Observación 3.4 . Si una fórmula de cuadratura tiene grado de precisión m , entonces la fórmula proporciona el valor exacto de la integral (3.1) para todo polinomio f de grado ≤m . Denición 3.5. Se dice que una fórmula de cuadratura es de tipo interpolatorio polinómico o de cuadratura interpolatoria si es una fórmula de cuadratura de la forma In(f) := Zb a Pn(x)ω(x)dx (3.3) siendo Pn(x) el polinomio de interpolación de Lagrange de f relativo a los nodos xn,0, xn,1, . . . , xn,n . Veamos que (3.3) es una fórmula de cuadratura. Representamos Pn(x) mediante la fórmula de interpolación de Lagrange: Pn(x) = n X i=0 f(xn,i)φn,i(x), donde los φn,i(x) , i= 0,1, . . . , n son los polinomios fundamentales de Lagrange relativos a los puntos xn,i , i= 0,1, . . . , n . Recordemos que φn,i(x) es el único polinomio de grado ≤n que cumple φn,i(xn,j) = δij, j = 0,1, . . . , n. Por lo tanto, Zb a Pn(x)ω(x)dx =Zb a n X i=0 f(xn,i)φn,i(x)ω(x)dx = n X i=0 Zb a φn,i(x)ω(x)dxf(xn,i) = n X i=0 wn,if(xn,i). Concluímos así que (3.3) es en efecto una fórmula de cuadratura y que la fórmula de cuadratura (3.2) es de tipo interpolatorio polinómico si y sólo si wn,i =Zb a φn,i(x)ω(x)dx, para i= 0,1, . . . , n. 3.1. Nociones básicas sobre fórmulas de cuadratura 27 Teorema 3.6. Una fórmula de cuadratura con n+ 1 nodos es de tipo interpolatorio polinómico si y sólo si es exacta en Pn . Demostración. Consideremos una fórmula de cuadratura con n+ 1 nodos Zb a f(x)ω(x)dx ≃In(f) = n X i=0 wn,if(xn,i). Supongamos que dicha fórmula de cuadratura es de tipo interpolatorio polinómico, es decir, wn,i =Zb a φn,i(x)ω(x)dx, i = 0,1, . . . , n y veamos que es exacta en Pn . Sea p∈Pn ; representamos p de la forma: p(x) = n X i=0 p(xn,i)φn,i(x); por tanto, Zb a p(x)ω(x)dx =Zb a n X i=0 p(xn,i)φn,i(x)ω(x)dx = n X i=0 p(xn,i)Zb a φn,i(x)ω(x)dx = = n X i=0 wn,ip(xn,i) = In(p). Recíprocamente, supongamos que la fórmula de cuadratura es exacta en Pn . Dado que para todo i= 0,1, . . . , n , φn,i ∈Pn , se cumple que Zb a φn,i(x)ω(x)dx =In(φn,i) = n X j=0 wn,jφn,i(xn,j) = n X i,j=0 wn,jδij =wn,i. Por consiguiente la fórmula de cuadratura es de tipo interpolatorio polinómico. Proposición 3.7. Si los nodos de cuadratura están simétricamente distribuídos respecto al punto a+b 2 , esto es, xn,i −a+b 2=−xn,n−i−a+b 2, i = 0,1,...,jn 2k, (3.4) y el peso es simétrico respecto de a+b 2 , es decir, ωt+a+b 2=ωa+b 2−t,∀t∈0,b−a 2, (3.5) entonces los coecientes de la fórmula de cuadratura de tipo interpolatorio polinómico con nodos xn,i , i= 0,1, . . . , n , satisfacen wn,i =wn,n−i, i = 0,1,...,jn 2k. 28 3. Fórmulas de cuadratura de tipo gaussiano Demostración. Tenemos que demostrar primero que hay una relación de simetría entre los polinomios fundamentales φn,i , i=,1, . . . , n . Recordemos que φn,i es el polinomio de grado ≤n tal que φn,i(xn,j) = δij para todo j= 0,1, . . . , n . La relación (3.4) se puede escribir de la forma xn,i =a+b−xn,n−i, i = 0,1, . . . , n. Consideramos ψi(x) := φn,i(a+b−x) , que es un polinomio de grado ≤n ; tenemos que ψi(xn,j) := φn,i(a+b−xn,j) = φn,i(xn,n−j) = δi,n−j=δn−i,j,0≤j≤n. Por la unicidad del polinomio de interpolación tenemos que ψi(x) = φn,n−i(x),∀x∈R,∀i= 0,1, . . . , n es decir, φn,i(a+b−x) = φn,n−i(x),∀x∈R,∀i= 0,1, . . . , n. Por otra parte, haciendo el cambio de variable τ=t+a+b 2 , podemos escribir la relación (3.5) como ω(τ) = ω(a+b−τ),∀τ∈(a, b). Por lo tanto, si hacemos el cambio de variable x=a+b−t , resulta que wn,i =Zb a φn,i(x)ω(x)dx =−Zb a φn,i(a+b−t)ω(a+b−t)dt =Zb a φn,n−i(t)ω(t)dt =wn,n−i, obteniendo así la simetría de los coecientes como queriamos demostrar. 3.2. Fórmulas de cuadratura de Gauss, Gauss-Radau y GaussLobatto Ahora nos centraremos en el estudio de las fórmulas de Gauss, Gauss-Radau y Gauss-Lobatto. Describiremos su construcción, estudiaremos cuál es el grado de exactitud que alcanzan y la positividad de sus coecientes. 3.2.1. Fórmulas de Gauss Sea (a, b) un intervalo (nito o innito) de R y ω una función peso en (a, b) . Sea n≥0 un entero. Nos planteamos el problema de encontrar una fórmula de cuadratura con n+ 1 nodos Zb a f(x)ω(x)dx ≃ n X i=0 wn,if(xn,i) (3.6) 3.2. Fórmulas de cuadratura de Gauss, Gauss-Radau y Gauss-Lobatto 29 cuyo grado de exactitud sea el máximo posible. Veremos que ese máximo es 2n+ 1 y que existe una única fórmula de cuadratura de n+ 1 nodos con grado de exactitud 2n+ 1 . Dicha fórmula se llama fórmula de Gauss de n+ 1 nodos. Observación 3.8 . Si una fórmula de cuadratura con n+ 1 nodos es exacta en P2n+1 , entonces es una fórmula de tipo interpolatorio polinómico. En efecto, como obviamente n≤2n+ 1 , Pn⊂P2n+1 ; por lo tanto, si la fórmula de cuadratura es exacta en P2n+1 , también lo es en Pn y por consiguiente es de tipo interpolatorio polinómico. Proposición 3.9. Si la fórmula de cuadratura (3.6) es exacta en P2n+1 , entonces sus n+ 1 nodos, que denotamos por xn,i , i= 0,1, . . . , n son las raíces de pn+1(x) , el polinomio ortogonal de grado n+ 1 relativo a la función peso ω . Demostración. Sea ψ(x) = Qn i=0 (x−xn,i) ; el grado de ψ es n+1 . Probaremos que ψ es ortogonal a Pn . Sea s∈Pn arbitario. El producto ψ(x)s(x) tiene grado a lo sumo 2n+ 1 , es decir, ψs ∈P2n+1 . Utilizando la hipótesis de que la fórmula de cuadratura es exacta en P2n+1 , tenemos que Zb a ψ(x)s(x)ω(x)dx = n X j=0 wn,jψ(xn,j)s(xn,j) = n X j=0 wn,js(xn,j) n Y i=0 (xn,j −xn,i) = 0. Por lo tanto, para todo s∈Pn , Rb aψ(x)s(x)ω(x)dx = 0 , es decir, ψ⊥Pn . Por la unicidad salvo una constante multiplicativa del polinomio ortogonal de grado n+1 , tenemos que pn+1(x) = cnψ(x) , con cn∈R , cn6= 0 . Si consideramos pn+1(x) mónico, entonces cn= 1 . De la Observación 3.8 y la Proposición 3.9 se deduce que si existe una fórmula de cuadratura con n+ 1 nodos exacta en P2n+1 es única. En efecto, sus nodos xn,i , i= 0,1, . . . , n , han de ser las raíces de pn+1(x) y como además la fórmula de cuadratura ha de ser de tipo interpolatorio polinómico, una vez determinados los xn,i , i= 0,1, . . . , n , los coecientes wn,i , i= 0,1, . . . , n vendrán dados por wn,i =Zb a φn,i(x)ω(x)dx, i = 0,1, . . . , n (3.7) siendo φn,i , i= 0,1, . . . , n los polinomios fundamentales de Lagrange relativos a los puntos xn,0, . . . , xn,n (es decir, φn,i ∈Pn , φn,i(xn,j) = δij , j= 0,1, . . . , n ). Esto establece la unicidad de la fórmula de n+ 1 puntos exacta en P2n+1 . Veremos ahora la existencia de dicha fórmula. Proposición 3.10. Sean xn,0,...xn,n las raíces de pn+1(x) (recordemos que son todas distintas y que están en (a,b)) y sean wn,0, . . . , wn,n los coecientes dados por la expresión (3.7). Entonces 36 3. Fórmulas de cuadratura de tipo gaussiano ya que ψ(xn,i) = 0 , i= 0,1, . . . , n . Por tanto, ψ es ortogonal a Pn−2 , y como además ψ∈Pn+1 lo podemos escribir como ψ(x) = αpn+1(x) + βpn(x) + γpn−1(x) con α, β, γ ∈R. Como ψ y pn+1 son mónicos, α= 1 . Además, xn,0=−1 y xn,n = 1 ; por tanto ψ(−1) = 0 y ψ(1) = 0 , y entonces β=a y γ=b . Con esto obtenemos que ψ(x) = pn+1(x) + apn(x) + bpn−1(x) = q(x). Luego los nodos xn,i , i= 1, . . . , n −1 con xn,0=−1 y xn,n = 1 son las raíces de q(x) como queriamos demostrar. Con la Observación 3.19 y la Proposición 3.20 queda determinada la unicidad de la fórmula. En efecto, los n+ 1 nodos, dos de ellos el −1 y el 1 , son las raíces de q(x) y como la fórmula de cuadratura es de tipo interpolatorio polinómico, quedan determinados sus coecientes por (3.7). Determinaremos ahora la existencia de la fórmula. Proposición 3.21. Sean −1 = xn,0< xn,1< . . . < xn,n−1< xn,n = 1 las n+ 1 raíces de q y sean wn,i , i= 0,1, . . . , n los coecientes dados por (3.7). Entonces la fórmula de cuadratura Z1 −1 f(x)ω(x)dx ≃ n X i=0 wn,if(xn,i) es exacta en P2n−1 . Demostración. Sea p∈P2n−1 arbitrario. Haciendo la división euclídea de p entre q , tenemos que existen s y r únicos con deg(s)≤n−2 y deg(r)≤n tales que p=qs +r. (3.20) Tenemos que probar que Z1 −1 p(x)ω(x)dx = n X i=0 wn,ip(xn,i). (3.21) Tenemos que Z1 −1 p(x)ω(x)dx =Z1 −1 q(x)s(x)ω(x)dx +Z1 −1 r(x)ω(x)dx =Z1 −1 r(x)ω(x)dx (3.22) ya que R1 −1q(x)s(x)ω(x)dx = 0 , porque q⊥Pn−2 y s∈Pn−2 . Como la fórmula de cuadratura (3.8) es interpolatoria polinómica, es exacta en Pn . Como además r∈Pn , tenemos que Z1 −1 r(x)ω(x)dx = n X i=0 wn,ir(xn,i). (3.23) 3.2. Fórmulas de cuadratura de Gauss, Gauss-Radau y Gauss-Lobatto 37 Dado que los xn,i son las raíces de q(x) , tenemos que q(xn,i)=0, i = 0,1, . . . , n lo que, junto con (3.20), implica p(xn,i) = r(xn,i), i = 0,1, . . . , n. (3.24) Finalmente, usando (3.22), (3.23) y (3.24) obtenemos (3.21). Concluímos así que la fórmula de cuadratura es exacta en P2n , como queríamos demostrar. Estudiaremos ahora la optimalidad de la fórmula. Proposición 3.22. No existe ninguna fórmula de cuadratura de n+ 1 nodos con xn,0=−1 y xn,n = 1 que sea exacta en P2n . Demostración. Sean −1 = xn,0< xn,1< . . . < xn,n−1< xn,n = 1 los nodos y wn,i , i= 0,1, . . . , n los coecientes de la fórmula. Sea s(x) = (1 + x)(1 −x)Qn−1 j=1 (x−xn,j)2 ; deg(s) = 2n . Probemos que la fórmula de cuadratura no es exacta para s . Tenemos que s(x)6≡ 0 y s(x)≥0 , ∀x∈[−1,1] , entonces por la condición (c) de la denición de función peso obtenemos que Z1 −1 s(x)ω(x)dx > 0. Por otro lado, como s(xn,i)=0 , i= 0,1, . . . , n , entonces n X i=0 wn,is(xn,i)=0. Por lo tanto, la fórmula de cuadratura no es exacta en P2n . Proposición 3.23. wn,i >0 , i= 0,1, . . . , n . Demostración. Para i6= 0 , consideramos R1 −1 (φn,i(x))2 1+xω(x)dx ; nótese que ψi(x) := (φn,i(x))2 1+x es un polinomio de grado 2n−1 y su factorización contiene el factor 1 + x , que es positivo en (−1,1) ; aplicando la condición (c) de la denición de función peso, tenemos que Z1 −1 (φn,i(x))2 1 + xω(x)dx = n X j=0 wn,jψi(xn,j) = wn,i 1 + xn,i >0. Como 1 + xn,i >0 , concluímos que wn,i >0 , para todo i= 1, . . . , n . 38 3. Fórmulas de cuadratura de tipo gaussiano Para i= 0 , consideramos R1 −1 (φn,0(x))2 1−xω(x)dx ; análogamente ψ0(x) := (φn,0(x))2 1−x es un polinomio de grado 2n−1 y ψ0(x)≥0 para todo x∈[−1,1] ; ahora tenemos que Z1 −1 (φn,0(x))2 1−xω(x)dx = n X j=0 wn,jψi(xn,j) = wn,0 1−xn,0 =wn,0 2>0. En el caso especial de los pesos de Jacobi, hay una caracterización alternativa de los nodos de la fórmula de cuadratura de Gauss-Lobatto de n+ 1 puntos (n≥1) : son −1 , 1 y las raíces de p0 n(x) . Nótese que, dado que pn(x) tiene n raíces distintas en (−1,1) , aplicando el Teorema de Rolle se deduce que las n−1 raíces de p0 n(x) son todas distintas y están en (−1,1) . La proposición siguiente justica la caracterización alternativa mencionada. Proposición 3.24. Sea n≥1 entero y sean xn,0=−1 , xn,n = 1 y xn,i , i= 1, . . . , n −1 las raíces de p0 n(x) . La fórmula de cuadratura interpolatoria polinómica Z1 −1 f(x)(1 −x)α(1 + x)βdx ≃ n X i=0 wn,if(xn,i) (3.25) es exacta en P2n−1 . Demostración. Sea p∈P2n−1 arbitrario. Haciendo la división euclídea de p entre el polinomio (1 −x2)p0 n(x) , obtenemos que existen polinomios s y r únicos con deg(s)≤n−2 y deg(r)≤n tales que p= (1 −x2)p0 ns+r. (3.26) Tenemos que ver que Z1 −1 p(x)ω(x)dx = n X i=0 wn,ip(xn,i). (3.27) Tenemos que Z1 −1 p(x)ω(x)dx =Z1 −1 (1 −x2)p0 n(x)s(x)ω(x)dx +Z1 −1 r(x)ω(x)dx. (3.28) Por una parte, integrando por parte Z1 −1 (1 −x2)p0 n(x)s(x)ω(x)dx = = (1 −x2)s(x)pn(x)ω(x) 1 −1−Z1 −1d dx (1 −x2)s(x)ω(x) + (1 −x2)s(x)ω0(x)pn(x)dx = =−Z1 −1 pn(x)d dx (1 −x2)s(x)ω(x)dx −Z1 −1 pn(x)s(x)(1 −x2)ω0(x) ω(x)ω(x)dx (3.29) 3.3. Ejemplos de fórmulas de cuadratura de tipo gaussiano 39 ya que (1−x2)ω(x) = (1−x)α+1(1+x)β+1 se anula en ±1 (ya que α+1 >0 y β+1 >0 ). Además en la primera integral, tenemos que deg((1−x2)s(x)) ≤n , con lo cual deg d dx (1 −x2)s(x)≤n−1 y como pn⊥Pn−1 , Z1 −1 pn(x)d dx (1 −x2)s(x)ω(x)dx = 0. (3.30) Consideremos ahora la segunda integral. Recordemos que ω(x) = (1 −x)α(1 + x)β ; es fácil comprobar que (1 −x2)ω0(x) ω(x)=−α(1 + x) + β(1 −x), con lo cual deg (1 −x2)ω0(x) ω(x)s(x)≤n−1 y como pn⊥Pn−1 , obtenemos que Z1 −1 pn(x)s(x)(1 −x2)ω0(x) ω(x)ω(x)dx = 0. (3.31) Usando (3.28)-(3.31), deducimos que Z1 −1 p(x)ω(x)dx =Z1 −1 r(x)ω(x)dx. (3.32) Por otra parte, como la fórmula de cuadratura es interpolatoria polinómica, es exacta en Pn . Como además r∈Pn , tenemos que Z1 −1 r(x)ω(x)dx = n X i=0 wn,ir(xn,i). (3.33) Dado que los xn,i , i= 1, . . . , n −1 son las raíces de p0 n(x) , tenemos que p0 n(xn,i) = 0, i = 1, . . . , n −1 lo que, junto con xn,0=−1 , xn,n = 1 y (3.26), implica p(xn,i) = r(xn,i), i = 0,1, . . . , n. (3.34) Finalmente, usando (3.32)-(3.34) obtenemos (3.27) como queríamos demostrar. Dado que la fórmula de cuadratura (3.25) tiene n+ 1 nodos, incluyendo −1 y 1 , y es exacta en P2n−1 , dicha fórmula es la fórmula de cuadratura de Gauss-Lobatto de n+ 1 puntos. 3.3. Ejemplos de fórmulas de cuadratura de tipo gaussiano En esta sección veremos algunas cuadraturas de tipo gaussiano para los pesos de Legendre y de Tchebyche. 40 3. Fórmulas de cuadratura de tipo gaussiano LEGENDRE-GAUSS Veamos cuál es la fórmula de cuadratura de Gauss de 3 nodos en [−1,1] para el peso de Legendre. Tenemos que Z1 −1 f(x)dx ≃ 2 X i=0 w2,if(x2,i) = w2,0f(x2,0) + w2,1f(x2,1) + w2,2f(x2,2). (3.35) Los nodos son las raíces del polinomio P3(x) = 1 25x3−3x . P3(x)=0 ⇐⇒ 1 25x3−3x= 0 ⇐⇒ x= 0 y x=±r3 5 Por tanto, los nodos de la fórmula son: x2,0=−r3 5, x2,1= 0, x2,2= +r3 5 Ahora tenemos que calcular los pesos correspondientes y para eso imponemos exactitud en P2 tomando como base el conjunto 1, x, x2 . Para P0(x)=1 w2,0+w2,1+w2,2=Z1 −1 1dx = 2. Para P1(x) = x −r3 5w2,0+0+r3 5w2,2=Z1 −1 xdx = 0. Para P2(x) = x2 3 5w2,0+0+3 5w2,2=Z1 −1 x2dx =2 3. Resolviendo el sistema obtenemos que w2,0=w2,2=5 9 y w2,1=8 9. Finalmente substituyendo en (3.35) los nodos y pesos obtenidos tenemos que Z1 −1 f(x)dx ≃5 9f(−r3 5) + 8 9f(0) + 5 9f(r3 5) es la fórmula de cuadratura de Legendre-Gauss de 3 nodos en [−1,1] . 3.3. Ejemplos de fórmulas de cuadratura de tipo gaussiano 41 LEGENDRE-GAUSS-LOBATTO Veamos cuál es la fórmula de cuadratura de Gauss-Lobatto de 4 nodos en [−1,1] para el peso de Legendre. Tenemos que Z1 −1 f(x)dx ≃ 3 X i=0 w3,if(x3,i) = w3,0f(x3,0) + w3,1f(x3,1) + w3,2f(x3,2) + w3,3f(x3,3). (3.36) Como la fórmula es de tipo Gauss-Lobatto tenemos que dos nodos son los extremos y los otros dos nodos están en el interior del intervalo [−1,1] . Como además los polinomios de Legendre son un caso particular de los polinomios de Jacobi, los nodos interiores son las raíces de P0 3(x) = 15 2x2−3 2 . P0 3(x) = 0 ⇐⇒ 15 2x2−3 2= 0 ⇐⇒ x2=1 5⇐⇒ x=±r1 5 Por tanto, los nodos de la fórmula son: x3,0=−1, x3,1=−r1 5, x3,2= +r1 5, x3,3= 1 Ahora tenemos que calcular los pesos correspondientes y para eso imponemos exactitud en P3 tomando como base el conjunto 1, x, x2, x3 . Para P0(x)=1 w3,0+w3,1+w3,2+w3,3=Z1 −1 1dx = 2. Para P1(x) = x −w3,0−r1 5w3,1+r1 5w3,2+w3,3=Z1 −1 xdx = 0. Para P2(x) = x2 w3,0+1 5w3,1+1 5w3,2+w3,3=Z1 −1 x2dx =2 3. Para P3(x) = x3 −w3,0−1 5√5w3,1+1 5√5w3,2+w3,3=Z1 −1 x3dx = 0. Además, como los nodos están simétricamente distribuídos y el peso es simétrico respecto del origen, aplicando la proposición 3.7, tenemos simetría en los pesos, de modo que w3,0=w3,3 y w3,1=w3,2 42 3. Fórmulas de cuadratura de tipo gaussiano y así basta con resolver el sistema 2w3,0+ 2w3,1= 2 2w3,0+2 5w3,1=2 3) de donde obtenemos que w3,0=w3,3=1 6 y w3,1=w3,2=5 6. Finalmente substituyendo en (3.36) los nodos y pesos obtenidos tenemos que la fórmula de cuadratura de Legendre-Gauss-Lobatto de 4 nodos en [−1,1] es: Z1 −1 f(x)dx ≃1 6f(−1) + 5 6f(−r1 5) + 5 6f(r1 5) + 1 6f(1). En este epígrafe y en los dos siguientes, para abreviar la notación suprimimos el índice n en los nodos y coecientes de las fórmulas de cuadratura, de modo que consideraremos xi y wi, i = 0,1, . . . , n. TCHEBYCHEFF-GAUSS Consideramos n+ 1 nodos, donde n es el grado de la interpolación polinómica subyacente. Veamos que Z1 −1 f(x)1 √1−x2dx ≃In(f) = π n+ 1 n X j=0 fcos (2j+ 1)π 2n+ 2  es la fórmula de Tchebyche-Gauss de n+ 1 nodos. Los n+ 1 nodos son las raíces de Tn+1(x) . Recordemos que Tn+1(x) = cos((n+ 1) arc cos(x)) . Haciendo el cambio de variable θ= arc cos(x) , tenemos que x= cos(θ) y Tn+1(x) = cos((n+ 1)θ). Por lo tanto, Tn+1(x)=0 ⇐⇒ cos((n+ 1)θ) = 0 ⇐⇒ (n+ 1)θ=π 2+jπ, j ∈Z. (3.37) Además como x∈[−1,1] , tenemos que θ∈[0, π] , con lo cual (n+ 1)θ∈[0,(n+ 1)π] . De esto, junto con (3.37), deducimos que (n+ 1)θ=π 2+jπ, j = 0,1, . . . , n. 3.3. Ejemplos de fórmulas de cuadratura de tipo gaussiano 43 Introducimos θj , que son los valores de θ correspondientes: θj= π 2+jπ n+ 1 =(2j+ 1)π 2n+ 2 , j = 0,1, . . . , n. Deshaciendo el cambio de variable tenemos que las raíces de Tn+1(x) son xj= cos(θj) = cos (2j+ 1)π 2n+ 2 , para j= 0,1, . . . , n. Además como la función coseno en [0, π] decrece estrictamente, las raíces son distintas y están ordenadas de forma decreciente, es decir, 1> x0> x1> . . . > xn>−1. Comprobemos ahora que los pesos son wj=π n+1 para j= 0,1, . . . , n ; para esto basta comprobar que la fórmula es exacta en P2n+1 . Tomamos como base de P2n+1 el conjunto {T0, T1, . . . , T2n+1} . Nótese que como T0= 1 , en virtud de (2.20) tenemos que Z1 −1 Tk(x)1 √1−x2dx =Z1 −1 Tk(x)T0(x)1 √1−x2dx =       0para k ≥1 π para k = 0 Ahora calcularemos el resultado de aplicar la fórmula de cuadratura a los Tk . Para k= 0 , In(T0) = In(1) = π n+ 1(n+ 1) = π. Para 1≤k≤2n+ 1 , In=π n+ 1 n X j=0 Tkcos (2j+ 1)π 2n+ 2 . (3.38) Tenemos que n X j=0 Tkcos (2j+ 1)π 2n+ 2 = n X j=0 cos (2j+ 1)kπ 2n+ 2 . (3.39) Como los argumentos del coseno están en progresión aritmética, transformamos la suma en una progresión geométrica, cuya razón es ekπi n+1 6= 1 ; tenemos que n X j=0 cos (2j+ 1)kπ 2n+ 2 =Re   n X j=0 e(2j+1)kπi 2n+2  =Re   e(2n+1)kπi 2n+2 ·e2kπi 2n+2 −ekπi 2n+2 e2kπi 2n+2 −1 = =Re   e(2n+3)kπi 2n+2 −ekπi 2n+2 e2kπi 2n+2 −1 =Re ekπi ·ekπi 2n+2 −ekπi 2n+2 e2kπi 2n+2 −1!= =Re (−1)k−1ekπi 2n+2 e2kπi 2n+2 −1! 44 3. Fórmulas de cuadratura de tipo gaussiano y multiplicando numerador y denominador por e −kπi 2n+2 obtenemos que Re (−1)k−1ekπi 2n+2 e2kπi 2n+2 −1!=Re (−1)k−1 ekπi 2n+2 −e −kπi 2n+2 !=Re   (−1)k−1 2isin kπ 2n+2 = 0. Esto, junto con (3.38) y (3.39) implica que In(Tk)=0 para todo k= 1,...,2n+ 1 . Luego la fórmula es exacta en P2n+1 y como tiene n+ 1 nodos, en efecto es la fórmula de Tchebyche-Gauss de n+ 1 nodos. TCHEBYCHEFF-GAUSS-LOBATTO Consideramos n+ 1 nodos. Veamos que Z1 −1 f(x)1 √1−x2dx ≃In(f) = π 2n[f(−1) + f(1)] + π n n−1 X j=1 fcos jπ n es la fórmula de Tchebyche-Gauss-Lobatto de n+ 1 nodos. Como la fórmula es de tipo Gauss-Lobatto tenemos dos nodos que corresponden a los extremos y los n−1 nodos restantes están en el interior del intervalo [−1,1] . Como además los polinomios de Tchebyche son un caso particular de los polinomios de Jacobi, los nodos interiores son las raíces de T0 n(x) . Tenemos que T0 n(x) = n √1−x2sin (narc cos(x)) y haciendo el cambio de variable θ= arc cos(x) , tenemos que x= cos(θ) y T0 n(x) = nsin(nθ) sin(θ). Por lo tanto, T0 n(x) = 0 ⇐⇒ nsin(nθ) sin(θ)= 0 =⇒sin(nθ) = 0 ⇐⇒ nθ =jπ, j ∈Z. (3.40) Además como x∈(−1,1) , tenemos que θ∈(0, π) , con lo cual nθ ∈(0, nπ) . De esto, junto con (3.40), deducimos que nθ =jπ, j = 1, . . . , n −1. Introducimos θj , que son los valores de θ correspondientes: θj=jπ n, j = 1, . . . , n −1. 3.3. Ejemplos de fórmulas de cuadratura de tipo gaussiano 45 Deshaciendo el cambio de variable tenemos que las raíces de T0 n(x) son xj= cos(θj) = cos jπ n, para j= 1, . . . , n −1. Por lo tanto, los nodos de la fórmula son xj= cos jπ n , j= 0, . . . , n . Comprobemos ahora que los pesos son wj=       π 2npara j = 0, n π npara j = 1, . . . , n −1 Para esto basta comprobar que la fórmula de cuadratura es exacta en P2n−1 . Tomamos como base de P2n−1 el conjunto {T0, T1, . . . , T2n−1} . Recuérdese que Z1 −1 Tk(x)1 √1−x2dx =       0para k ≥1 π para k = 0 Ahora calcularemos el resultado de aplicar la fórmula de cuadratura a los Tk . Para k= 0 , In(T0) = In(1) = π 2n2 + π n(n−1) = π. Para 1≤k≤2n−1 , tenemos que In(Tk) = π 2n[Tk(−1) + Tk(1)] + π nS, (3.41) donde S= n−1 X j=1 Tkcos jπ n= n−1 X j=1 cos kjπ n. Como los argumentos del coseno están en progresión aritmética, transformamos la suma en una progresión geométrica, cuya razón es ekπi n6= 1 ; tenemos que S= n−1 X j=1 cos kjπ n=Re   n−1 X j=1 ekjπi n =Re e(n−1)kπi n·ekπi n−ekπi n ekπi n−1!= =Re ekπi −ekπi n ekπi n−1!=Re (−1)k−ekπi n ekπi n−1! y multiplicando numerador y denominador por e−kπi 2n obtenemos que S=Re (−1)ke−kπi 2n−ekπi 2n ekπi 2n−e−kπi 2n!=Re (−1)k−isin kπ 2n−isin kπ 2n isin kπ 2n+isin kπ 2n!= =Re (−1)k+1 −1isin kπ 2n 2isin kπ 2n!=Re (−1)k+1 −1 2=(−1)k+1 −1 2. (3.42)