Rory Winston wrote:
>
> I am currently (for pedagogical purposes) writing a simple numerical
> analysis library in R. I have come unstuck when writing a simple
> Newton-Raphson implementation, that looks like this:
>
> f <- function(x) { 2*cos(x)^2 + 3*sin(x) + 0.5 }
>
> root <- newton(f, tol=0.0001, N=20, a=1)
>
> My issue is calculating the symbolic derivative of f() inside the
> newton() function.
>
If it's pedagogical, maybe returning to basics could help.
What is f'(x)?
It's the limit of (f(x + h) - f(x)) / h when h tends to zero.
So, do it numerically: take a sufficiently small h and
compute the limit. h must be small enough
that h^2 f''(x) is much smaller than h f'(x), but big
enough that f(x+h) is not f(x)
numerical.derivative <- function(f, x, h = 0.0001)
{
# test something
(f(x + h) - f(x)) / h
}
Ex:
numerical.derivative(cos, pi) = 5e-05 # should be -sin(pi) = 0
numerical.derivative(sin, pi) = -1 # ok
numerical.derivative(exp, 0) = 1.00005 # close enough
numerical.derivative(sqrt, 0) = 100 # should be Inf
Alberto Monteiro
______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
