Notes on Difference Equations (1)

(This part mainly discusses dynamics of first-order difference equations)

1.1 Introduction

$x(n+1) = f(x(n)) \tag{1.1.1}$

(Time-invariant) difference equations and discrete dynamical systems represents two sides of the same coin.

1.2 Linear First-Order Difference Equations

Homogeneous FOE:
$x(n+1) = a(n) x(n) \tag{1.2.1}$
Non-homogeneous FOE:
$x(n+1) = a(n) x(n) + g(n) \tag{1.2.2}$
Important Cases
$y(n+1) = ay(n) + g(n), \quad y(0) = y_0 \tag{1.2.3}$
The solution is given by
$y(n) = a^n y_0 + b[\frac{a^n-1}{a-1}], \quad if \quad a\neq1 \tag{1.2.4}$
Notice that the solution of nonhomogeneous differential equation
$\frac{dy}{dt} = ay(t) + g(t), \quad y(0) = y_0 \tag{1.2.5}$
is given by
$y(t) = e^{at} y_0 + \int_0^t e^{a(t-s)}g(s) ds \tag{1.2.6}$

1.3 Equilibrium Points

What is the equilibrium point?

Definition 1: (Equilibrium Points) A point $x^{*}$ in the domain of $f$ is said to be an equilibrium point if it is a fixed point of $f$ , i.e., $f(x^*) = x^*$ .

Definition 2: (Eventually Equilibrium Points) Let $x$ be a point in the domain of $f$ . If there exists a positive integer $r$ and an equilibrium point $x^{*}$ , such that $f^r(x) = x^{*}$ , $f^{r-1}(x) \neq x^{*}$ , then $x$ is an eventually equilibrium (fixed) point.

Stability of the equilibrium point

Definition 3

(a) The EP $x^{*}$ is stable if given $\varepsilon > 0$ there exists $\delta > 0$ such that $|x_0 - x^*| < \delta$ implies $|f^n(x_0) - x^*| < \varepsilon$ for all $n > 0$ .

(b) The point $x^*$ is said to be attracting if there exists $\eta > 0$ such that

$|x(0) - x^*| < \eta$ implies $\lim_{n \to \infty}(x_n) = x^{*}$ .

If $\eta = \infty$ , $x^*$ is called global attractor.

If $\eta = \infty$ , $x^*$ is said to be globally asymptotically stable.

To determine the stability of an EP from these definitions may prove to be a mission impossible in many cases. This is due to the fact that we may not be able to find the solution in a closed form for most equations.

The simplest bur most powerful tools are graph techniques.

1.4 Criterion for the Asymptotic Stability of Equilibrium Points

Theorem 1: Let $x^{*}$ be an EP of the difference equation
$x(n+1) = f(x(n)) \tag{1.4.1}$
where $f$ is continuously differentiable at $x^{*}$ . The following statements then hold true:

(i) If $|f'(x^*)| < 1$ , then $x^*$ is asymptotically stable.

(ii) If $|f'(x^*)| > 1$ , then $x^*$ is unstable.

Theorem 2: Suppose $f'(x^*) = 1$ . The following statements then hold:

(i) If $f''(x^*) \neq 0$ , then $x^*$ is unstable.

(ii) If $f''(x^*) = 0$ and $f'''(x^*) > 0$ , then $x^*$ is unstable.

(iii) If $f''(x^*) = 0$ and $f'''(x^*) < 0$ , then $x^*$ is asymptotically stable.

1.5 Periodic Points and Cycles

Definition 4: Let $b$ in the domain of $f$ . Then:

(i) $b$ is called a periodic point of $f$ if for some positive integer $k$ , $f^k(b) = b$ . Hence a point is k-periodic if it is a fixed point of $f^k$ , that is, if it is an EP of the difference equation
$x(n+1) = g(x(n)) \tag{1.5.1}$
where $g = f^k$ .

The periodic orbit of $b$ , $O(b) = \{b, f(b), f^2(b), ...,f^{k-1}(b)\}$ is called a $k-cycle$ .

(ii) $b$ is called eventually k-periodic if for some positive integer $m$ , $f^m(b)$ is a k-periodic point. In other words, b is eventually k-periodic if
$f^{m+k}(b) = f^m(b)$
Stability of Periodic Points

Deﬁnition 1.20 Let $b$ be a k-period point of $f$ . Then $b$ is:

(i) stable if it is a stable ﬁxed point of $f^k$

(ii) asymptotically stable if it is an asymptotically stable ﬁxed point of $f^k$

(iii) unstable if it is an unstable ﬁxed point of $f^k$

REFERENCE

Elaydi S N. An Introduction to Difference Equations[M]. 2011.

要对经济动力系统进行研究，就离不开对差分方程问题的讨论。特别是经济中所出现的经济周期等现象，更有可能与差分方程系统的性质有着密切联系，因此，我选择Elaydi (2011)进行阅读与研究，为差分方程方面的技术打基础。

上面是对Elaydi (2011) 第一章的简要总结，当然还有关于周期解的稳定性的讨论以及虫口模型的内容还没有加入进来，因为定理的内容比较复杂，所以等以后有时间再加入整理也好。关于一阶差分方程，这一章讨论的主要问题就是不动点(Equilibrium Point)与周期性，而周期性问题又可以拆解为不动点问题。学了这一章之后，大概了解了：

差分方程的定义
一阶差分方程与微分方程之间的联系
不动点的定义
不动点的稳定性
差分方程的周期性
周期点的稳定性

为了更直观的了解差分方程的性质，我们可以利用程序模拟差分方程的解。

import numpy as np
import matplotlib.pyplot as plt

首先，我们定义差分方程函数，这个函数应当给定其函数形式与初始值之后，可以生成orbits.

# Basic Blocks of Difference Equations
def iter_func(var):
    """ Define the iterative function in difference equations.
        Given args of the function, it returns the value of the calculus.
    """
    if 0 <= var <= (1/2):
        return 2 * var
    elif (1/2) < var <= 1:
        return 2 * (1 - var)

def diff_eq(iter_func, init_val=0,  num=100):
    """ Define difference equations.
        Given the function form and initial value of the equation,
        it returns sereis of solutions of the equation.
        Length of the series is controled by num.
    """
    orb = list(range(num))
    orb[0] = init_val
    for i in range(1, num):
        orb[i] = iter_func(orb[i-1])
    return orb

我们定义一个函数来绘制我们所生成的orbits

def draw_series(series):
    """ Draw the series of orbits
    """
    n = range(len(series))
    fig = plt.figure(1)
    ax = plt.subplot(1,1,1)
    ax.plot(n, series)
    plt.xlabel('n')
    plt.ylabel('x(n)')
    plt.show()

x = diff_eq(iter_eq, 0.2)

draw_series(x)

可以看出，这个差分方程在给定初始值0.2的情况下，经过不断地震荡，最终收敛到0.

output_6_0.png

下面我们通过一个函数绘制图像，帮我们找到一个差分方程的不动点。

#Find Equilibrium Point
def find_EP(iter_eq):
    """ Given the form of the equation,
    it helps to find EP by the graph.
    """
    x = np.linspace(0,1,100)
    viter_eq = np.vectorize(iter_eq)
    y = viter_eq(x)
    
    fig = plt.figure(1)
    ax = plt.subplot(1,1,1)
    ax.plot(x,x)
    ax.plot(x,y)
    plt.show()

find_EP(iter_eq)

通过图像可以看出，不动点有两个，一个是0，另一个在0.6到0.8之间（2/3）。

output_9_0.png

下面我们通过图像来进行相图分析，来观察不动点的稳定性。

# Stair-step Analysis
def stair_step(iter_func):
    """ Presenting the stair step diagram of the
    difference equation.
    """
    x = np.linspace(0,1,100)
    viter_eq = np.vectorize(iter_func)
    y = viter_eq(x)
    
    ssd = diff_eq(iter_func, 0.2)
    h = np.zeros(2 * len(ssd) - 1)
    v = np.zeros(2 * len(ssd) - 1)
    v[0] = 0
    for i in range(0, len(h)-1):
        h[i] = ssd[i // 2]
    for i in range(1, len(v)-1):
        v[i] = ssd[(i+1) // 2]
        
    fig = plt.figure(1)
    ax = plt.subplot(1,1,1)
    ax.plot(h, v, '^--')
    ax.plot(x,x)
    ax.plot(x,y)
    plt.show()

stair_step(iter_eq)

output_12_0.png

可以看出，当 $x$ 靠近0.5的时候很容易被发散到1，进而导致迭代值为0，这与我们之前的直觉也是一致的。