# Unbounded solution but bounded Euler discretization

Is there an ordinary differential equation in $$\mathbb{R}^d$$ induced by a gradient vector field with unbounded solutions, for which the difference equations obtained by using the forward Euler method have all solutions bounded?

I read the paper Guy, Richard T.; Misiurewicz, Michał, Euler approximations can destroy unbounded solutions, Fasc. Math. 44, 43-52 (2010). ZBL1215.65124 (pdf also available here). In this paper, the authors construct examples of differential equations with the form $$\dot{x}=F(x)$$ that has an unbounded (in forward time) solution $$x(t)$$, but its Euler approximations for sufficiently small $$\epsilon>0$$ $$x_n=x_{n−1}+\epsilon F(x_{n−1})$$ have only bounded solutions. However, all examples they propose are not induced by a gradient vector field, i.e., we cannot write $$F(x)=\nabla G(x)$$ for some $$G(x)$$, where $$\nabla$$ stands for the gradient operator. It makes me think whether such phenomenon could happen even for a gradient vector field.

My guess is yes, but it seems hard to find an example that can be verified rigorously by analysis, not just by numerical experiments. Can anyone give a hint on it? It would be better if the example is a polynomial.

• That’s really interesting. By conservative system, it seems like you mean “gradient system”. Conservative system typically means that the dynamic conserves something like a first integral. Commented Sep 3, 2022 at 16:10
• @NawafBou-Rabee Thanks for clarifying, I am confused with those two concepts, I have changed it to gradient system. Commented Sep 3, 2022 at 17:51
• Here is a classical reference that clarifies these terms: books.google.com/books/about/…. See Chapter 2 Commented Sep 3, 2022 at 17:58
• @NawafBou-Rabee Thanks for the reference. Commented Sep 3, 2022 at 20:38

The central problem with forward Euler is that it doesn’t continuously update the vector field, and therefore can tunnel through fine features of the vector field. The gradient case is no exception. For concreteness, consider the 1D potential $$G(x) =(1/2) x^2 + \sum_{i=1}^{\infty} \int_0^x \eta_{2^{-2-i}}(s - 1 + 2^{-i}) ds$$ where $$\eta_{\delta}(\cdot)$$ is a smooth nascent delta function with support in $$(-\delta/2, \delta/2)$$. Note that the supports of these nascent delta functions are disjoint. Starting say from $$1/3$$, the corresponding continuous solution is infinitely expansive in finite time and hence unbounded; whereas, for any step size $$\epsilon>0$$, the corresponding Euler solution will tunnel through an infinite # of nascent delta functions and hence remain bounded over a finite time window.

In the comments, the OP shared a potential counterexample, which was brought up in a recent MSE post on August 25th 2022 that remains unanswered. In this "counterexample", the underlying potential/loss/objective/cost function is: $$G(x,y) = \frac{1}{2} x^2 y^2 - x - y \;.$$ As noted in the post, based on numerical evidence, it looks like the forward Euler solution is bounded, but the exact solution seems unbounded. These numerical observations can be validated (or refuted) using a standard stability analysis of the continuous/discrete dynamics outlined below.

First, note that the exact gradient dynamics has only one fixed point at $$(1,1)$$. The eigenvalues of the Hessian evaluated at this fixed point $$\nabla^2 G(1,1)$$ are $$3$$ and $$-1$$ with corresponding eigenvectors $$(1,1)$$ and $$(-1,1)$$, therefore the fixed point is linearly unstable, and more precisely the instability is in the antidiagonal direction.

To determine the nonlinear stability of the fixed point in the anti diagonal direction, consider the Lyapunov function $$\mathcal{V}(x,y) = (x-y)^2 \;,$$ which measures the distance to the diagonal in the $$xy$$-plane. A direct computation shows that $$\mathcal{L} V(x,y) = 2 x y \mathcal{V}(x,y) \;,$$ where $$\mathcal{L}$$ is the differential operator of the exact gradient dynamics. This shows that the distance from the diagonal is non-decreasing if $$x$$ and $$y$$ have the same sign, and in particular, increasing in a neighborhood of the fixed point $$(1,1)$$, until either \begin{align*} &y = 1/x^2 \implies \dot{y}=0 ~~\text{and}~~ \dot{x}>0 ~~ \text{or} \\ & x=1/y^2 \implies \dot{x}=0 ~~ \text{and} ~~ \dot{y}>0 \;. \end{align*} In other words, the exact dynamics becomes (nonlinearly) unstable due to the gravity-like terms in $$G(x,y)$$, but this instability is weak (non-exponential).

On the other hand, the situation is a bit different along the forward Euler solution. The fixed point is still unstable with eigenvalues of the linearized (discrete) dynamics being $$1-3 \epsilon$$ and $$1+\epsilon$$ with the same eigenvectors. To assess the nonlinear stability in the anti diagonal direction, we use the same Lyapunov function $$\mathcal{V}(x_{n+1},y_{n+1}) = (1+\epsilon x_n y_n) \mathcal{V}(x_n, y_n)$$ Thus, the distance to the diagonal is increasing in the region $$\{ x, y : |1 + \epsilon x y | > 1 \}$$ which includes the region where $$x y >0$$ but also $$x y < -2 /\epsilon$$. At a fixed step size, for small time interval the Euler solution essentially behaves like the exact solution, but for larger time intervals instead of asymptoting to one of the following trajectories \begin{align*} &y = 1/x^2 \implies \dot{y}=0 ~~\text{and}~~ \dot{x}>0 ~~ \text{or} \\ & x=1/y^2 \implies \dot{x}=0 ~~ \text{and} ~~ \dot{y}>0 \;, \end{align*} it can catastrophically lose stability and blow up because the component that is small can change sign and the dynamics in that component becomes unstable, and in turn, the dynamics in the other component grows without bound.

This is borne out numerically. Here is a a trajectory of forward Euler at $$\epsilon = 1/6$$ in the $$xy$$-plane over the time interval $$[0,10]$$ initiated at $$(1.6,1.5)$$. The size of the dots increase with time. At first, the Euler solution is essentially following the continuous dynamics and everything looks fine. However, note that when the $$y$$-component becomes small, and because forward Euler does not continuously update the gradient, the $$y$$-component can change sign, which leads to an unexpected/dramatic/artificial instability.

• Nice example, but it would be better if the function is more elementary, like a polynomial (maybe 2-d, because 1-d polynomial seems doesn't work). Commented Sep 7, 2022 at 14:46
• It might be possible to construct a 2D polynomial example, but the failure might be less spectacular. Importantly, the construction will follow the same principle: forward Euler can fail to "see" important features of the vector field if the step size is too big. Commented Sep 7, 2022 at 15:18
• Thanks. Now I understand your point. I would like an example so that forward Euler fail for all sufficiently small step size, i.e., I fix a function, and there exists $\alpha_0>0$ such that for any step size $\alpha\in(0,\alpha_0]$, forward Euler is bounded even if the continuous solution is unbounded (but the bound could depend on $\alpha$). Your example seems to say that I first choose any small step size and then find a corresponding function. Commented Sep 7, 2022 at 19:03
• This example seems to satisfy that requirement, because for any $\epsilon>0$ forward Euler will only see a finite number of the nascent delta functions; whereas the true solution will see an infinite number of them. Commented Sep 7, 2022 at 19:33
• Your numerical experiment is unbounded because you take too large step size for the initial point you choose. Since the gradient of $G(x)$ does not have a global Lipschitz constant, we need to select small enough step size for each different initial point (that is what I mean sufficiently small, not a single small step size for all initial point, but one step size for one initial point). I used your initial point, but reduced $\epsilon$ to be $\epsilon=0.01$, and run $10000$ iterations, I saw all iterates lie in a bounded region and some periodic pattern occurred. Commented Sep 14, 2022 at 4:03