Review for test 3 in MTH 229

[ Get this as an ipynb notebook here ]

Test 3 is 12/19/2016.

This review assumes you have done:

using MTH229

Newton's method

Newton's method is an iterative method to solve an equation of the form $f(x) = 0$. A guess, $x_i$ is updated by the formula $x_{i+1} = x_i - f(x_i)/f'(x_i)$ (which comes from finding the zero of the tangent line to the graph of $f$ at $((x_i), f(x_i))$). In the projects we saw we can write one step as:

x = 2
f(x) = x^3 - 3x + 1
x = x - f(x)/f'(x)

1.6666666666666667

Newton's method basically approximates the answer by the intersection point of the tangent line:

The guess of 2 is improved in one step to 1.6666666666666667. Doing another step, this guess is is improved to 1.548611111111111. We can visualize that doing more steps will further improve the guess, but the improvements will just be refinements on what we can reasonably read from a graph.

We can run many such steps by copy and paste (avoiding @take5):

x = x - f(x)/f'(x)
x = x - f(x)/f'(x)
x = x - f(x)/f'(x)
x = x - f(x)/f'(x)

1.532088886237968

That means starting from 2 we have run 5 iterations. The value of x and the value of f(x) are now:

x, f(x)

(1.532088886237968,4.884981308350689e-14)

Which we see has f(x) at the level of $10^{-14}$, not quite 0. Another step gets us as close as we can get:

x = x - f(x)/f'(x)
x, f(x)

(1.532088886237956,-4.440892098500626e-16)

In the project we saw also the newton function that makes implementing Newton's method easy. As well, the related fzero function does the same thing – find a zero of $f$ from an initial guess – but is a bit more robust.

newton(f, 2), fzero(f, 2)

(1.532088886237956,1.532088886237956)

You might be interested in verbose output, which is found through newton(f, 2, verbose=true).

Extrema

We consider problems like:

What inscribed rectangle in a half circle of radius 1 has maximum area? This figure shows 3 such inscribed rectangles.

In a math class, we might have these steps, where $(x,y)$ is the point where the rectangle intersects the half-circle in the first quadrant:

$A = (2x) \cdot y$ – the area of a rectangle
$r^2 = x^2 + y^2$ the point is on the unit circle (a constraint)
$y = (r^2 - x^2)^{1/2}$ solve for $y$
$A = (2x) \cdot (1 - x^2)^{1/2}$ substitute

Then we find the critical points of $A(x)$ and compare to the endpoints over $[0,r]$.

In julia we have with r=1:

y(x) = sqrt(1 - x^2)
A(x) = 2x * y(x)
plot(A, 0, 1)         ## see that max occurs at a critical point and then find that
x0 = fzero(D(A), 0.6) ## find the critical point -- and assign it a name (x0)
x0, y(x0), A(x0)      ## The three main values to answer a question

(0.7071067811865476,0.7071067811865475,1.0)

Setting up the problem – finding $A$ and the constraint – is what makes these problems difficult. Between the project and the notes there are several examples.

Integration

In the integration project we saw how to compute the integral of a function many ways:

Using Riemann integration. With this function

function riemann(f::Function, a::Real, b::Real, n::Int; method="right")
  if method == "right"
     meth(f,l,r) = f(r) * (r-l)
  elseif method == "left"
     meth(f,l,r) = f(l) * (r-l)
  elseif method == "trapezoid"
     meth(f,l,r) = (1/2) * (f(l) + f(r)) * (r-l)
  elseif method == "simpsons"
     meth(f,l,r) = (1/6) * (f(l) + 4*(f((l+r)/2)) + f(r)) * (r-l)
  end

  xs = a + (0:n) * (b-a)/n
  as = [meth(f, l, r) for (l,r) in zip(xs[1:end-1], xs[2:end])]
  sum(as)
  end

riemann (generic function with 1 method)

We could integrate most functions using Riemann sums:

riemann(sin, 0, pi, 10000)                     ## 10,000 riemann sums using right endpoints
riemann(sin, 0, pi, 10000, method="left")      ## 10,000 riemann sums using left endpoints
riemann(sin, 0, pi, 100,   method="simpsons")  ## 100 sums using Simpson's method

2.0000000006764704

All are approximations to the actual answer of 2.

We saw that julia ships with a powerful built in integration function quadgk:

a, err = quadgk(sin, 0, pi)

(2.0,1.7905676941154525e-12)

We don't specify the size of $n$, as quadgk adaptively finds its answer. The function also gives an estimate on the error, which is not the case for our simple-minded riemann function.

The one application of the integral in the assignment is the use of the integral to represent the volume of a figure of revolution: $V = \int_0^b \pi r(h)^2 dh$. For example, this red solo cup has radius as a function of $h$ given approximately by

r(h) = 2.85 + 0.15*h  ## 0 <= h <= 10.8

r (generic function with 1 method)

This will give a total volume when filled to the top of:

v(h) = pi * r(h)^2
a, err = quadgk(v, 0, 10.8)  ## 461 ml is about 15.58 ozs

(461.9223165286893,0.0)

To solve for what height gives 355 ml (about 12oz) we can guess:

quadgk(v, 0, 8) ## too small
quadgk(v, 0, 8.5) ## still too small
quadgk(v, 0, 9)  ## close enough (fzero(b -> quadgk(v, 0, b)[1]-355, 9) is exact)

(355.6204344047306,0.0)