JuliaMath / JuliaMath/Calculus.jl
problems with finite difference
- Dominant language
- Julia
- Stars
- 285
- Forks
- 78
- PR merge metrics
- No merged PRs in 30d
Description
I think there is a problem with the `central` difference code.
Here is a simply example:
```
julia> Calculus.derivative(x -> x/(x+ 1.4424183196362515e-9),2e-8,:central)
-39.33717713979761
```
```
julia> Calculus.derivative(x -> x/(x+ 1.4424183196362515e-9),2e-8,:forward)
1.8509290259770826e6
```
We know that the analytical derivative is
```
1.4424183196362515e-9/(x+1.4424183196362515e-9)
```
which at `x=2e-8` is `3.137210795286552e6`, similar to the `forward difference` result (at least on the same order of magnitude) but there is a 5 orders of magnitude difference compared to the `central difference` result, which is supposed to be `2nd` order accurate!
I looked through the source code, and I think the problem is with the macro `@centralrule`
```
julia> @centralrule 2e-8 epsilon
6.0554544523933395e-6
```
in comparison
```
julia> @forwardrule 2e-8 epsilon
1.4901161193847656e-8
```
The stepping given by `@centralrule` is clearly unacceptably large (2 orders of magnitude greater than `x`!. I don't quite understand why the stepping in the `@centralrule` is taken as the cubic square root of `eps()` while the `@forwardrule` uses the square root of `eps()`?
```
macro forwardrule(x, e)
x, e = esc(x), esc(e)
quote
$e = sqrt(eps(eltype($x))) * max(one(eltype($x)), abs($x))
end
end
macro centralrule(x, e)
x, e = esc(x), esc(e)
quote
$e = cbrt(eps(eltype($x))) * max(one(eltype($x)), abs($x))
end
end
```
In this case an appropriate stepping show be smaller than the order of `1e-10`, which nether of these two methods gives. Why not just use `eps()` as the step size?
Contributor guide
No contributing guide indexed for this repository
Assessment
This issue has not been assessed yet.