A modification to the regula falsi method uses two starting estimates that must bracket the root. The algorithm is as follows:
1.a.Given two starting estimates, x1 and x2 that bracket the root, compute the mid point, xm = (x1 + x2)/2.
b.Compute the corresponding function evaluations for these three points: y1 = f(x1), y2 = f(x2) and ym = f(xm).
c. Compute s = If s = 0, then set x = xm and stop. Otherwise continue.
d. Update the estimate
e. Tighten the bracket by keeping xn and one of x1, x2 or xm with:
If xm and xn bracket the root, set x1 = xm and x2 = xn
Else, if x1 and xn bracket the root, set x2 = xn
Else, it must be the case that xn and x2 bracket the root, so set x1 = xn.
f. Continue by returning to step 1.
Write a Matlab root finding routine using the following algorithm. Test your algorithm on a variety of nonlinear functions. How does it compare to other algorithms?
2. A simple model of the spread of disease gives
where P(t) represents the number of individuals in the population who are infected, and C is the constant size of the total population.
(a) Show that the solution of this differential equation is the logistic function.
(b) Suppose that the parameter k fluctuates (perhaps because the population is more susceptible during certain seasons). Solve the modified problem
with P(0) = 10, k = 2, C = 200.
3. Solve the Bessel equation of order one,
starting with y(1) = 1 and y 0 (1) = 0 on 1 ≤ x ≤ 4
4. The motion of one body around another, such as a comet orbiting around the sun, can be described by the system
where the distance r =
With distance in AU (1 AU = 1.496 × 1011 m) and time in years, we have K ≈ 40 for an object rotating around the sun.
Take as initial conditions x(0) = 1, x0 (0) = 0 and y(0) = 0, y0 (0) = 2, and solve for 0 ≤ t ≤ 4.
Investigate the effects of different values of y 0 (0).
5. The equations for the deflection y and rotation z or a simply supported beam with a uniformly distributed load of intensity 9 · 103 N/ft and bending moment M(x) = 10x−x² can be expressed as
where E is the modulus of elasticity, and I is the moment of inertia of the cross section of the beam.
Taking EI = 1.6 · 107 N and y(0) = 0 and z(0) = −0.02, find y and z for 0 < x < 10.
6. Numerically compare the following 4th order Runga-Kutta-Gill method with the classical Runge-Kutta method.