0/11
Computational Electromagnetics :
Numerical Integration
Uday Khankhoje
Electrical Engineering, IIT Madras
1/11
Topics in this module
1 Simple Numerical Integration
2 Interpolating a Function
3 Advanced Numerical Integration
1/11
Table of Contents
1 Simple Numerical Integration
2 Interpolating a Function
3 Advanced Numerical Integration
2/11
Simple numerical integration
Rb
Take a function f (x) that has no analytical integral, want: a f (x) dx
Rules: Rectangle Trapezoidal Simpson’s
I≈ hf ( a+b
2 )
h
2 [f (a)+ f (b)] h
+ 4f ( a+b
3 [f (a) 2 ) + f (b)]
3/11
Simple Numerical Integration (contd.)
What was the basis of these simple rules? Taylor’s theorem:
(x−x0 )2 00
f (x) = f (x0 ) + (x − x0 )f 0 (x0 ) + 2 f (x − x0 ) + · · ·
Rh
Now, 0 f (x) dx =
Rectangle
Trapezoidal
[called Newton-Cotes formulas]
4/11
Extended rules – idea of ‘quadrature rule’
Rh R xN
Extend trapezoidal rule 0 f (x) dx ≈ h2 [f (0) + f (h)] for x1 f (x) dx
Quadrature rule – a formula in terms of:
∗ weights,
∗ function evaluation points
4/11
Table of Contents
1 Simple Numerical Integration
2 Interpolating a Function
3 Advanced Numerical Integration
5/11
Interpolating a function
Many functions we want to integrate can’t be done analytically.
• What can be done? Weierstrass Approximation Theorem:
If f is a continuous real-valued function on [a, b] and if any > 0 is given, then
there exists a polynomial p on [a, b] such that |f (x) − p(x)| < for all x ∈ [a, b].
• Fn known at few points – interpolate a polynomial fn: 2, 3, . . . , N
Lagrange polynomials
6/11
Integrating this interpolated function
R xN R xN
x1 f (x) dx ≈ x1 fN −1 (x) dx
Another kind of quadrature rule
Accurate to polynomial order N − 1, needing N points
6/11
Table of Contents
1 Simple Numerical Integration
2 Interpolating a Function
3 Advanced Numerical Integration
7/11
Can we do better? Gaussian quadrature
Summary: Clever math gives poly accuracy of order 2N − 1 using N points
0. We know the inner product of vectors . . . but functions?
Rb
1. Construct a polynomial pN (x) (order N ) s.t. a xk pN (x) dx = 0, k ∈ [0, N − 1]
Gram-Schmidt → Legendre Polynomials
8/11
Gaussian quadrature (contd.)
2. Roots of pN (x)?
3. Approximate f (x) to poly order 2N − 1 using pN (x) (Euclidean division)
f2N −1 (x)
pN (x) =
3x3 −2x2 +4x−3 28x+30
e.g. x2 +3x+3
= (3x − 11) + x2 +3x+3
=⇒ f2N −1 (x) = q(x)PN (x) + r(x)
9/11
Gaussian quadrature (contd.)
4. Now integrate on both sides of: f2N −1 (x) = q(x)PN (x) + r(x)
Lagrange polynomials for r(x)
5. If the N points are chosen as roots of pN (x)?
Called Gauss-Lengedre quadrature rule,
accurate to order 2N − 1.
10/10
Gaussian quadrature (contd.)
Take an example, f (x) = (x + 1)3 over [−1, 1]
p0 (x) = 1, p1 (x) = x, p2 (x) = (3x2 − 1)/2
2-pt Gauss-Lengendre quadrature
Exact calculation
3-pt trapezoidal rule
10/11
Topics that were covered in this module
1 Simple Numerical Integration
2 Interpolating a Function
3 Advanced Numerical Integration
Reference: Chapter 4 of Numerical recipes in C++ - Brian P. Flannery, Saul Teukolsky,
William H. Press, and William T. Vetterling