4.11 Gaussian Quadrature
\[ \int_a^b f(x)dx \approx \sum_{i=1}^{n} w_i f(x_i) \] \[ \int_{-1}^1 f(x)dx \approx \sum w_i f(x_i) \]
| N | \(\pm x_i\) | \(w_i\) |
|---|---|---|
| 2 | 0.5773502692 | 1.0000000000 |
| 3 | 0.0000000000 0.7745966692 | 0.8888888889 0.5555555556 |
| 4 | 0.3399810436 0.8611363116 | 0.6521451549 0.3478548451 |
| 5 | 0.0000000000 0.5384693101 0.9061798459 | 0.5688888889 0.4786286705 0.2369268850 |
| 8 | 0.1834346425 0.5255324099 0.7966664774 0.9602898565 | 0.3626837834 0.3137066459 0.2223810345 0.1012285363 |
| 10 | 0.1488743389 0.4333953941 0.6794095682 0.8650633666 0.9739065285 | 0.2955242247 0.2692667193 0.2190863625 0.1494513491 0.0666713443 |
| 12 | 0.1252334085 0.3678314990 0.5873179543 0.7699026742 0.9041172564 0.9815606342 | 0.2491470458 0.2334925365 0.2031674267 0.1600783285 0.1069393260 0.0471753364 |
For general case, we need to transform varaible
\[ t = \frac{b-a}{2} x_i + \frac{a + b}{2} \]
\[ \int_a^b f(t)dt \approx \frac{b - a}{2} \sum w_i f \left( \frac{b - a}{2} x_i + \frac{a + b}{2} \right) \]
Practice script
GQuad8## function (fx, a, b)
## {
## xi = c(0.1834346425, 0.5255324099, 0.7966664774, 0.9602898565)
## wi = c(0.3626837834, 0.3137066459, 0.2223810345, 0.1012285363)
## xa = (b - a)/2
## xb = (a + b)/2
## Ar = xa * sum(wi * (fx(xi * xa + xb) + fx(-xi * xa + xb)))
## return(Ar)
## }
## <bytecode: 0x000001e2a72690c8>
## <environment: namespace:math>fx = function(x) exp(-x/100)
format(GQuad8(fx, 0, 24), digits=22)## [1] "21.33721389547881130966"format(integrate(fx, 0, 24), digits=22)## value abs.error
## "21.33721389334466067567" "2.368906614574302115298e-13"
## subdivisions message
## "1" "OK"
## call
## "integrate(f = fx, lower = 0, upper = 24)"fx2 = function(x) 3 * x * sin(x) - log(x)
format(GQuad8(fx2, 0.8, 2.6), digits=22)## [1] "6.887419623706937166219"format(integrate(fx2, 0.8, 2.6), digits=22)## value
## "6.887419623354108288993"
## abs.error
## "7.646571846102854117574e-14"
## subdivisions
## "1"
## message
## "OK"
## call
## "integrate(f = fx2, lower = 0.8, upper = 2.6)"