Overview
A common problem is the following:
Given a continuous function \(f\) defined on the real numbers, find values such that \(f(x)=0\).
This is called the root finding problem. It comes up, for example, in finding the equilibrium or steady-state solutions for a one-dimensional autonomous system. Note that the roots of a function \(y=f(x)\) coincide exactly with the \(x\)-intercepts of the graph of \(y=f(x)\). In general, when \(f\) is a nonlinear function, the root finding problem is intractable in terms of finding exact, closed-form solutions. However, algorithms for finding approximate solutions to the root finding problem is a major theme in numerical analysis.
In this post, we discuss practical numerical root-finding by way of the R package rootSolve. The documentation for this package is available here.
Practical Root-Finding
Suppose we would like to find the equilibria for the one-dimensional autonomous system
\(\frac{dx}{dt}=f(x) = x^2 - 1\)
While it’s obvious what the answer is, let us address the problem computationally for illustrative purposes. We begin by plotting the function \(y=f(x)\) in order to check for the existence of roots “by inspection.”
The graph tells us that there are two roots in the interval \([-2,2]\). Let’s call the function
uniroot.all from rootSolve:
uniroot.all(function(x){x^2-1}, c(-2, 2))
[1] -1 1Notice that this in fact returns the roots for \(f(x)=x^2-1\).
What uniroot.all attempts to do is to find all of the
roots in the specified interval. Thus, this provides a handy tool for
finding equillbria for one-dimensional autonomous systems.
The rootSolve package is by no means the only option for
numerical root-finding in R. We choose to highlight it as it is
developed by the group that developed the deSolve package
for the numerical approximation for differential equations.