[ Get this as an ipynb notebook here ]
Test 3 is May 23, 2016. It will cover the 8th, 9th and 10th projects on Newton's method; max and min problems; and integration.
To begin, you should "run all" so that these packages are loaded into your Julia session.
using Plots using Roots Base.ctranspose(f::Function) = D(f) ## make f' work in place of D(f)
ctranspose (generic function with 38 methods)
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)
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:
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.
In the integration project we saw how to compute the integral of a function many ways:
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)
Use these cells or add some more to answer the questions on the exam.
