oce (version 0.1-82)

oce.bisect: Find the root of a function by bisection

Description

Find the root of a function by bisection

Usage

oce.bisect(f, xleft, xright, dx, debug=getOption("oce.debug"))

Arguments

f

a function of a single variable, returning a single variable.

xleft

a value of the independent variable, to the left of a root.

xright

a value of the independent variable, to the right of a root.

dx

optional value indicating the desired tolerance of the independent variable, near the root. If not provided, this is taken as (xlight-xleft)/1e5.

debug

a debugging parameter, set nonzero to make oce.bisect print a trail of its choices of independent variable, as it searches for the root.

Value

• A value of the independent variable that is within eps of a root.

Details

The root is found by a classical bisection search. That is, the first step involves defining xmiddle as the mid-point between xleft and xright. Then oce.bisect checks whether the sign of the function at xmiddle matches that at xleft. If the sign of the function at xmiddle is the reverse of the sign at xleft, then a root must be between these two points, and so xmiddle is taken as the new xright, and the process is repeated. (The other case is handled similarly.) The procedure continues until the width of the bracketting range gets smaller than dx.

The process is fast, because every step reduces the search range by a factor of 2. Note that functions with multiple roots are not well-handled by oce.bisect (or any bisection algorithm), because only one root will be found.

Examples

hfx.sun.angle <- function(t) sun.angle(t, lat=44+39/60, lon=-(63+36/60))$elevation
sunrise <- oce.bisect(hfx.sun.angle, 
                      as.numeric(as.POSIXct("2011-07-01 05:00:00", tz="UTC")), 
                      as.numeric(as.POSIXct("2011-07-01 10:00:00", tz="UTC")))
cat("Sunrise on Canada Day, 2011 is at", format(number.as.POSIXct(sunrise)), "UTC 
")