HP Forums

Full Version: (SR-52) In defense of linear quadrature rules
You're currently viewing a stripped down version of our content. View the full version with proper formatting.
An extract from In Defense of Linear Quadrature Rules, William Squire (Aerospace Engineering Dept., West Virginia University), Comp. & Maths with Apple, Vol. 7, pp. 147.-149, Pergamon Press Ltd., 1981

Abstract--lt is shown that appropiate linear quadrature rules can handle integrands with singularities at or near the end points more effectively than the nonlinear methods proposed by Werner and Wuytack. A special 10 point Gauss rule gives good results. A method with exponential convergence gives high accuracy with a moderate number of nodes. Both procedures were implemented on a programmable hand calculator.
The purpose of this note is to demonstrate that:
(1) a special 10 point Gauss rule for integrands with singularities at or near the endpoints proposed by Harris and Evans [2] will give results comparing favourably to any other procedure using a comparable number of nodes.
(2) A quadrature rule, which Stenger [3] has shown to have exponential convergence, gives very accurate results for such integrands with a moderate number of function evaluations.
Both these procedures were implemented on an SR 52 programmable hand calculator.

The SR-52 implementation is given in Appendix A …

The method was implemented on an SR-52 as described in Appendix B …

Harris-Evans 10 point rule …

SR-52 program for Stenger quadrature (equation 1) …

Thank you, SlideRule.
I tried this modified point / weight quadrature on \(\int _0^6 e^{x^3}\; dx ≈ 5.963938092 × 10^{91}\)

gp/gw from Gaussian Quadrature Weights and Abscissae (n=10, rounded to 10 digits)
mp/mw is modified points/weights from In Defense of Linear Quadrature Rules, table 1

lua> gp = {0.1488743390, 0.4333953941, 0.6794095683, 0.8650633667, 0.9739065285}
lua> gw = {0.2955242247, 0.2692667193, 0.2190863625, 0.1494513492, 0.06667134431}
lua> mp = {0.2295037173, 0.6364758401, 0.9015072053, 0.9928383122, 0.9999843443}
lua> mw = {0.4501100825, 0.3483026852, 0.1744679776, 0.02696299772, 0.0001562579734}

lua> function integ(f,a,b,p,w) -- 10 points quadrature
:         local t, m, c = 0, (b-a)/2, (b+a)/2
:         for i=1,5 do t = t + (f(-m*p[i]+c) + f(m*p[i]+c)) * w[i] end
:         return m * t
:     end

lua> function f(x) return math.exp(x^3) end
lua> integ(f, 0, 6, gp, gw)     → 3.052910317e+089
lua> integ(f, 0, 6, mp, mw)   → 5.444304730e+091

Modified point/wieght looks much closer.
However, if we integrate only the dominant part, plain guassian points is better.

lua> integ(f, 5.50, 6, gp, gw)     → 5.942395811e+091
lua> integ(f, 5.50, 6, mp, mw)   → 5.577982336e+091

lua> integ(f, 5.75, 6, gp, gw)     → 5.963893713e+091
lua> integ(f, 5.75, 6, mp, mw)   → 5.907966288e+091
Reference URL's