Monthly Archives: March 2017

Simpson’s Rule (Riemann Sums with Parabolas)

`int_a^b f(x) dx ~~ (h/3)((y_0+4y_1+y_2)+(y_2+4y_3+y_4)+…+(y_(n-2)+4y_(n-1)+y_n))`

`int_a^b f(x) dx ~~ (h/3)(y_0+4y_1+2y_2+4y_3+2y_4+…+2y_(n-2)+4y_(n-1)+y_n)`

`int_a^b f(x) dx ~~ sum_(j=1)^(n/2) (h/3)(f(x_(2j-2))+4f(x_(2j-1))+f(x_(2j)))`

where `h=(b-a)/n` and the error is `|E_S|<=(M(b-a)^5)/(180*n^4)` and `M` is an upper bound for `|f^((4))|` on `[a,b]`.

Ln with parabola

Example 1:

Evaluate `int_5^7 5x dx` using Simpson’s rule with `n=4`

`x_0` through `x_4` are `{5, 5.5, 6, 6.5, 7}` and `h=2/4=1/2`

Area `~~ 1/6*(5(5)+4*5(5.5)+2*5(6)+4*5(6.5)+5(7))=60`

Since `f'(x)=5` and `f”(x)=0`, then `f^((4))(x)=0` and the error is zero here.

Example 2:

Consider `int_-1^1 (x^2+5) dx` with `n=4`

`h=2/4=1/2` so `x_k=-1+k/2`, which is the set `{-1,-1/2,0,1/2,1}`

Area `~~1/6(((-1)^2+5)+4((-1/2)^2+5)+2((0)^2+5)+4((1/2)^2+5)+((1)^2+5))`

Area`=10 2/3`

Since `f^((4))(x)=0` the error here is also zero.

Example 3:

`int_0^2 (5t^3+8t) dt` with `n=4`

`h=2/4=1/2` and `x_k=0+k/2`, so the x-values are `{0,1/2,1,1 1/2,2}`

Area `~~ 1/6((5(0)^3+8(0))+4(5(1/2)^3+8(1/2))+2(5(1)^3+8(1))+4(5(1 1/2)^3+8(1 1/2))+(5(2)^3+8(2)))`

Area `=36`

Since `f^((4))(x)=0`, the error is zero here.

Example 4:

`int_3^9 4/s^2 ds` with `n=4`

`h=6/4=3/2` and `x_k=3+3k/2`, so the x-values are `{3, 4.5, 6, 7.5, 9}`

Area `~~1/2*(4/3^2+4*4/4.5^2+2*4/6^2+4*4/7.5^2+4/9^2)`

Example 5:

How many intervals are needed to estimate `int_1^2 1/x dx` to within 0.0001 square units of the actual area?

Recall the error is `|E_S|<=(M(b-a)^5)/(180*n^4)` and `M` is an upper bound for `|f^((4))|` on `[a,b]`.

`f(x)=x^(-1)`

`f'(x)=-1*x^(-2)`

`f”(x)=2*x^(-3)`

`f”(x)=-6*x^(-4)`

`f^((4))(x)=24*x^(-5)`

`f^((5))(x)=-144*x^(-6)`

I found the 5th derivative to show that `f^((4))(x)` is decreasing, since its derivative is negative on the interval `[1,2]`. So, `f^((4))(1)=24` is the maximum needed for the error estimation.

Plugging into `|E_S|<=(M(b-a)^5)/(180*n^4)`:

`0.0001<=(24*(1)^5)/(180*n^4)`

Solving for `n` we get `n <=((24)/(180*0.0001))^(1/4)~~6.04`

`A_(S_6) = 1/3*1/6*(1/1+4*1/(1 1/6)+2*1/(1 2/6)+4*1/(1 3/6)+2*1/(1 4/6)+4*1/(1 5/6)+1/2)~~0.693169`

`A_(S_4) = 1/3*1/4*(1/1+4*1/(1 1/4)+2*1/(1 2/4)+4*1/(1 3/4)+1/2)~~0.693254`

`A_(S_2) = 1/3*1/2*(1+4*1/1.5+1/2)~~.694444`

`ln(2)~~0.69314718`

The actual error from the 6th Simpson’s sum is about 0.00002.

`A_(S_20) = sum_(n=1)^10 1/(3*20)*(1/(1 + 1/20 (2 n – 2)) + 4/(1 + 1/20 (2 n – 1)) + 1/(1 + (2 n)/20)) = 5555158368718531/8014397185594800~~0.693147375`

For comparison, consider 100 terms with right endpoints: `sum_(k=1)^100 1/(100 (k/100 + 1))`
`~~0.690653430481824215252268721472608478928428119821791431288…`

Or the Trapezoidal Rule with 50 intervals:
`A_(T_50)=sum_(k=1)^50 1/(2*50)*(1/(1+(K-1)/50)+1/(1+k/50))~~0.693172179`

Example 6:

Here is a picture of a parabola intersecting the reciprocal function at `x=1,2,3`

Ln with parabola

The parabola, `y=1/6*x^2-1*x+11/6`, is shown in red.
The beautiful thing about Simpson’s rule is that we don’t actually need to calculate the coefficients of any of these parabolas. The weighted sum does all of the work for us.

Relatively Prime Linear Combinations

If there exist x and y for which ax+by=GCD(a,b) then GCD(x,y)=1.


Let g=GCD(a,b) and x and y solutions to ax+by=g.

Since g|a and g|b, there exist integers n and m such that a=gn and b=gm.

ax+by=gnx+gmy=g

So, nx+my=1

The GCD(n,m)=1 since g is the GCD(a,b).

Let d be a divisor of x. If d also divides my, then d must divide 1, which isn’t possible. So, GCD(x,y)=1.