# multiroot.1D: Solves for n roots of n (nonlinear) equations, created by... In rootSolve: Nonlinear Root Finding, Equilibrium and Steady-State Analysis of Ordinary Differential Equations

## Description

multiroot.1D finds the solution to boundary value problems of ordinary differential equations, which are approximated using the method-of-lines approach.

Assumes a banded Jacobian matrix, uses the Newton-Raphson method.

## Usage

 ```1 2 3 4``` ```multiroot.1D(f, start, maxiter = 100, rtol = 1e-6, atol = 1e-8, ctol = 1e-8, nspec = NULL, dimens = NULL, verbose = FALSE, positive = FALSE, names = NULL, parms = NULL, ...) ```

## Arguments

 `f ` function for which the root is sought; it must return a vector with as many values as the length of `start`. It is called either as `f(x, ...)` if `parms = NULL` or as `f(x, parms, ...)` if `parms` is not `NULL`. `start ` vector containing initial guesses for the unknown x; if `start` has a name attribute, the names will be used to label the output vector. `maxiter ` maximal number of iterations allowed. `rtol ` relative error tolerance, either a scalar or a vector, one value for each element in the unknown x. `atol ` absolute error tolerance, either a scalar or a vector, one value for each element in x. `ctol ` a scalar. If between two iterations, the maximal change in the variable values is less than this amount, then it is assumed that the root is found. `nspec ` the number of *species* (components) in the model. If `NULL`, then `dimens` should be specified. `dimens` the number of *boxes* in the model. If NULL, then `nspec` should be specified. `verbose ` if `TRUE`: full output to the screen, e.g. will output the steady-state settings. `positive` if `TRUE`, the unknowns y are forced to be non-negative numbers. `names ` the names of the components; used to label the output, which will be written as a matrix. `parms ` vector or list of parameters used in `f`. `... ` additional arguments passed to function `f`.

## Details

`multiroot.1D` is similar to `steady.1D`, except for the function specification which is simpler in `multiroot.1D`.

It is to be used to solve (simple) boundary value problems of differential equations.

The following differential equation:

0=f(x,y,y',y'')

with boundary conditions

y(x=a) = ya, at the start and y(x=b)=yb at the end of the integration interval [a,b] is approximated

as follows:

1. First the integration interval x is discretized, for instance:

`dx <- 0.01`

`x <- seq(a,b,by=dx)`

where `dx` should be small enough to keep numerical errors small.

2. Then the first- and second-order derivatives are differenced on this numerical grid. R's `diff` function is very efficient in taking numerical differences, so it is used to approximate the first-, and second-order derivates as follows.

A first-order derivative y' can be approximated either as:

y'=`diff(c(ya,y))/dx`

if only the initial condition ya is prescribed,

y'=`diff(c(y,yb))/dx`

if only the final condition, yb is prescribed,

y'=`0.5*(diff(c(ya,y))/dx+diff(c(y,yb))/dx)`

if initial, ya, and final condition, yb are prescribed.

The latter (centered differences) is to be preferred.

A second-order derivative y” can be approximated by differencing twice.

y”=`diff(diff(c(ya,y,yb))/dx)/dx`

3. Finally, function `multiroot.1D` is used to locate the root.

See the examples

## Value

a list containing:

 `root ` the values of the root. `f.root ` the value of the function evaluated at the `root`. `iter ` the number of iterations used. `estim.precis ` the estimated precision for `root`. It is defined as the mean of the absolute function values (mean(abs(`f.root`))).

## Note

`multiroot.1D` makes use of function `steady.1D`.

It is NOT guaranteed that the method will converge to the root.

## Author(s)

Karline Soetaert <[email protected]>

## See Also

`stode`, which uses a different function call.

`uniroot.all`, to solve for all roots of one (nonlinear) equation

`steady`, `steady.band`, `steady.1D`, `steady.2D`, `steady.3D`, steady-state solvers, which find the roots of ODEs or PDEs. The function call differs from `multiroot`.

`jacobian.full`, `jacobian.band`, estimates the Jacobian matrix assuming a full or banded structure.

`gradient`, `hessian`, estimates the gradient matrix or the Hessian.

## Examples

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64``` ```## ======================================================================= ## Example 1: simple standard problem ## solve the BVP ODE: ## d2y/dt^2=-3py/(p+t^2)^2 ## y(t= -0.1)=-0.1/sqrt(p+0.01) ## y(t= 0.1)= 0.1/sqrt(p+0.01) ## where p = 1e-5 ## ## analytical solution y(t) = t/sqrt(p + t^2). ## ## ======================================================================= bvp <- function(y) { dy2 <- diff(diff(c(ya, y, yb))/dx)/dx return(dy2 + 3*p*y/(p+x^2)^2) } dx <- 0.0001 x <- seq(-0.1, 0.1, by = dx) p <- 1e-5 ya <- -0.1/sqrt(p+0.01) yb <- 0.1/sqrt(p+0.01) print(system.time( y <- multiroot.1D(start = runif(length(x)), f = bvp, nspec = 1) )) plot(x, y\$root, ylab = "y", main = "BVP test problem") # add analytical solution curve(x/sqrt(p+x*x), add = TRUE, type = "l", col = "red") ## ======================================================================= ## Example 2: bvp test problem 28 ## solve: ## xi*y'' + y*y' - y=0 ## with boundary conditions: ## y0=1 ## y1=3/2 ## ======================================================================= prob28 <-function(y, xi) { dy2 <- diff(diff(c(ya, y, yb))/dx)/dx # y'' dy <- 0.5*(diff(c(ya, y)) +diff(c(y, yb)))/dx # y' - centered differences xi*dy2 +dy*y-y } ya <- 1 yb <- 3/2 dx <- 0.001 x <- seq(0, 1, by = dx) N <- length(x) print(system.time( Y1 <- multiroot.1D(f = prob28, start = runif(N), nspec = 1, xi = 0.1) )) Y2<- multiroot.1D(f = prob28, start = runif(N), nspec = 1, xi = 0.01) Y3<- multiroot.1D(f = prob28, start = runif(N), nspec = 1, xi = 0.001) plot(x, Y3\$root, type = "l", col = "green", main = "bvp test problem 28") lines(x, Y2\$root, col = "red") lines(x, Y1\$root, col = "blue") ```

rootSolve documentation built on May 29, 2017, 3:28 p.m.