# How to take the second derivative of a complicated function using limits [closed]

I have a very complicated function of one variable consisting of exponential and error functions combined in various ways. Mathematica can't take the derivative of this function in the usual way (using prime or D) because it makes use of the function Im. Following the discussion here, I decided to use the limit definition of the derivative instead of prime or D:

dn[\[Delta]_] = Assuming[Element[\[Delta], Reals], Limit[(Realn[(\[Delta] + \[CurlyEpsilon])] - Realn[\[Delta]])/\[CurlyEpsilon], \[CurlyEpsilon] -> 0]];


This works for the first derivative, but returns a function which is infinity everywhere for the second derivative:

ddn[\[Delta]_] =
Assuming[Element[\[Delta], Reals],
Limit[(Realn[(\[Delta] + \[Xi])] - 2*Realn[\[Delta]] +
Realn[(\[Delta] - \[Xi])])/\[Xi]^2, \[Xi] -> 0]];


at any point. While trying to debug this problem, I discovered that multiplying the numerator and the denominator in the definition of the first derivative by 2 causes it too to give a function which is infinity at all points:

dn[\[Delta]_] = Assuming[Element[\[Delta], Reals], Limit[(2 Realn[(\[Delta] + \[CurlyEpsilon])] - 2 Realn[\[Delta]])/(2 \[CurlyEpsilon]), \[CurlyEpsilon] -> 0]];


and which I am utterly unable to explain. I can do both derivatives numerically using a finite number of points, and the results look reasonably well-behaved (i.e finite and smooth) within the region where the function is varying (it approaches zero at both infinities). Besides, the function made by summing terms which are all infinitely differentiable, so I don't see why the analytical method doesn't work...

• I can't reproduce your problem when setting for example Realn[x_]:=Re[Gamma[I x]]. Could you add an example that shows the issue? – Jens Aug 10 '15 at 23:19
• I don't know how to produce a function that reproduces the error and that can be shared here. Even copying straight from mathematica's output doesn't work, and produces a function that mathematica won't compute the limit of at all. I suspect that this is due to needing arbitrary precision in some of the coefficients for everything to work. (In fact mathematica won't even plot the function copied from the output.) – Mat6554 Aug 11 '15 at 15:22
• If it works for Jens's example and not for yours, I've got to think that the problem has to do with something specific in your function. Any chance of narrowing it down to an example that can be shared? – Michael E2 Aug 11 '15 at 17:44