Does the limit of $\frac{\partial f}{\partial u}$ at $u=0$ exist?

For an optimization routine I needed to compute the derivative of the right-hand side $\: f_u(x_k, u_k)$ of a discrete-time system $x_{k+1} = f(x_k, u_k)$. Since $\: f_u(x_k, u_k)$ includes terms that are divided by $u_k$, it is not possible to evaluate $f(x_k, u_k=0)$. From the physics, it's pretty clear that a limit must exist. The numerics poses a problem here, though. Any idea, how I could circumvent it?

To illustrate the behavior, I added a double logarithmic plot showing $f_u(x_k, u_k)$ as $u \rightarrow 0$. It doesn't make any sense that $f_u$ escapes here. What the plot does not show is that there are multiple sign changes of $f_u$ as $u$ drops below $10^{-4}$, let's say. The function $f_u(x_k, u_k)$ is defined as follows:

$$f_u(x, u) = \partial_u \int_0^T \int_0^t \sin(x_5 + x_6 \tau + u \frac{\tau^2}{2}) \: d\tau.$$

I had a CAS do the job of analytically integrating and derivating the sine term. As already pointed out this produces some division by $u$ which causes the problems when it comes to evaluating at $u=0$.

• Do you have the analytical form of $f$? If so, can you post it here? Jan 24 '16 at 15:41
• I would segregate out the part that involves division by $u$ (which might be smaller than the 2800+ characters long). To that part I would try to apply l'Hopital's Rule to get the limit as $u\to 0$. Jan 24 '16 at 17:21
• I don't think I understand -- is $f(x,u=0)$ finite, or does it go to infinity? In the former case, you can apply l'Hopital's rule, as suggested. In the latter case, however, the function is not continuous at $u=0$, and consequently also not differentiable. Jan 25 '16 at 0:47
• Can you evaluate $f(x,u)$ for small $u$, or is the error just as bad as for $f_u$? Is the function holomorphic? I feel the question would be much clearer if you said more about the function. Even at 3k characters, can you at least put it into something like a gist (gist.github.com)? Jan 25 '16 at 18:00
• @MaxHerrmann I got the full series for the integrand of $f_u$ in $u$ as $$\sum_{n\geq1} \frac{(-1)^{\lfloor n/2\rfloor}}{2^n (n-1)!}\tau^{2n}Q_n u^{n-1},$$ where $Q_n$ is $\cos(x_5+x_6\tau)$ (odd $n$), or $\sin(x_5+x_6\tau)$ (even $n$). I expect it's more accurate to take more terms than just the first one. Incidentally, if you've answered your own question with this, please submit a full answer — [SE encourages people to answer their own questions.] Jan 26 '16 at 14:40

Given that $u \frac{\tau^2}{2} \ll 1$, one way of tackling the numerical oscillations, even before the actual term emerges, is a Taylor approximation in $u$ of the sine term (Thanks to Kirill for the full series):

$$\sin(\underbrace{x_5 + x_6 \tau}_{u_0} + \underbrace{u \frac{\tau^2}{2}}_{\Delta u}) = \sin(x_5 + x_6 \tau) + \sum_{n=1}^N \frac{(-1)^{\lfloor n/2 \rfloor}}{2^n(n-1)!} \tau^{2n}Q_nu^{n} + o(\Delta u^{N}),$$

where

$$Q_n = \begin{cases} \cos(x_5 + x_6 \tau) & n \text{ odd} \\ \sin(x_5 + x_6 \tau) & \text{else.} \end{cases}$$

Introducing this approximation into the original equation and exchanging order of integration and differentiation gives

$$\begin{array}{ll} f_u(x,u) & = \int_0^T \int_0^t \partial_u \left( \sin(x_5 + x_6 \tau) + \sum_{n=1}^N \frac{(-1)^{\lfloor n/2 \rfloor}}{2^n(n-1)!} \tau^{2n}Q_nu^{n} + o(\Delta u^{N}) \right) \: d\tau \\ & = \int_0^T \int_0^t \sum_{n=1}^N \frac{(-1)^{\lfloor n/2 \rfloor}}{2^n(n-1)!} \tau^{2n}Q_nu^{n-1} + o(\Delta u^{N-1}) \: d\tau. \end{array}$$

A first-order approximation ($N=1$) is then given as

$$\begin{array}{ll} f_u(x,u) & \approx \int_0^T \int_0^t \cos(x_5 + x_6 \tau) \frac{\tau^2}{2} \: d\tau \\ & = \frac{-3 \cos(x_5)}{x_6^4} + \frac{(6-T^2x_6^2)\cos(x_5 + Tx_6) + 2Tx_6(\sin(x_5) + 2 \sin(x_5 + Tx_6))}{2x_6^4}, \end{array}$$

which does not depend on $u$, naturally.