Gaussian Quadrature September 25, 2011 Interpolation Approximation - PowerPoint PPT Presentation

Interpolation Gaussian Quadrature September 25, 2011

Interpolation Approximation of integrals Approximation of integrals by quadrature ◮ Many definite integrals cannot be computed in closed form, and must be approximated numerically. ◮ Basic building block � 1 N � u ( ξ ) w ( ξ ) d ξ ≈ u ( ξ i ) ω i (1) − 1 i = 0 ◮ ξ i are the quadrature points ◮ ω i are the quadrature weights ◮ w ( ξ ) > 0 is a weight function ◮ Integration can be done with function evaluations at ξ i ◮ How do I pick ξ i and ω i to minimum errors and/or maximize efficiency

Interpolation Approximation of integrals Example Quadrature: Trapezoidal integration ◮ w = 1 ◮ Quadrature points are equally spaced: ξ i = 2 i N − 1 ◮ Quadrature weights are � 2 for i = 0 , N 2 N ω i = (2) 2 for 1 ≤ i < N N ◮ Formula takes the form, with u i = u ( ξ i ) , � 1 u ( ξ ) d ξ ≈ 2 � u 0 2 + u 1 + · · · + u N − 1 + u N � (3) N 2 − 1 ◮ Can we do better?

Interpolation Approximation of integrals Possible route ◮ Choose ξ i and vary ω i to maximize order of exact integration. ◮ I can impose ( N + 1 ) constraints on ω i , e.g. maximize the degree of polynomials that quadrature approximates exactly: � 1 � ξ k d ξ, 0 ≤ k ≤ N ξ k i ω k = (4) − 1 i = 0 ◮ Gauss quadrature adjust ξ i so that formula above applies for k > ( N + 1 ) . Relies on properties of orthogonal polynomials.

Interpolation Approximation of integrals Gauss quadrature Let p m ( x ) be the set of orthogonal polynomials on the interval a ≤ x ≤ b w.r.t. weight function w ( x ) and of degree m � b p m ( x ) p n ( x ) w ( x ) d x = δ m , n � p m � 2 (5) a Let the collocation points be the roots of p N + 1 ( x i ) = 0 and the weights: � b � x k x k w ( x ) d x , 0 ≤ k ≤ N i ω k = (6) a i = 0 then...

Interpolation Approximation of integrals Gauss quadrature ◮ The weights are all positive ω i > 0. ◮ Gauss quadrature is exact for all polynomials, q , of degree less or equal to 2 N + 1 � b � q ( x i ) ω k = q ( x ) w ( x ) d x (7) a i = 0 ◮ It is not possible to find a x i , ω i combination where the integration is exact for polynomials of degree 2 N + 2.

Interpolation Approximation of integrals Gauss quadrature flavors Gauss quadrature comes in different flavors depending on whether the quadrature points include none, one, or both end points. ◮ Gauss quadrature: all quadrature points are strictly inside the interval: a < x i < b ◮ Gauss-Lobatto quadrature: x 0 = a and x N = b . Useful for imposing BC as we will see later. ◮ Gauss-Radau quadrature: x 0 = a or x N = b but not both. Useful for integrals or PDEs on semi-infinite intervals

Interpolation Approximation of integrals Gauss-Lobatto quadrature ◮ Formula is exact for polynomials of degree ≤ 2 N − 1. ◮ For Jacobi type polynomials, quadrature points are the roots of ( 1 − x 2 ) p ′ N ( x ) = 0 Jacobi polynomials include Chebyshev, Legendre and associated Legendre polynomials.

Interpolation Examples Example: Chebyshev Gauss Quadrature ◮ Roots are known analytically since T N + 1 ( x ) = cos ( N + 1 ) θ , with x = cos θ . ◮ The roots are then ξ i = cos θ i = ( 2 i + 1 ) π 2 ( N + 1 ) , 0 ≤ i ≤ N (8) ◮ Quadrature weights are π ω i = (9) N + 1 ◮ The formula is exact for polynomial of degree 2 N + 1: � 1 N q ( x ) � √ 1 − x 2 d x = q ( x i ) ω i (10) − 1 i = 0

Interpolation Examples Example: Chebyshev Gauss-Lobatto Quadrature N ( x ) = N sin N θ ◮ Roots are known analytically since T ′ sin θ ◮ The roots of ( 1 − x 2 ) p ′ N ( x ) = sin θ sin N θ ξ i = i π N , 0 ≤ i ≤ N (11) ◮ Quadrature weights are ω i = π 2 N for i = 0 , N , and ω i = π N , 1 ≤ i < N (12) ◮ The formula is exact for polynomial of degree 2 N − 1: � 1 N q ( x ) � √ 1 − x 2 d x = q ( x i ) ω i (13) − 1 i = 0

Interpolation Examples Example: Legendre Gauss Quadrature ◮ Roots must be computed numerically L N + 1 ( x i ) = 0 ◮ Quadrature weights are 2 ω i = (14) ( 1 − x i ) 2 L ′ N + 1 ( x i ) ◮ The formula is exact for polynomial of degree 2 N + 1: � 1 N � q ( x ) d x = q ( x i ) ω i (15) − 1 i = 0

Interpolation Examples Example: Legendre Gauss-Lobatto Quadrature ◮ Roots must be computed numerically ( 1 − x 2 ) L ′ N ( x ) = 0 ◮ Quadrature weights are 2 ω i = (16) N ( N + 1 )[ L N ( ξ )] 2 ◮ The formula is exact for polynomial of degree 2 N − 1: � 1 N � q ( x ) d x = q ( x i ) ω i (17) − 1 i = 0

Interpolation Some Math Proofs Some mathematical proofs for Gaussian Quadrature ◮ How do we know that the orthogonal polynomial p N + 1 has N + 1 roots? ◮ How can we prove that formula is exact for polynomials of degree 2 N + 1? ◮ Starting point is that p m is the set of orthogonal polynomials: � b p m ( x ) p n ( x ) w ( x ) d x = δ m , n � p m � 2 (18) a

Interpolation Some Math Proofs Multiplicity and location of roots ◮ Assume p N has only k roots in [ a , b ] : p N can be written as p N = ( x − r 1 )( x − r 2 ) . . . ( x − r k ) h ( x ) (19) where h k is a polynomial of degree N − k that must be single-signed in [ a , b ] . If it changes sign, then there must be another root r k + 1 ◮ Consider the following polynomial: z ( x ) = ( x − r 1 )( x − r 2 ) . . . ( x − r k ) p N ( x ) (20) [( x − r 1 )( x − r 2 ) . . . ( x − r k )] 2 h ( x ) z ( x ) = (21) ◮ z ( x ) must also be single-signed.

Interpolation Some Math Proofs Multiplicity and location of roots p N is orthogonal to all polynomials of degree k < N then � b ( x − r 1 )( x − r 2 ) . . . ( x − r k ) p N ( x ) w ( x ) d x = δ N , k a � �� � polynomial of degree k � b [( x − r 1 )( x − r 2 ) . . . ( x − r k )] 2 h ( x ) w ( x ) d x = δ N , k a Since z is single signed this is possible only if k = N , then p N has N roots in [ a , b ] Since p N is of degree N and has N roots, the roots must be isolated. Likewise p ′ N ( x ) must have N − 1 roots corresponding to the extrema of p N ( x ) .

Interpolation Some Math Proofs Exactness of quadrature for polynomials Let f ( x ) be a polynomial of degree 2 N + 1, it can always be written as f ( x ) = p N + 1 ( x ) q ( x ) + r ( x ) (22) where q ( x ) , and r ( x ) are polynomials of degree at most N . � b � b � b f ( x ) w ( x ) d x = p N + 1 ( x ) q ( x ) w ( x ) d x + r ( x ) w ( x ) d x a a a � �� � � �� � 0 by orthogonality Exact quadrature N N � � = p N + 1 ( x i ) q ( x i ) ω i + r ( x i ) ω i � �� � i = 0 i = 0 0 N N � � = [ p N + 1 ( x i ) q ( x i ) + r ( x i )] ω i = f ( x i ) ω i � �� � i = 0 i = 0 f ( x i )