Lesson 20: Newton-Cotes Rules restart; with(Student[Calculus1]):

I want to look at the errors in Newton-Cotes rules with different orders, all using the same n, for some functions on the interval 0 .. 1. I'll take n to be 36, so the order k can be any factor of 36. These are the possibilities. K := [1,2,3,4,6,9,12,18,36]; First I'll use our function LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2J0YrLUYjNiYtRiw2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y6USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRC8lKXN0cmV0Y2h5R0ZELyUqc3ltbWV0cmljR0ZELyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRkQvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZTLUkobWZlbmNlZEdGJDYkLUYsNiVRInhGJ0Y2RjlGQEZALUY9Ni1RKSZlcXVhbHM7RidGQEZCRkVGR0ZJRktGTUZPL0ZSUSwwLjI3Nzc3NzhlbUYnL0ZVRmpuLUkmbWZyYWNHRiQ2KC1JI21uR0YkNiRRIjFGJ0ZALUYjNiZGX28tRj02LVEnJnBsdXM7RidGQEZCRkVGR0ZJRktGTUZPL0ZSUSwwLjIyMjIyMjJlbUYnL0ZVRmlvRllGQC8lLmxpbmV0aGlja25lc3NHRmJvLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRl9wLyUpYmV2ZWxsZWRHRkRGQEYrRkA=. Digits:= 30: J:= int(1/(1+x),x=0..1): seq(evalf(J - ApproximateInt(1/(1+x),x=0..1,method=newtoncotes[K[j]], partition=36/K[j])), j=1..9); For this function, the higher-order rules turned out to be better. But if we take an f whose higher derivatives grow faster, that might not be true. f:= x -> 1/(x^2 + 1/400): J:= int(f(x),x=0..1): seq(evalf(J - ApproximateInt(f(x),x=0..1,method=newtoncotes[K[j]], partition=36/K[j])), j=1..9); Here the best answer was obtained with k = 1, i.e. the Trapezoid rule.

Here's an innocent-looking function where, taking partition = 1 (so n = k), the errors tend to get worse as n and k increase. f:= x -> 1/(x^2+1); J:= int(f(x),x=-5..5); seq(evalf(J - ApproximateInt(f(x),x=-5..5,method=newtoncotes[k], partition=1)), k=1..30); The Newton-Cotes rule has a close connection to interpolation by polynomials. Consider the Newton-Cotes rule of order k with n = k. The result depends on the value of your function f at the LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEia0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RIitGJy9GOlEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkQvJSlzdHJldGNoeUdGRC8lKnN5bW1ldHJpY0dGRC8lKGxhcmdlb3BHRkQvJS5tb3ZhYmxlbGltaXRzR0ZELyUnYWNjZW50R0ZELyUnbHNwYWNlR1EsMC4yMjIyMjIyZW1GJy8lJ3JzcGFjZUdGUy1JI21uR0YkNiRRIjFGJ0ZARkBGK0ZARitGQA== equally spaced points LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUkjbW5HRiQ2JFEiMEYnL0Y2USdub3JtYWxGJ0Y+LyUvc3Vic2NyaXB0c2hpZnRHRj1GPg== ,..., LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRImtGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPQ==. If f was a polynomial of degree at most k, this result would be correct. So what the rule gives you is the integral for a polynomial of degree at most k that agrees with f at those points. We say that polynomial interpolates the values of f at LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYpLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUkjbW5HRiQ2JFEiMEYnL0Y2USdub3JtYWxGJ0Y+LyUvc3Vic2NyaXB0c2hpZnRHRj0tSSNtb0dGJDYtUSIsRidGPi8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGNC8lKXN0cmV0Y2h5R0ZILyUqc3ltbWV0cmljR0ZILyUobGFyZ2VvcEdGSC8lLm1vdmFibGVsaW1pdHNHRkgvJSdhY2NlbnRHRkgvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR1EsMC4zMzMzMzMzZW1GJy1GQzYtUSMuLkYnRj5GRi9GSkZIRktGTUZPRlFGUy9GVlEsMC4yMjIyMjIyZW1GJy9GWUZXLUZDNi1RIi5GJ0Y+RkZGaG5GS0ZNRk9GUUZTRlVGW29GQi1GLDYlRi4tRiM2JC1GLzYlUSJrRidGMkY1Rj5GQEY+. If partition > 1, the Newton-Cotes rule does this on each partition of k intervals. ApproximateInt with the option output=plot shows you your function f and the interpolating polynomials on each interval. ApproximateInt(1/(x^2+1),x=-5..5,method=newtoncotes[4],partition=1,output=plot); Here's an animation of this for different k with partition=1. For technical reasons, the animate command doesn't work here. However, as we saw in Lesson 15 there's another way to produce an animation, using the display command: you give it a list of plots (one for each frame), and the option insequence=true. with(plots): display([seq(ApproximateInt(1/(x^2+1),x=-5..5,method=newtoncotes[k],partition=1, output=plot),k=1..20)],insequence=true); display([seq(ApproximateInt(1/(x+1),x=0..1,method=newtoncotes[k],partition=1, output=plot),k=1..20)],insequence=true);

Suppose LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= is some quantity you want to calculate, and you have available some approximations LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnRj1GPQ==. Often you know something about how well LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnRj1GPQ== approximates LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=, e.g. LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYtLUkjbWlHRiQ2I1EhRictSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuMGVtRicvJSZkZXB0aEdGNC8lKmxpbmVicmVha0dRKG5ld2xpbmVGJy1GMDYmRjJGNUY4L0Y7USVhdXRvRictRiw2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RIj1GJy9GSFEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRlIvJSlzdHJldGNoeUdGUi8lKnN5bW1ldHJpY0dGUi8lKGxhcmdlb3BHRlIvJS5tb3ZhYmxlbGltaXRzR0ZSLyUnYWNjZW50R0ZSLyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGW28tRiw2JVEiQUYnRkRGRy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0ZERkdGTkZOLUZLNi1RIitGJ0ZORlBGU0ZVRldGWUZlbkZnbi9Gam5RLDAuMjIyMjIyMmVtRicvRl1vRl1wLUYsNiVRIk9GJy9GRUZSRk4tRmJvNiQtRiM2JC1JJW1zdXBHRiQ2JUZmby1GIzYlLUZLNi1RKiZ1bWludXMwO0YnRk5GUEZTRlVGV0ZZRmVuRmduRlxwRl5wLUYsNiVRInBGJ0ZERkdGTi8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGTkZORk4=. That is, the error in approximating J by LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y6USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRC8lKXN0cmV0Y2h5R0ZELyUqc3ltbWV0cmljR0ZELyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRkQvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZTLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEibkYnRjZGOUZARkBGQEYrRkBGK0ZA is less than some constant times LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JC1JJW1zdXBHRiQ2JS1GLDYlUSJuRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUYjNiUtSSNtb0dGJDYtUSomdW1pbnVzMDtGJy9GO1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkcvJSlzdHJldGNoeUdGRy8lKnN5bW1ldHJpY0dGRy8lKGxhcmdlb3BHRkcvJS5tb3ZhYmxlbGltaXRzR0ZHLyUnYWNjZW50R0ZHLyUnbHNwYWNlR1EsMC4yMjIyMjIyZW1GJy8lJ3JzcGFjZUdGVi1GLDYlUSJwRidGN0Y6RkMvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRkNGK0ZD. But let's suppose you have more, say LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Li1GLDYlUSJKRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEiPUYnL0Y4USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQi8lKXN0cmV0Y2h5R0ZCLyUqc3ltbWV0cmljR0ZCLyUobGFyZ2VvcEdGQi8lLm1vdmFibGVsaW1pdHNHRkIvJSdhY2NlbnRHRkIvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZRLUYsNiVRIkFGJ0Y0RjctSShtZmVuY2VkR0YkNiQtRiM2JC1GLDYlUSJuRidGNEY3Rj5GPi1GOzYtUSIrRidGPkZARkNGRUZHRklGS0ZNL0ZQUSwwLjIyMjIyMjJlbUYnL0ZTRl1vLUYsNiVRIkNGJ0Y0RjctRjs2LVExJkludmlzaWJsZVRpbWVzO0YnRj5GQEZDRkVGR0ZJRktGTS9GUFEmMC4wZW1GJy9GU0Zmby1JJW1zdXBHRiQ2JUZmbi1GIzYlLUY7Ni1RKiZ1bWludXMwO0YnRj5GQEZDRkVGR0ZJRktGTUZcb0Zeby1GLDYlUSJwRidGNEY3Rj4vJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRmluLUYsNiVRIk9GJy9GNUZCRj4tRlg2JC1GIzYkLUZpbzYlRmZuLUYjNidGXXBGYHAtRjs2LVEoJm1pbnVzO0YnRj5GQEZDRkVGR0ZJRktGTUZcb0Zeby1GLDYlUS0mdmFyZXBzaWxvbjtGJ0ZpcEY+Rj5GY3BGPkY+Rj5GK0Y+ for some LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEtJnZhcmVwc2lsb247RicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJ0Yy > 0 (and typically LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEtJnZhcmVwc2lsb247RicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJ0Yy = 1 or 2). Thus the error is approximately some constant times LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYlLUkjbW9HRiQ2LVEqJnVtaW51czA7RicvRjZRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZCLyUpc3RyZXRjaHlHRkIvJSpzeW1tZXRyaWNHRkIvJShsYXJnZW9wR0ZCLyUubW92YWJsZWxpbWl0c0dGQi8lJ2FjY2VudEdGQi8lJ2xzcGFjZUdRLDAuMjIyMjIyMmVtRicvJSdyc3BhY2VHRlEtRi82JVEicEYnRjJGNUY+LyUxc3VwZXJzY3JpcHRzaGlmdEdRIjBGJ0Y+. Unfortunately you don't know the constant LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiQ0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=. If you did know it, you could make a better approximation by using LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYpLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnRj0tSSNtb0dGJDYtUSIrRidGPS8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRS8lKXN0cmV0Y2h5R0ZFLyUqc3ltbWV0cmljR0ZFLyUobGFyZ2VvcEdGRS8lLm1vdmFibGVsaW1pdHNHRkUvJSdhY2NlbnRHRkUvJSdsc3BhY2VHUSwwLjIyMjIyMjJlbUYnLyUncnNwYWNlR0ZULUYsNiVRIkNGJ0YvRjItRkA2LVExJkludmlzaWJsZVRpbWVzO0YnRj1GQ0ZGRkhGSkZMRk5GUC9GU1EmMC4wZW1GJy9GVkZobi1JJW1zdXBHRiQ2JUY6LUYjNiUtRkA2LVEqJnVtaW51czA7RidGPUZDRkZGSEZKRkxGTkZQRlJGVS1GLDYlUSJwRidGL0YyRj0vJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRj0= instead of LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnRj1GPQ==. Richardson extrapolation remedies this difficulty by looking at two different LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y6USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRC8lKXN0cmV0Y2h5R0ZELyUqc3ltbWV0cmljR0ZELyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRkQvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZTLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEibkYnRjZGOUZARkBGQEYrRkBGK0ZA. We'll get both a better approximation for LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= and some idea of the error in LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnRj1GPQ==. Suppose we calculate LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnRj1GPQ== and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JKG1mZW5jZWRHRiQ2JC1GIzYkLUkmbWZyYWNHRiQ2KC1GIzYkLUYsNiVRIm5GJ0YvRjIvRjNRJ25vcm1hbEYnLUYjNiQtSSNtbkdGJDYkUSIyRidGQkZCLyUubGluZXRoaWNrbmVzc0dRIjFGJy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0ZPLyUpYmV2ZWxsZWRHUSZmYWxzZUYnRkJGQkZC. J:= 'J': e1:= J = A(n) + C*n^(-p) + O(n^(-p-epsilon)); e2:= eval(e1,n=n/2); Think of these as two equations in the two unknowns J and C. S:= solve({e1,e2},{J, C}); If we neglect the O terms: SR := simplify(eval(S,O=0)); CR:= eval(C,SR); JR:= eval(J,SR); simplify(CR-eval(C,S)); So the difference between the Richardson value LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiQ0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRIlJGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPQ== and the true C is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJPRicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjcvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkobWZlbmNlZEdGJDYkLUYjNiQtSSVtc3VwR0YkNiUtRiw2JVEibkYnL0Y1USV0cnVlRicvRjhRJ2l0YWxpY0YnLUYjNiUtRjs2LVEqJnVtaW51czA7RidGN0Y+RkBGQkZERkZGSEZKL0ZNUSwwLjIyMjIyMjJlbUYnL0ZQRmBvLUYsNiVRLSZ2YXJlcHNpbG9uO0YnRjRGN0Y3LyUxc3VwZXJzY3JpcHRzaGlmdEdRIjBGJ0Y3RjdGN0YrRjc=. simplify(JR-eval(J,S)); The difference between the Richardson value LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRIlJGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPQ== and the true J is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJPRicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjcvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkobWZlbmNlZEdGJDYkLUYjNiQtSSVtc3VwR0YkNiUtRiw2JVEibkYnL0Y1USV0cnVlRicvRjhRJ2l0YWxpY0YnLUYjNictRjs2LVEqJnVtaW51czA7RidGN0Y+RkBGQkZERkZGSEZKL0ZNUSwwLjIyMjIyMjJlbUYnL0ZQRmBvLUYsNiVRInBGJ0ZmbkZobi1GOzYtUSgmbWludXM7RidGN0Y+RkBGQkZERkZGSEZKRl9vRmFvLUYsNiVRLSZ2YXJlcHNpbG9uO0YnRjRGN0Y3LyUxc3VwZXJzY3JpcHRzaGlmdEdRIjBGJ0Y3RjdGN0YrRjc=. When n is large, the main contribution to the error in LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y6USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRC8lKXN0cmV0Y2h5R0ZELyUqc3ltbWV0cmljR0ZELyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRkQvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZTLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEibkYnRjZGOUZARkBGQEYrRkBGK0ZA is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2JVEiQ0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMSZJbnZpc2libGVUaW1lcztGJy9GM1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRj0vJSlzdHJldGNoeUdGPS8lKnN5bW1ldHJpY0dGPS8lKGxhcmdlb3BHRj0vJS5tb3ZhYmxlbGltaXRzR0Y9LyUnYWNjZW50R0Y9LyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGTC1JJW1zdXBHRiQ2JS1GLDYlUSJuRidGL0YyLUYjNiUtRjY2LVEqJnVtaW51czA7RidGOUY7Rj5GQEZCRkRGRkZIL0ZLUSwwLjIyMjIyMjJlbUYnL0ZORmVuLUYsNiVRInBGJ0YvRjJGOS8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGOQ==. If we approximate that error as LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiQ0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRIlJGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRictSSNtb0dGJDYtUTEmSW52aXNpYmxlVGltZXM7RidGPS8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGSC8lKXN0cmV0Y2h5R0ZILyUqc3ltbWV0cmljR0ZILyUobGFyZ2VvcEdGSC8lLm1vdmFibGVsaW1pdHNHRkgvJSdhY2NlbnRHRkgvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZXLUklbXN1cEdGJDYlLUYvNiVRIm5GJ0YyRjUtRiM2JS1GQzYtUSomdW1pbnVzMDtGJ0Y9RkZGSUZLRk1GT0ZRRlMvRlZRLDAuMjIyMjIyMmVtRicvRllGYG8tRi82JVEicEYnRjJGNUY9LyUxc3VwZXJzY3JpcHRzaGlmdEdGQUY9, how far off are we (i.e. what is the error in our approximation of the error in our approximation)? simplify(eval(J - A(n) - CR*n^(-p),S)); This is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEiT0YnLyUnaXRhbGljR1EmZmFsc2VGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRictSSNtb0dGJDYtUTAmQXBwbHlGdW5jdGlvbjtGJ0Y5LyUmZmVuY2VHRjgvJSpzZXBhcmF0b3JHRjgvJSlzdHJldGNoeUdGOC8lKnN5bW1ldHJpY0dGOC8lKGxhcmdlb3BHRjgvJS5tb3ZhYmxlbGltaXRzR0Y4LyUnYWNjZW50R0Y4LyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGUC1JKG1mZW5jZWRHRiQ2JC1GIzYkLUklbXN1cEdGJDYlLUYsNiVRIm5GJy9GN1EldHJ1ZUYnL0Y6USdpdGFsaWNGJy1GIzYnLUY9Ni1RKiZ1bWludXMwO0YnRjlGQEZCRkRGRkZIRkpGTC9GT1EsMC4yMjIyMjIyZW1GJy9GUkZiby1GLDYlUSJwRidGaG5Gam4tRj02LVEoJm1pbnVzO0YnRjlGQEZCRkRGRkZIRkpGTEZhb0Zjby1GLDYlUScmIzk0OTtGJ0Y2RjlGOS8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGOUY5RjlGK0Y5RitGOQ==. As long as LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2JVEiQ0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RJSZuZTtGJy9GM1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRj0vJSlzdHJldGNoeUdGPS8lKnN5bW1ldHJpY0dGPS8lKGxhcmdlb3BHRj0vJS5tb3ZhYmxlbGltaXRzR0Y9LyUnYWNjZW50R0Y9LyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGTC1JI21uR0YkNiRRIjBGJ0Y5Rjk=, that's much smaller than the actual error when LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= is large. So this should be a good approximation for the error in LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEiQUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y6USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRC8lKXN0cmV0Y2h5R0ZELyUqc3ltbWV0cmljR0ZELyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRkQvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZTLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEibkYnRjZGOUZARkBGQEYrRkBGK0ZA. simplify(CR*n^(-p)); We don't know a good approximation for the error in our improved approximation LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRIlJGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPQ==, only that it is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJPRicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjcvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkobWZlbmNlZEdGJDYkLUYjNiQtSSVtc3VwR0YkNiUtRiw2JVEibkYnL0Y1USV0cnVlRicvRjhRJ2l0YWxpY0YnLUYjNictRjs2LVEqJnVtaW51czA7RidGN0Y+RkBGQkZERkZGSEZKL0ZNUSwwLjIyMjIyMjJlbUYnL0ZQRmBvLUYsNiVRInBGJ0ZmbkZobi1GOzYtUSgmbWludXM7RidGN0Y+RkBGQkZERkZGSEZKRl9vRmFvLUYsNiVRLSZ2YXJlcHNpbG9uO0YnRjRGN0Y3LyUxc3VwZXJzY3JpcHRzaGlmdEdRIjBGJ0Y3RjdGN0YrRjc=. But LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiQ0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRIlJGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRictSSNtb0dGJDYtUTEmSW52aXNpYmxlVGltZXM7RidGPS8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGSC8lKXN0cmV0Y2h5R0ZILyUqc3ltbWV0cmljR0ZILyUobGFyZ2VvcEdGSC8lLm1vdmFibGVsaW1pdHNHRkgvJSdhY2NlbnRHRkgvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZXLUklbXN1cEdGJDYlLUYvNiVRIm5GJ0YyRjUtRiM2JS1GQzYtUSomdW1pbnVzMDtGJ0Y9RkZGSUZLRk1GT0ZRRlMvRlZRLDAuMjIyMjIyMmVtRicvRllGYG8tRi82JVEicEYnRjJGNUY9LyUxc3VwZXJzY3JpcHRzaGlmdEdGQUY9 should be a fairly conservative estimate for it, at least if LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= is large.

It turns out that for any LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiTkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= the error in LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJURicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEwJkFwcGx5RnVuY3Rpb247RicvRjhRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZCLyUpc3RyZXRjaHlHRkIvJSpzeW1tZXRyaWNHRkIvJShsYXJnZW9wR0ZCLyUubW92YWJsZWxpbWl0c0dGQi8lJ2FjY2VudEdGQi8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRlEtSShtZmVuY2VkR0YkNiQtRiM2JC1GLDYlUSJuRidGNEY3Rj5GPkY+RitGPg== for smooth functions LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= can be expressed in the formLUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KEYrLUYjNictSSttdW5kZXJvdmVyR0YkNictSSNtb0dGJDYtUSYmU3VtO0YnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGPy8lKXN0cmV0Y2h5R1EldHJ1ZUYnLyUqc3ltbWV0cmljR0Y/LyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRj8vJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR1EsMC4xNjY2NjY3ZW1GJy1GIzYmLUYsNiVRImtGJy8lJ2l0YWxpY0dGRC9GO1EnaXRhbGljRictRjc2LVEiPUYnRjpGPUZAL0ZDRj9GRS9GSEY/L0ZKRj9GSy9GTlEsMC4yNzc3Nzc4ZW1GJy9GUUZdby1JI21uR0YkNiRRIjFGJ0Y6RjotRiw2JVEiTkYnRlhGWkZLLyUsYWNjZW50dW5kZXJHRj9GKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EkNS4wRicvJSZkZXB0aEdGXXAvJSpsaW5lYnJlYWtHUSVhdXRvRictSSZtZnJhY0dGJDYoLUYjNiQtSSVtc3ViR0YkNiUtRiw2JVEiY0YnRlhGWi1GIzYmRistRiM2Ji1GYG82JFEiMkYnRjotRjc2LVExJkludmlzaWJsZVRpbWVzO0YnRjpGPUZARmluRkVGam5GW29GS0ZNL0ZRRk9GVUY6RitGOi8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRjotRiM2JC1JJW1zdXBHRiQ2JS1GLDYlUSJuRidGWEZaRmNxLyUxc3VwZXJzY3JpcHRzaGlmdEdGXnJGOi8lLmxpbmV0aGlja25lc3NHRmJvLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRl1zLyUpYmV2ZWxsZWRHRj9GOi1GNzYtUSIrRidGOkY9RkBGaW5GRUZqbkZbb0ZLL0ZOUSwwLjIyMjIyMjJlbUYnL0ZRRmZzLUYjNiYtRiw2JVEiT0YnL0ZZRj9GOi1GNzYtUTAmQXBwbHlGdW5jdGlvbjtGJ0Y6Rj1GQEZpbkZFRmpuRltvRktGTUZbci1JKG1mZW5jZWRHRiQ2JC1GIzYkLUZncDYoLUYjNiRGX29GOi1GIzYkLUZicjYlRmRyLUYjNidGKy1GIzYmRmVxRmhxRmNvRjpGYnNGZXFGOkZnckY6RmlyRltzRl5zRmBzRjpGOkY6RitGOkYrRjo= for some constants LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiY0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYmLUYvNiNRIUYnLUYjNiYtSSNtbkdGJDYkUSIyRicvRjZRJ25vcm1hbEYnLUkjbW9HRiQ2LVExJkludmlzaWJsZVRpbWVzO0YnRkMvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRksvJSlzdHJldGNoeUdGSy8lKnN5bW1ldHJpY0dGSy8lKGxhcmdlb3BHRksvJS5tb3ZhYmxlbGltaXRzR0ZLLyUnYWNjZW50R0ZLLyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGWi1GLzYlUSJrRidGMkY1RkNGOkZDLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGQw== depending on LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYlLUkobWZlbmNlZEdGJDYkLUYjNiktSSNtbkdGJDYkUSIyRicvRjZRJ25vcm1hbEYnLUkjbW9HRiQ2LVEifkYnRkMvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRksvJSlzdHJldGNoeUdGSy8lKnN5bW1ldHJpY0dGSy8lKGxhcmdlb3BHRksvJS5tb3ZhYmxlbGltaXRzR0ZLLyUnYWNjZW50R0ZLLyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGWi1GLzYlUSJrRidGMkY1LUZGNi1RKCZtaW51cztGJ0ZDRklGTEZORlBGUkZURlYvRllRLDAuMjIyMjIyMmVtRicvRmZuRl5vLUZANiRRIjFGJ0ZDRjJGNUZDRjJGNS8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGQw==, LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiYUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiYkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= (but not on LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiTkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=). I want to (somewhat) justify this statement. Consider the integral of LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJmRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEwJkFwcGx5RnVuY3Rpb247RicvRjhRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZCLyUpc3RyZXRjaHlHRkIvJSpzeW1tZXRyaWNHRkIvJShsYXJnZW9wR0ZCLyUubW92YWJsZWxpbWl0c0dGQi8lJ2FjY2VudEdGQi8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRlEtSShtZmVuY2VkR0YkNiQtRiM2JC1GLDYlUSJ4RidGNEY3Rj5GPkY+RitGPg== on one of our subintervals LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkobWZlbmNlZEdGJDYmLUYjNiYtSSVtc3ViR0YkNiUtSSNtaUdGJDYlUSJ4RicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUYjNiYtRjQ2I1EhRictRiM2Ji1GNDYlUSJrRidGN0Y6LUkjbW9HRiQ2LVEoJm1pbnVzO0YnL0Y7USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGTy8lKXN0cmV0Y2h5R0ZPLyUqc3ltbWV0cmljR0ZPLyUobGFyZ2VvcEdGTy8lLm1vdmFibGVsaW1pdHNHRk8vJSdhY2NlbnRHRk8vJSdsc3BhY2VHUSwwLjIyMjIyMjJlbUYnLyUncnNwYWNlR0Zobi1JI21uR0YkNiRRIjFGJ0ZLRktGP0ZLLyUvc3Vic2NyaXB0c2hpZnRHUSIwRictRkg2LVEiLEYnRktGTS9GUUY5RlJGVEZWRlhGWi9GZ25RJjAuMGVtRicvRmpuUSwwLjMzMzMzMzNlbUYnLUYxNiVGMy1GIzYkRkRGS0Zfb0ZLRksvJSVvcGVuR1EiW0YnLyUmY2xvc2VHUSJdRidGSw==, which I'll write as LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkobWZlbmNlZEdGJDYmLUYjNigtSSNtaUdGJDYjUSFGJy1GIzYmLUYxNiVRIm1GJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictSSNtb0dGJDYtUSgmbWludXM7RicvRj1RJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZHLyUpc3RyZXRjaHlHRkcvJSpzeW1tZXRyaWNHRkcvJShsYXJnZW9wR0ZHLyUubW92YWJsZWxpbWl0c0dGRy8lJ2FjY2VudEdGRy8lJ2xzcGFjZUdRLDAuMjIyMjIyMmVtRicvJSdyc3BhY2VHRlYtRjE2JVEickYnRjlGPEZDLUZANi1RIixGJ0ZDRkUvRklGO0ZKRkxGTkZQRlIvRlVRJjAuMGVtRicvRlhRLDAuMzMzMzMzM2VtRictRiM2JkY2LUZANi1RIitGJ0ZDRkVGSEZKRkxGTkZQRlJGVEZXRllGQ0YwRkNGQy8lJW9wZW5HUSJbRicvJSZjbG9zZUdRIl1GJ0ZD, where LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJyRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEiPUYnL0Y4USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQi8lKXN0cmV0Y2h5R0ZCLyUqc3ltbWV0cmljR0ZCLyUobGFyZ2VvcEdGQi8lLm1vdmFibGVsaW1pdHNHRkIvJSdhY2NlbnRHRkIvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZRLUkmbWZyYWNHRiQ2KC1GIzYmRistRiM2Ji1GLDYlUSJiRidGNEY3LUY7Ni1RKCZtaW51cztGJ0Y+RkBGQ0ZFRkdGSUZLRk0vRlBRLDAuMjIyMjIyMmVtRicvRlNGXG8tRiw2JVEiYUYnRjRGN0Y+RitGPi1GIzYmRistRiM2Ji1JI21uR0YkNiRRIjJGJ0Y+LUY7Ni1RMSZJbnZpc2libGVUaW1lcztGJ0Y+RkBGQ0ZFRkdGSUZLRk0vRlBRJjAuMGVtRicvRlNGXXAtRiw2JVEibkYnRjRGN0Y+RitGPi8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGZ3AvJSliZXZlbGxlZEdGQkY+RitGPg== and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJtRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEpJmVxdWFscztGJy9GOFEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkIvJSlzdHJldGNoeUdGQi8lKnN5bW1ldHJpY0dGQi8lKGxhcmdlb3BHRkIvJS5tb3ZhYmxlbGltaXRzR0ZCLyUnYWNjZW50R0ZCLyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGUS1JJm1mcmFjR0YkNigtSShtZmVuY2VkR0YkNiQtRiM2Ji1JJW1zdWJHRiQ2JS1GLDYlUSJ4RidGNEY3LUYjNiYtRiw2JVEia0YnRjRGNy1GOzYtUSgmbWludXM7RidGPkZARkNGRUZHRklGS0ZNL0ZQUSwwLjIyMjIyMjJlbUYnL0ZTRmVvLUkjbW5HRiQ2JFEiMUYnRj5GPi8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnLUY7Ni1RJyZwbHVzO0YnRj5GQEZDRkVGR0ZJRktGTUZkb0Zmby1GZ242JUZpbkZeb0ZbcEY+Rj4tRmhvNiRRIjJGJ0Y+LyUubGluZXRoaWNrbmVzc0dGam8vJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGanAvJSliZXZlbGxlZEdGQkY+RitGPg==. f:= 'f': Errorterm := Int(f(x),x=m-r .. m+r) - r*(f(m-r)+f(m+r)); J:= Int(f(x),x=m-r .. m+r); The main tool in analyzing this is Integration by Parts. We've seen that the Student[Calculus1] package can do integration step-by-step. In particular, you can integrate by parts using Rule[parts, u, v](J) where u is the part to be differentiated and v is an antiderivative of the other part (where LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNictRiw2JVEiSkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RIj1GJy9GOlEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkQvJSlzdHJldGNoeUdGRC8lKnN5bW1ldHJpY0dGRC8lKGxhcmdlb3BHRkQvJS5tb3ZhYmxlbGltaXRzR0ZELyUnYWNjZW50R0ZELyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGUy1GIzYoLUY9Ni1RKCYjODc0NztGJ0ZARkJGRUZHRklGS0ZNRk8vRlJRJjAuMGVtRicvRlVGZm4tRiw2JVEidUYnRjZGOS1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0Zgby8lKmxpbmVicmVha0dRJWF1dG9GJy1GPTYtUTAmRGlmZmVyZW50aWFsRDtGJ0ZARkJGRUZHRklGS0ZNRk9GZW5GZ24tRiw2JVEidkYnRjZGOUZARitGQEYrRkBGK0ZA). In this case, I want LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNictRiw2JVEidUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RIj1GJy9GOlEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkQvJSlzdHJldGNoeUdGRC8lKnN5bW1ldHJpY0dGRC8lKGxhcmdlb3BHRkQvJS5tb3ZhYmxlbGltaXRzR0ZELyUnYWNjZW50R0ZELyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGUy1GIzYmLUYsNiVRImZGJ0Y2RjktRj02LVEwJkFwcGx5RnVuY3Rpb247RidGQEZCRkVGR0ZJRktGTUZPL0ZSUSYwLjBlbUYnL0ZVRmluLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEieEYnRjZGOUZARkBGQEYrRkBGK0ZARitGQA== and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEjZHZGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictSSNtb0dGJDYtUSI9RicvRjpRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZELyUpc3RyZXRjaHlHRkQvJSpzeW1tZXRyaWNHRkQvJShsYXJnZW9wR0ZELyUubW92YWJsZWxpbWl0c0dGRC8lJ2FjY2VudEdGRC8lJ2xzcGFjZUdRLDAuMjc3Nzc3OGVtRicvJSdyc3BhY2VHRlMtRiw2JVEjZHhGJ0Y2RjlGQEYrRkBGK0ZA. For v I'll use the form that best expresses the symmetry of our interval: LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RKCZtaW51cztGJy9GOlEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkQvJSlzdHJldGNoeUdGRC8lKnN5bW1ldHJpY0dGRC8lKGxhcmdlb3BHRkQvJS5tb3ZhYmxlbGltaXRzR0ZELyUnYWNjZW50R0ZELyUnbHNwYWNlR1EsMC4yMjIyMjIyZW1GJy8lJ3JzcGFjZUdGUy1GLDYlUSJtRidGNkY5RkBGK0ZARitGQA==. with(Student[Calculus1]): Rule[parts,f(x),x-m](J); The first two terms are this interval's contribution to the Trapezoid Rule. So the rest is the error term. Errorterm := op(3,rhs(%)); It's more convenient to use the D form of the derivative here. Errorterm := convert(%,D); Now in my next integration by parts, I'll differentiate LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNidGKy1GIzYmLUYsNiVRIkRGJy8lJ2l0YWxpY0dRJmZhbHNlRicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLUkjbW9HRiQ2LVEwJkFwcGx5RnVuY3Rpb247RidGOy8lJmZlbmNlR0Y6LyUqc2VwYXJhdG9yR0Y6LyUpc3RyZXRjaHlHRjovJSpzeW1tZXRyaWNHRjovJShsYXJnZW9wR0Y6LyUubW92YWJsZWxpbWl0c0dGOi8lJ2FjY2VudEdGOi8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRlItSShtZmVuY2VkR0YkNiQtRiM2JC1GLDYlUSJmRicvRjlRJXRydWVGJy9GPFEnaXRhbGljRidGO0Y7RjtGPi1GVjYkLUYjNiQtRiw2JVEieEYnRmduRmluRjtGO0Y7RitGO0YrRjs= and integrate LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RKCZtaW51cztGJy9GOlEnbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkQvJSlzdHJldGNoeUdGRC8lKnN5bW1ldHJpY0dGRC8lKGxhcmdlb3BHRkQvJS5tb3ZhYmxlbGltaXRzR0ZELyUnYWNjZW50R0ZELyUnbHNwYWNlR1EsMC4yMjIyMjIyZW1GJy8lJ3JzcGFjZUdGUy1GLDYlUSJtRidGNkY5RkBGK0ZARitGQA==. It's not so obvious what to use for the constant in this one, so for now I'll just put in a c. Rule[parts,D(f)(x),int(t-m,t=m..x)+c](Errorterm); As we'll see, the D(f) terms here will give us the LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYkLUkjbW5HRiQ2JFEiMUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJ0Y0LUYjNiYtSSNtaUdGJDYjUSFGJy1GIzYkLUklbXN1cEdGJDYlLUY6NiVRIm5GJy8lJ2l0YWxpY0dRJXRydWVGJy9GNVEnaXRhbGljRictRjE2JFEiMkYnRjQvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRjRGOUY0LyUubGluZXRoaWNrbmVzc0dGMy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0ZULyUpYmV2ZWxsZWRHUSZmYWxzZUYnRjQ= term in the expansion. For that to be the whole LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JC1JJm1mcmFjR0YkNigtRiM2JC1JI21uR0YkNiRRIjFGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRidGOi1GIzYmRistRiM2JC1JJW1zdXBHRiQ2JS1GLDYlUSJuRicvJSdpdGFsaWNHUSV0cnVlRicvRjtRJ2l0YWxpY0YnLUY3NiRRIjJGJ0Y6LyUxc3VwZXJzY3JpcHRzaGlmdEdRIjBGJ0Y6RitGOi8lLmxpbmV0aGlja25lc3NHRjkvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGVi8lKWJldmVsbGVkR1EmZmFsc2VGJ0Y6RitGOg== term, we'll need the integral on the right to be 0 when LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEkZicnRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnL0YzUSdub3JtYWxGJw== is constant (otherwise it would give us another LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYkLUkjbW5HRiQ2JFEiMUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJ0Y0LUYjNiYtSSNtaUdGJDYjUSFGJy1GIzYkLUklbXN1cEdGJDYlLUY6NiVRIm5GJy8lJ2l0YWxpY0dRJXRydWVGJy9GNVEnaXRhbGljRictRjE2JFEiMkYnRjQvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRjRGOUY0LyUubGluZXRoaWNrbmVzc0dGMy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0ZULyUpYmV2ZWxsZWRHUSZmYWxzZUYnRjQ= term). So I want int(1/2*(x^2+m^2-2*m*x+2*c),x=m-r .. m+r) = 0; That can be simplified. normal(%); solve(%,c); So this is what I want to take as v in this integration by parts. v2:= unapply(normal(int(t-m,t=m..x)-r^2/6), x); Rule[parts,D(f)(x),v2(x)](Errorterm); Errorterm:= simplify(%); The sum of the terms LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYmLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEickYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21uR0YkNiRRIjJGJy9GO1Enbm9ybWFsRicvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnLUkjbW9HRiQ2LVExJkludmlzaWJsZVRpbWVzO0YnRkEvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkwvJSlzdHJldGNoeUdGTC8lKnN5bW1ldHJpY0dGTC8lKGxhcmdlb3BHRkwvJS5tb3ZhYmxlbGltaXRzR0ZMLyUnYWNjZW50R0ZMLyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGZW4tSShtZmVuY2VkR0YkNiQtRiM2KC1GNDYjUSFGJy1GIzYmLUY0NiVRI2YnRidGN0Y6LUZHNi1RMCZBcHBseUZ1bmN0aW9uO0YnRkFGSkZNRk9GUUZTRlVGV0ZZRmZuLUZpbjYkLUYjNiZGXW8tRiM2Ji1GNDYlUSJtRidGN0Y6LUZHNi1RKCZtaW51cztGJ0ZBRkpGTUZPRlFGU0ZVRlcvRlpRLDAuMjIyMjIyMmVtRicvRmduRmVwRjNGQUZdb0ZBRkFGQUZhcC1GIzYmRmJvRmVvLUZpbjYkLUYjNiZGXW8tRiM2JkZecC1GRzYtUSIrRidGQUZKRk1GT0ZRRlNGVUZXRmRwRmZwRjNGQUZdb0ZBRkFGQUZdb0ZBRkFGQS1GIzYkLUY+NiRRIjNGJ0ZBRkEvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRlxyLyUpYmV2ZWxsZWRHRkxGQQ== over all subintervals is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYmLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEickYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21uR0YkNiRRIjJGJy9GO1Enbm9ybWFsRicvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnLUkjbW9HRiQ2LVExJkludmlzaWJsZVRpbWVzO0YnRkEvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkwvJSlzdHJldGNoeUdGTC8lKnN5bW1ldHJpY0dGTC8lKGxhcmdlb3BHRkwvJS5tb3ZhYmxlbGltaXRzR0ZMLyUnYWNjZW50R0ZMLyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGZW4tSShtZmVuY2VkR0YkNiQtRiM2KC1GNDYjUSFGJy1GIzYmLUY0NiVRI2YnRidGN0Y6LUZHNi1RMCZBcHBseUZ1bmN0aW9uO0YnRkFGSkZNRk9GUUZTRlVGV0ZZRmZuLUZpbjYkLUYjNiQtRjQ2JVEiYUYnRjdGOkZBRkFGQS1GRzYtUSgmbWludXM7RidGQUZKRk1GT0ZRRlNGVUZXL0ZaUSwwLjIyMjIyMjJlbUYnL0ZnbkZjcC1GIzYmRmJvRmVvLUZpbjYkLUYjNiQtRjQ2JVEiYkYnRjdGOkZBRkFGQUZdb0ZBRkFGQS1GIzYkLUY+NiRRIjNGJ0ZBRkEvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRmhxLyUpYmV2ZWxsZWRHRkxGQQ==. Since LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJyRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEiPUYnL0Y4USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQi8lKXN0cmV0Y2h5R0ZCLyUqc3ltbWV0cmljR0ZCLyUobGFyZ2VvcEdGQi8lLm1vdmFibGVsaW1pdHNHRkIvJSdhY2NlbnRHRkIvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZRLUkmbWZyYWNHRiQ2KC1GIzYmRistRiM2Ji1GLDYlUSJiRidGNEY3LUY7Ni1RKCZtaW51cztGJ0Y+RkBGQ0ZFRkdGSUZLRk0vRlBRLDAuMjIyMjIyMmVtRicvRlNGXG8tRiw2JVEiYUYnRjRGN0Y+RitGPi1GIzYmRistRiM2Ji1JI21uR0YkNiRRIjJGJ0Y+LUY7Ni1RMSZJbnZpc2libGVUaW1lcztGJ0Y+RkBGQ0ZFRkdGSUZLRk0vRlBRJjAuMGVtRicvRlNGXXAtRiw2JVEibkYnRjRGN0Y+RitGPi8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGZ3AvJSliZXZlbGxlZEdGQkY+RitGPg==, this gives us a term LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1JJm1mcmFjR0YkNigtRiM2JC1JJW1zdWJHRiQ2JS1GLDYlUSJjRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUYjNiQtSSNtbkdGJDYkUSIyRicvRkBRJ25vcm1hbEYnRkgvJS9zdWJzY3JpcHRzaGlmdEdRIjBGJ0ZILUYjNiZGKy1GIzYkLUklbXN1cEdGJDYlLUYsNiVRIm5GJ0Y8Rj9GRC8lMXN1cGVyc2NyaXB0c2hpZnRHRkxGSEYrRkgvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRmhuLyUpYmV2ZWxsZWRHUSZmYWxzZUYnLUkjbW9HRiQ2LVEiPUYnRkgvJSZmZW5jZUdGXW8vJSpzZXBhcmF0b3JHRl1vLyUpc3RyZXRjaHlHRl1vLyUqc3ltbWV0cmljR0Zdby8lKGxhcmdlb3BHRl1vLyUubW92YWJsZWxpbWl0c0dGXW8vJSdhY2NlbnRHRl1vLyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGYnAtRjI2KC1GIzYnLUkobWZlbmNlZEdGJDYkLUYjNihGKy1GIzYmLUYsNiVRI2YnRidGPEY/LUZfbzYtUTAmQXBwbHlGdW5jdGlvbjtGJ0ZIRmJvRmRvRmZvRmhvRmpvRlxwRl5wL0ZhcFEmMC4wZW1GJy9GZHBGZ3EtRmpwNiQtRiM2JC1GLDYlUSJhRidGPEY/RkhGSEZILUZfbzYtUSgmbWludXM7RidGSEZib0Zkb0Zmb0Zob0Zqb0ZccEZecC9GYXBRLDAuMjIyMjIyMmVtRicvRmRwRmRyLUYjNiZGYHFGY3EtRmpwNiQtRiM2JC1GLDYlUSJiRidGPEY/RkhGSEZIRitGSEZILUZfbzYtUTEmSW52aXNpYmxlVGltZXM7RidGSEZib0Zkb0Zmb0Zob0Zqb0ZccEZecEZmcUZocS1GIzYkLUZSNiUtRmpwNiQtRiM2JkZcc0ZgckZdckZIRkhGREZXRkhGK0ZILUYjNiZGKy1GIzYnLUZFNiRRIzEyRidGSEZfc0ZPRitGSEYrRkhGWUZmbkZpbkZbb0ZIRitGSA== in the Trapezoid Rule error. What's left? whatsleft := op(3,rhs(Errorterm)); The first factor of the integrand (divided by 6) is what we called v2(t) last time. v2:= unapply(normal(int(t-m,t=m..x)-r^2/6), x); We'll do another integration by parts. v3:= unapply(simplify(int(6*v2(t),t=m..x)),x); whatsleft := Rule[parts,(D@@2)(f)(x),v3(x)](whatsleft); whatsleft := simplify(rhs(whatsleft)); So far so good. We'll integrate by parts again. v4 := unapply(simplify(int(-v3(t),t=m..x)+c4),x); Again the constant c4 will be chosen so the integral of v4(x) over our interval is 0. normal(int(v4(x),x=m-r..m+r)); c4:= solve(%=0,c4); simplify(rhs(Rule[parts,(D@@3)(f)(x),v4(x)](whatsleft))); The sum of the terms LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYmLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEickYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21uR0YkNiRRIjRGJy9GO1Enbm9ybWFsRicvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnLUkjbW9HRiQ2LVExJkludmlzaWJsZVRpbWVzO0YnRkEvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkwvJSlzdHJldGNoeUdGTC8lKnN5bW1ldHJpY0dGTC8lKGxhcmdlb3BHRkwvJS5tb3ZhYmxlbGltaXRzR0ZMLyUnYWNjZW50R0ZMLyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdGZW4tSShtZmVuY2VkR0YkNiQtRiM2KC1GNDYjUSFGJy1GIzYmLUY0NiVRJWYnJydGJ0Y3RjotRkc2LVEwJkFwcGx5RnVuY3Rpb247RidGQUZKRk1GT0ZRRlNGVUZXRllGZm4tRmluNiQtRiM2JkZdby1GIzYmLUY0NiVRIm1GJ0Y3RjotRkc2LVEiK0YnRkFGSkZNRk9GUUZTRlVGVy9GWlEsMC4yMjIyMjIyZW1GJy9GZ25GZXBGM0ZBRl1vRkFGQUZBLUZHNi1RKCZtaW51cztGJ0ZBRkpGTUZPRlFGU0ZVRldGZHBGZnAtRiM2JkZib0Zlby1GaW42JC1GIzYmRl1vLUYjNiZGXnBGZ3BGM0ZBRl1vRkFGQUZBRl1vRkFGQUZBLUYjNiQtRj42JFEjNDVGJ0ZBRkEvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRlxyLyUpYmV2ZWxsZWRHRkxGQQ== over all subintervals is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1JJm1mcmFjR0YkNigtRiM2Ji1JJW1zdXBHRiQ2JS1GLDYlUSJyRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW5HRiQ2JFEiNEYnL0ZAUSdub3JtYWxGJy8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRictSSNtb0dGJDYtUTEmSW52aXNpYmxlVGltZXM7RidGRi8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGUS8lKXN0cmV0Y2h5R0ZRLyUqc3ltbWV0cmljR0ZRLyUobGFyZ2VvcEdGUS8lLm1vdmFibGVsaW1pdHNHRlEvJSdhY2NlbnRHRlEvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0Zqbi1JKG1mZW5jZWRHRiQ2JC1GIzYoRistRiM2Ji1GLDYlUSVmJycnRidGPEY/LUZMNi1RMCZBcHBseUZ1bmN0aW9uO0YnRkZGT0ZSRlRGVkZYRlpGZm5GaG5GW28tRl5vNiQtRiM2JC1GLDYlUSJiRidGPEY/RkZGRkZGLUZMNi1RKCZtaW51cztGJ0ZGRk9GUkZURlZGWEZaRmZuL0ZpblEsMC4yMjIyMjIyZW1GJy9GXG9GZXAtRiM2JkZkb0Znby1GXm82JC1GIzYkLUYsNiVRImFGJ0Y8Rj9GRkZGRkZGK0ZGRkZGRi1GIzYkLUZDNiRRIzQ1RidGRkZGLyUubGluZXRoaWNrbmVzc0dRIjFGJy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0ZqcS8lKWJldmVsbGVkR0ZRLUZMNi1RIj1GJ0ZGRk9GUkZURlZGWEZaRmZuL0ZpblEsMC4yNzc3Nzc4ZW1GJy9GXG9GY3ItRjI2KC1GIzYmLUY3NiUtRl5vNiQtRiM2JkZecEZhcEZdcUZGRkZGQkZIRktGXW9GRi1GIzYmRistRiM2Jy1GQzYkUSQ3MjBGJ0ZGRkstRiM2JC1GNzYlLUYsNiVRIm5GJ0Y8Rj9GQkZIRkZGK0ZGRitGRkZlcUZocUZbckZdckZGRitGRg==, which is our next term. And the pattern continues...

I want to apply Richardson extrapolation to the Trapezoid rule. h := n -> (b-a)/n: X := (k,n) -> a + k*h(n): a:= 0: b:= 1: T := n -> add((f(X(k-1,n)) + f(X(k,n)))/2 * h(n), k=1..n); J := int(f(x),x=a..b); JSFH The Trapezoid Rule LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNiYtRiw2JVEiVEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y6USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGRC8lKXN0cmV0Y2h5R0ZELyUqc3ltbWV0cmljR0ZELyUobGFyZ2VvcEdGRC8lLm1vdmFibGVsaW1pdHNHRkQvJSdhY2NlbnRHRkQvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZTLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEibkYnRjZGOUZARkBGQEYrRkBGK0ZA has error LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JkYrLUYjNictSSZtZnJhY0dGJDYoLUYjNiQtSSVtc3ViR0YkNiUtRiw2JVEiY0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUkjbW5HRiQ2JFEiMkYnL0ZCUSdub3JtYWxGJ0ZKLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGSi1GIzYmRistRiM2JC1JJW1zdXBHRiQ2JS1GLDYlUSJuRidGPkZBRkYvJTFzdXBlcnNjcmlwdHNoaWZ0R0ZORkpGK0ZKLyUubGluZXRoaWNrbmVzc0dRIjFGJy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0Zqbi8lKWJldmVsbGVkR1EmZmFsc2VGJy1JI21vR0YkNi1RIitGJ0ZKLyUmZmVuY2VHRl9vLyUqc2VwYXJhdG9yR0Zfby8lKXN0cmV0Y2h5R0Zfby8lKnN5bW1ldHJpY0dGX28vJShsYXJnZW9wR0Zfby8lLm1vdmFibGVsaW1pdHNHRl9vLyUnYWNjZW50R0Zfby8lJ2xzcGFjZUdRLDAuMjIyMjIyMmVtRicvJSdyc3BhY2VHRmRwLUYjNiYtRiw2JVEiT0YnL0Y/Rl9vRkotRmFvNi1RMCZBcHBseUZ1bmN0aW9uO0YnRkpGZG9GZm9GaG9Gam9GXHBGXnBGYHAvRmNwUSYwLjBlbUYnL0ZmcEZhcS1JKG1mZW5jZWRHRiQ2JC1GIzYkLUY0NigtRiM2JC1GRzYkRmduRkpGSi1GIzYmRistRiM2JC1GVDYlRlYtRkc2JFEiNEYnRkpGWUZKRitGSkZlbkZobkZbb0Zdb0ZKRkpGSkYrRkpGK0ZKRitGSg==, so the improved approximation using Richardson extrapolation would be TR[1] := n -> (2^2*T(n) - T(n/2))/(2^2-1); I'm calling it LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JC1JI21uR0YkNiRRIjFGJy9GNlEnbm9ybWFsRidGPi8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRj4= instead of just TR because, as we'll see, this will be the start of a sequence LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JC1GLzYlUSJrRidGMkY1L0Y2USdub3JtYWxGJy8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRj0= TR[1](2); That should look familiar. LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JC1JJW1zdWJHRiQ2JS1GLDYlUSNUUkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUkjbW5HRiQ2JFEiMUYnL0Y7USdub3JtYWxGJ0ZDLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGQ0YrRkM= is Simpson's rule. If T(n) has error LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYrLUkrbXVuZGVyb3ZlckdGJDYnLUkjbW9HRiQ2LVEmJlN1bTtGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRjcvJSlzdHJldGNoeUdRJXRydWVGJy8lKnN5bW1ldHJpY0dGNy8lKGxhcmdlb3BHRjwvJS5tb3ZhYmxlbGltaXRzR0Y8LyUnYWNjZW50R0Y3LyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdRLDAuMTY2NjY2N2VtRictRiM2Ji1JI21pR0YkNiVRImtGJy8lJ2l0YWxpY0dGPC9GM1EnaXRhbGljRictRi82LVEiPUYnRjJGNUY4L0Y7RjdGPS9GQEY3L0ZCRjdGQy9GRlEsMC4yNzc3Nzc4ZW1GJy9GSUZmbi1JI21uR0YkNiRRIjFGJ0YyRjItRk42JVEiTkYnRlFGU0ZDLyUsYWNjZW50dW5kZXJHRjctRk42I1EhRictSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJDUuMEYnLyUmZGVwdGhHRmlvLyUqbGluZWJyZWFrR1ElYXV0b0YnLUkmbWZyYWNHRiQ2KC1GIzYkLUklbXN1YkdGJDYlLUZONiVRImNGJ0ZRRlMtRiM2JkZhby1GIzYmLUZpbjYkUSIyRidGMi1GLzYtUTEmSW52aXNpYmxlVGltZXM7RidGMkY1RjhGWEY9RllGWkZDRkUvRklGR0ZNRjJGYW9GMi8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRjItRiM2JC1JJW1zdXBHRiQ2JS1GTjYlUSJuRidGUUZTRl9xLyUxc3VwZXJzY3JpcHRzaGlmdEdGanFGMi8lLmxpbmV0aGlja25lc3NHRltvLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRmlyLyUpYmV2ZWxsZWRHRjctRi82LVEiK0YnRjJGNUY4RlhGPUZZRlpGQy9GRlEsMC4yMjIyMjIyZW1GJy9GSUZicy1GTjYlUSJPRicvRlJGN0YyLUkobWZlbmNlZEdGJDYkLUYjNiQtRmNwNigtRiM2JEZobkYyLUYjNiQtRl5yNiVGYHItRiM2J0Zhby1GIzYmRmFxRmRxRlxvRjJGXnNGYXFGMkZjckYyRmVyRmdyRmpyRlxzRjJGMkZhb0Yy, what about LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JC1JJW1zdWJHRiQ2JS1GLDYlUSNUUkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUkjbW5HRiQ2JFEiMUYnL0Y7USdub3JtYWxGJ0ZDLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGQ0YrRkM=? It's not hard to see that this will have error LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYrLUkrbXVuZGVyb3ZlckdGJDYnLUkjbW9HRiQ2LVEmJlN1bTtGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRjcvJSlzdHJldGNoeUdRJXRydWVGJy8lKnN5bW1ldHJpY0dGNy8lKGxhcmdlb3BHRjwvJS5tb3ZhYmxlbGltaXRzR0Y8LyUnYWNjZW50R0Y3LyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdRLDAuMTY2NjY2N2VtRictRiM2Ji1JI21pR0YkNiVRImtGJy8lJ2l0YWxpY0dGPC9GM1EnaXRhbGljRictRi82LVEiPUYnRjJGNUY4L0Y7RjdGPS9GQEY3L0ZCRjdGQy9GRlEsMC4yNzc3Nzc4ZW1GJy9GSUZmbi1JI21uR0YkNiRRIjJGJ0YyRjItRk42JVEiTkYnRlFGU0ZDLyUsYWNjZW50dW5kZXJHRjctRk42I1EhRictSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJDUuMEYnLyUmZGVwdGhHRmlvLyUqbGluZWJyZWFrR1ElYXV0b0YnLUkmbWZyYWNHRiQ2KC1GIzYlLUZONiVRImNGJ0ZRRlMtSSVtc3ViR0YkNiUtRi82LVEiJ0YnRjJGNUY4RlhGPUZZRlpGQy9GRlEsMC4xMTExMTExZW1GJy9GSUZHLUYjNiZGYW8tRiM2JkZobi1GLzYtUTEmSW52aXNpYmxlVGltZXM7RidGMkY1RjhGWEY9RllGWkZDRkVGYnFGTUYyRmFvRjIvJS9zdWJzY3JpcHRzaGlmdEdRIjBGJ0YyLUYjNiQtSSVtc3VwR0YkNiUtRk42JVEibkYnRlFGU0ZlcS8lMXN1cGVyc2NyaXB0c2hpZnRHRlxyRjIvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRlxzLyUpYmV2ZWxsZWRHRjctRi82LVEiK0YnRjJGNUY4RlhGPUZZRlpGQy9GRlEsMC4yMjIyMjIyZW1GJy9GSUZlcy1GTjYlUSJPRicvRlJGN0YyLUkobWZlbmNlZEdGJDYkLUYjNiQtRmNwNigtRiM2JC1GaW42JEZpckYyRjItRiM2JC1GYHI2JUZici1GIzYnRmFvLUYjNiZGaG5GZ3FGXG9GMkZhc0ZobkYyRmVyRjJGZ3JGanJGXXNGX3NGMkYyLUYvNi1RIn5GJ0YyRjVGOEZYRj1GWUZaRkNGRUZicUYy(where LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYlLUkjbWlHRiQ2JVEiY0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JJW1zdWJHRiQ2JS1JI21vR0YkNi1RIidGJy9GM1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkAvJSlzdHJldGNoeUdGQC8lKnN5bW1ldHJpY0dGQC8lKGxhcmdlb3BHRkAvJS5tb3ZhYmxlbGltaXRzR0ZALyUnYWNjZW50R0ZALyUnbHNwYWNlR1EsMC4xMTExMTExZW1GJy8lJ3JzcGFjZUdRJjAuMGVtRictRiM2Ji1JI21uR0YkNiRRIjJGJ0Y8LUY5Ni1RIn5GJ0Y8Rj5GQUZDRkVGR0ZJRksvRk5GUkZQLUYsNiVRImtGJ0YvRjJGPC8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRjw= are other constants). So Richardson extrapolation improves LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JS1JI21uR0YkNiRRIjFGJy9GNlEnbm9ybWFsRidGMkY1LyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPg== to this: TR[2]:= n -> (2^4*TR[1](n)-TR[1](n/2))/(2^4-1); TR[2](4); This turns out to be the same as Newton-Cotes rule of order 4. with(Student[Calculus1]): ApproximateInt(f(x), x=a..b, partition=1, method=newtoncotes[4]); This should have error LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYrLUkrbXVuZGVyb3ZlckdGJDYnLUkjbW9HRiQ2LVEmJlN1bTtGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRicvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRjcvJSlzdHJldGNoeUdRJXRydWVGJy8lKnN5bW1ldHJpY0dGNy8lKGxhcmdlb3BHRjwvJS5tb3ZhYmxlbGltaXRzR0Y8LyUnYWNjZW50R0Y3LyUnbHNwYWNlR1EmMC4wZW1GJy8lJ3JzcGFjZUdRLDAuMTY2NjY2N2VtRictRiM2Ji1JI21pR0YkNiVRImtGJy8lJ2l0YWxpY0dGPC9GM1EnaXRhbGljRictRi82LVEiPUYnRjJGNUY4L0Y7RjdGPS9GQEY3L0ZCRjdGQy9GRlEsMC4yNzc3Nzc4ZW1GJy9GSUZmbi1JI21uR0YkNiRRIjNGJ0YyRjItRk42JVEiTkYnRlFGU0ZDLyUsYWNjZW50dW5kZXJHRjctRk42I1EhRictSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJDUuMEYnLyUmZGVwdGhHRmlvLyUqbGluZWJyZWFrR1ElYXV0b0YnLUkmbWZyYWNHRiQ2KC1GIzYmLUZONiVRImNGJ0ZRRlMtRi82LVEiJ0YnRjJGNUY4RlhGPUZZRlpGQy9GRlEsMC4xMTExMTExZW1GJy9GSUZHLUklbXN1YkdGJDYlRmpwLUYjNiZGYW8tRiM2Ji1GaW42JFEiMkYnRjItRi82LVExJkludmlzaWJsZVRpbWVzO0YnRjJGNUY4RlhGPUZZRlpGQ0ZFRl9xRk1GMkZhb0YyLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGMi1GIzYkLUklbXN1cEdGJDYlLUZONiVRIm5GJ0ZRRlNGZXEvJTFzdXBlcnNjcmlwdHNoaWZ0R0ZfckYyLyUubGluZXRoaWNrbmVzc0dRIjFGJy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0Zfcy8lKWJldmVsbGVkR0Y3LUYvNi1RIitGJ0YyRjVGOEZYRj1GWUZaRkMvRkZRLDAuMjIyMjIyMmVtRicvRklGaHMtRk42JVEiT0YnL0ZSRjdGMi1JKG1mZW5jZWRHRiQ2JC1GIzYkLUZjcDYoLUYjNiQtRmluNiRGXHNGMkYyLUYjNiQtRmNyNiVGZXItRiM2J0Zhby1GIzYmRmdxRmpxRlxvRjJGZHNGZ3FGMkZockYyRmpyRl1zRmBzRmJzRjJGMi1GLzYtUSIuRidGMkY1RjhGWEY9RllGWkZDRkVGX3FGMg== TR[3] := n -> (2^6*TR[2](n)-TR[2](n/2))/(2^6-1); TR[3](8); This one is not the same as the Newton-Cotes rule of order 8, although they evaluate f at the same points. ApproximateInt(f(x), x=a..b, partition=1, method=newtoncotes[8]); We could go farther with these "TR rules", but we won't. The correct name is "Romberg Integration".

LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JS1JI21uR0YkNiRRIjFGJy9GNlEnbm9ybWFsRidGMkY1LyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPg== was Simpson's Rule (the Newton-Cotes rule of order 2), and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JS1JI21uR0YkNiRRIjJGJy9GNlEnbm9ybWFsRidGMkY1LyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPg== was the Newton-Cotes rule of order 4, but LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JS1JI21uR0YkNiRRIjNGJy9GNlEnbm9ybWFsRidGMkY1LyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPg== is not a Newton-Cotes rule. Which is better, LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictSSNtbkdGJDYkUSIzRicvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPA== or the Newton-Cotes rule of order 8? On the one hand, LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JC1JI21uR0YkNiRRIjNGJy9GNlEnbm9ybWFsRidGPi8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRj4= should have error LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJPRicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjcvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkobWZlbmNlZEdGJDYkLUYjNiQtSSVtc3VwR0YkNiUtRiw2JVEibkYnL0Y1USV0cnVlRicvRjhRJ2l0YWxpY0YnLUYjNiUtRjs2LVEqJnVtaW51czA7RidGN0Y+RkBGQkZERkZGSEZKL0ZNUSwwLjIyMjIyMjJlbUYnL0ZQRmBvLUkjbW5HRiQ2JFEiOEYnRjdGNy8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGN0Y3RjdGK0Y3, while LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEsbmV3dG9uY290ZXNGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictSSNtbkdGJDYkUSI4RicvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPA== should have LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJPRicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjcvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkobWZlbmNlZEdGJDYkLUYjNiQtSSVtc3VwR0YkNiUtRiw2JVEibkYnL0Y1USV0cnVlRicvRjhRJ2l0YWxpY0YnLUYjNiUtRjs2LVEqJnVtaW51czA7RidGN0Y+RkBGQkZERkZGSEZKL0ZNUSwwLjIyMjIyMjJlbUYnL0ZQRmBvLUkjbW5HRiQ2JFEjMTBGJ0Y3RjcvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRjdGN0Y3RitGNw==. So LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEsbmV3dG9uY290ZXNGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictSSNtbkdGJDYkUSI4RicvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPA== should be better for large LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=. On the other hand, if LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic= is fairly small LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEjVFJGJy8lJ2l0YWxpY0dRJXRydWVGJy8lLG1hdGh2YXJpYW50R1EnaXRhbGljRictRiM2JC1JI21uR0YkNiRRIjNGJy9GNlEnbm9ybWFsRidGPi8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRj4= might be as good or better. Here is our function from last time that was bad for the Newton-Cotes rules with partition=1. f := x -> 1/((8*x-4)^2+1); evalf(J-TR[3](8)); evalf(J-ApproximateInt(f(x),x=0..1,method=newtoncotes[8],partition=1)); So in this particular case LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2JC1JJW1zdWJHRiQ2JS1GLDYlUSNUUkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUkjbW5HRiQ2JFEiM0YnL0Y7USdub3JtYWxGJ0ZDLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGQ0YrRkM= is much better. If we used a larger n, Newton-Cotes might win. seq([evalf(J-TR[3](8*k)), evalf(J-ApproximateInt(f(x),x=0..1,method=newtoncotes[8],partition=k))], k=1..10);

ApproximateInt(..., output=plot) in Student[Calculus1] package Rule[parts,...] in Student[Calculus1] package op

LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYjLUkjbWlHRiQ2I1EhRic=