Thank you for the comment, David!
I changed the source code as below and it is really fast indeed.
function hder(f, n, a)
t = taylor1_variable(n);
g(x) = f(x+t);
hod = get_coeff(g(a), n) * gamma(n+1);
return hod
end
f30(x) = hder(f, 30, x)
f40(x) = hder(f, 40, x)
f60(x) = hder(f, 60, x)
Now it takes about 0.000155 seconds to calculate 40th or 60th order
derivatives.
p.s. I learned Julia from your tutorial video at Scipy 2014 :)
Thanks for all of the help.
On Sunday, July 17, 2016 at 8:18:30 PM UTC+1, David P. Sanders wrote:
>
> Here is how to use TaylorSeries.jl. In the next version, the syntax
> taylor1_variable will change to just Taylor1:
>
> julia> t = taylor1_variable(40)
> 1.0 t + 𝒪(t⁴¹)
>
> julia> f(x) = sin(exp(x))
> f (generic function with 1 method)
>
> julia> f(t)
> 0.8414709848078965 + 0.5403023058681398 t - 0.15058433946987837 t² -
> 0.4207354924039483 t³ - 0.3229307265911699 t⁴ - 0.13861923299172246 t⁵ -
> 0.016963650188307994 t⁶ + 0.027151311649628966 t⁷ + 0.028082409487798127 t⁸
> + 0.016534675299332297 t⁹ + 0.006744265666861196 t¹⁰ +
> 0.0016199249536140723 t¹¹ - 0.00020257707395350685 t¹² -
> 0.0005026787959383626 t¹³ - 0.00034829060136926197 t¹⁴ -
> 0.00016780164241401782 t¹⁵ - 6.0803736936094484e-5 t¹⁶ -
> 1.4645171902851875e-5 t¹⁷ - 1.1346823483603049e-8 t¹⁸ +
> 2.515198118591517e-6 t¹⁹ + 1.817031683619697e-6 t²⁰ + 8.873827219817266e-7
> t²¹ + 3.383697539084375e-7 t²² + 9.90873953752995e-8 t²³ +
> 1.7331465080481204e-8 t²⁴ - 2.8377431603022023e-9 t²⁵ -
> 4.478125328289587e-9 t²⁶ - 2.633876550034189e-9 t²⁷ - 1.1444925030062194e-9
> t²⁸ - 4.013808037099447e-10 t²⁹ - 1.1102527569232767e-10 t³⁰ -
> 1.95623844080032e-11 t³¹ + 1.8017173073914597e-12 t³² +
> 3.715400951669623e-12 t³³ + 2.1778749843678503e-12 t³⁴ +
> 9.347841028230405e-13 t³⁵ + 3.279923193009391e-13 t³⁶ +
> 9.423306246867712e-14 t³⁷ + 1.997073344234444e-14 t³⁸ +
> 1.3114196619133462e-15 t³⁹ - 1.5394356646368331e-15 t⁴⁰ + 𝒪(t⁴¹)
>
> julia> get_coeff(f(t), 40)
> -1.5394356646368331e-15
>
> Note that this coefficient is the 40th derivative divided by 40!. (The
> nth Taylor coefficient is f^{(n)} / n!. )
>
>
>
> El domingo, 17 de julio de 2016, 21:15:33 (UTC+2), David P. Sanders
> escribió:
>>
>>
>>
>> El domingo, 17 de julio de 2016, 17:13:38 (UTC+2), Young Chun escribió:
>>>
>>> I have an optimization problem that involves series of high order
>>> derivatives of f(x).
>>> So, to get a decent value, I need to calculate higher order derivative
>>> at least to 40~50th order.
>>>
>>
>> I suggest you use TaylorSeries.jl.
>>
>> Note that ForwardDiff from version 0.2 removed ForwardDiff.derivative(f);
>> this must now be written as the anonymous function x ->
>> ForwardDiff.derivative(f, x).
>>
>>
>>>
>>> Initially, I tried to use ForwardDiff package with automatic
>>> differentiation but kept facing the following error message
>>>
>>> <div id="inner-editor"><br class="Apple-interchange-newline">using
>>> DualNumbers using ForwardDiff f(x::Vector) = sum(sin, x) + prod(tan, x) *
>>> sum(sqrt, x); x = rand(5) g = ForwardDiff.derivative(f);</div>
>>>
>>> using DualNumbers
>>>
>>> using ForwardDiff
>>>
>>> f(x::Vector) = sum(sin, x) + prod(tan, x) * sum(sqrt, x);
>>>
>>> x = rand(5)
>>>
>>> g = ForwardDiff.derivative(f);
>>>
>>> LoadError: MethodError: `derivative` has no method matching
>>> derivative(::Function)
>>> Closest candidates are:
>>> derivative(::Any, !Matched::Any)
>>> while loading In[28], in expression starting on line 4
>>>
>>>
>>>
>>> So, I moved on to D(f) instead from Root package and it works quite well
>>> for a small order (like 10th order derivative)
>>>
>>> <div id="inner-editor"><br class="Apple-interchange-newline">using
>>> DualNumbers using ForwardDiff using Roots using ForwardDiff using JuMP #
>>> Need to say it whenever we use JuMP using Gadfly function hder(f, n::Real)
>>> # subroutine to calculate higher order derivative temp = f; if n==0 return
>>> temp else for i=1:1:n temp = D(temp) end return temp end end</div>
>>>
>>> using DualNumbers
>>>
>>> using ForwardDiff
>>>
>>> using Roots
>>>
>>> using ForwardDiff
>>>
>>>
>>>
>>> using JuMP
>>>
>>> using Gadfly
>>>
>>>
>>>
>>> function hder(f, n::Real) # subroutine to calculate n-th order derivative
>>>
>>> temp = f;
>>>
>>> if n==0
>>>
>>> return temp
>>>
>>> else
>>>
>>> for i=1:1:n
>>>
>>> temp = D(temp)
>>>
>>> end
>>>
>>> return temp
>>>
>>> end
>>>
>>> end
>>>
>>> Out[7]:
>>>
>>> hder (generic function with 1 method)
>>>
>>> f(x) = (cos(x)) * exp(-1/5 * x)
>>> g0 = hder(f, 0)
>>> g1 = hder(f, 1)
>>> g2 = hder(f, 2)
>>> g3 = hder(f, 3)
>>> g4 = hder(f, 4)
>>> g5 = hder(f, 5)
>>> g6 = hder(f, 6)
>>> g7 = hder(f, 7)
>>> g8 = hder(f, 8)
>>>
>>> plot([g0, g1, g2, g3, g4, g5, g6, g7, g8], 0, 5pi)
>>>
>>>
>>> <https://lh3.googleusercontent.com/-RCrcnLVFPhE/V4tfcG9uSEI/AAAAAAAAGUU/unxjN7x60jgPaEcCmcRlykx8WpbMxJBnACLcB/s1600/Screen%2BShot%2B2016-07-17%2Bat%2B11.35.09.png>
>>> However, if I push the limit and ask to calculate, like 30th order
>>> derivative, the program never ends and keep calculating for hours.
>>> Is there a better way to do this type of task? either by using
>>> ForwardDiff package or by modifying my function?
>>> I just learned Julia last week, so please understand if my question
>>> sounds stupid.
>>>
>>>
>>>
>>>
