混沌系统在图像加密中的应用（三维连续Lorenz混沌系统及三种离散化method）

前言

本节内容将介绍3维连续混沌系统，随后读者可以根据第四和第五节的内容用连续混沌系统对图像进行加密，并和用离散混沌系统加密的图片进行对比。

一、连续混沌系统和离散混沌系统区别

连续混沌系统是指动力学方程在连续时间范围内描述系统演化的混沌系统。在连续混沌系统中，系统的状态变量（如位置、速度等）是连续变化的，时间是连续的，其动力学方程通常由一组非线性微分方程表示。这些方程描述了系统在连续时间下的行为，系统的状态可以通过连续时间的轨迹来表示。

离散混沌系统是指动力学方程在离散时间范围内描述系统演化的混沌系统。在离散混沌系统中，系统的状态变量和时间都是离散的，在每个时间步长上更新。系统的状态变量根据离散时间的迭代规则进行更新，而不是根据连续时间的微分方程。离散混沌系统的动力学方程通常由一系列非线性差分方程或迭代方程表示。这些方程描述了系统在离散时间下的行为，系统的状态可以通过离散时间的轨迹来表示。

二、Lorenz混沌系统

1.Lorenz混沌系统介绍

Lorenz系统是一种经典的混沌系统，由物理学家Edward Lorenz于1963年首次提出。它是一组非线性的三维常微分方程，用来描述大气层中对流运动的模型。Lorenz系统也被广泛研究和应用于混沌理论、动力系统和非线性科学等领域。

Lorenz系统的方程可以表示为：

其中，x、y、z是系统的状态变量，t是时间，σ、ρ、β是方程的参数。

Lorenz系统具有以下显著特点：

（1）非线性行为：方程中的非线性项使得系统的行为变得复杂且难以预测。
（2）混沌特性：当参数取特定值时，Lorenz系统表现出混沌行为，即具有无法重复的敏感依赖于初始条件的轨迹。这意味着微小的初始条件差异会导致系统最终的状态差异非常大。
（3）吸引子：Lorenz系统的吸引子是一个具有奇异吸引特性的结构，类似于蝴蝶的形状，因此也被称为"Lorenz吸引子"。吸引子上的轨迹形成一个涡旋状，表现出不断变形和蜿蜒回归的特点。

由于其简单而引人注目的行为，Lorenz系统成为了研究混沌现象和非线性科学的重要模型。它对于理解混沌的基本特性、敏感依赖于初始条件和确定性混沌等概念有着重要贡献。Lorenz系统也广泛应用于通信、密码学、天气预报、混沌加密和控制系统等领域，为实际应用提供了理论依据。

2.Lorenz混沌系统离散化

对于连续非线性系统而言，计算机是无法直接绘出微分方程的轨迹图，需要先对连续非线性系统进行离散化，才能搜集每个step的轨迹，有三种方法可以将其离散化

（一）欧拉（Euler）方法

欧拉法是一种常见的数值方法，用于近似求解常微分方程（ODEs）。该方法基于离散化时间的思想，将连续变量在每个时间步骤上进行逐步更新。以下是欧拉法的简要描述：

（1）初始条件：给定初始状态的值。
（2）离散化时间：选择适当的步长（即离散化时间间隔）。
（3）迭代更新：通过近似微分方程的导数来更新系统的状态。 公式如下：

（4）重复迭代：重复步骤（3），直到达到所需的时间点或时间步数。下面是仿真实验

欧拉法是一种简单易懂的数值方法，并且容易实现。然而，它的主要缺点是积累误差，因为它基于导数值的线性近似。随着步长的增大，误差会增加并导致结果的不准确性。为了获得更准确的结果，可能需要采用更高阶的数值方法，如改进的欧拉法、四阶龙格-库塔法等。

（二）改进的欧拉（Euler）方法

改进的欧拉方法（Improved Euler method），也称为梯形方法（Trapezoidal method）或Heun’s方法，是对欧拉法的一种改进和精确度提升的数值方法。改进的欧拉方法通过在每个时间步骤中使用两个导数估计来减小误差。

以下是改进的欧拉方法的简要描述：

（1）初始条件：给定初始状态的值。
（2）离散化时间：选择适当的步长（即离散化时间间隔）。
（3）第一次更新：使用欧拉法进行第一次近似更新。
（4）第二次更新：利用第一次更新得到的值进行第二次近似更新。公示如下：

（5）重复迭代：重复步骤4，直到达到所需的时间点或时间步数。下面是仿真实验

（三）四阶龙格-库塔法（Fourth-Order Runge-Kutta method）

四阶龙格-库塔法（Fourth-Order Runge-Kutta method）是一种常用的数值方法，用于解决常微分方程（ODEs）。它是改进的欧拉方法的一个高阶拓展，具有更高的精度和稳定性。

以下是四阶龙格-库塔法的简要描述：

（1）初始条件：给定初始状态的值。
（2）离散化时间：选择适当的步长（即离散化时间间隔）。
（3）第一次斜率估计：计算初始点的斜率估计。
（4）第二次斜率估计：使用第一次斜率估计得到的点进行更准确的斜率估计。
（5）第三次斜率估计：再次使用斜率估计。
（6）第四次斜率估计：最终斜率估计。公示如下：

其中，

（7）更新下一个时间步骤的值：根据四个斜率估计和步长更新下一个时间步骤的值。
（8）重复迭代：重复步骤3到7，直到达到所需的时间点或时间步数。。下面是仿真实验

四阶龙格-库塔法通过使用四次斜率估计的加权平均值来计算下一个时间步骤的值，从而提供了更高的数值精度。它通常比欧拉法和改进的欧拉方法更准确，尤其在更小的时间步长和具有快速变化的系统方程中。

三、python代码

1.欧拉（Euler）方法python代码

import numpy as np
import matplotlib.pyplot as plt


# 定义Lorenz系统的方程
def lorenz_equations(x, y, z, sigma, rho, beta):
    dx_dt = sigma * (y - x)
    dy_dt = x * (rho - z) - y
    dz_dt = x * y - beta * z
    return dx_dt, dy_dt, dz_dt


# 使用欧拉法进行离散化
def euler_discretization(x0, y0, z0, sigma, rho, beta, dt, num_steps):
    # 初始化变量数组
    x = np.zeros(num_steps + 1)
    y = np.zeros(num_steps + 1)
    z = np.zeros(num_steps + 1)

    # 初始条件
    x[0] = x0
    y[0] = y0
    z[0] = z0

    # 迭代计算
    for i in range(num_steps):
        dx_dt, dy_dt, dz_dt = lorenz_equations(x[i], y[i], z[i], sigma, rho, beta)
        x[i + 1] = x[i] + dx_dt * dt
        y[i + 1] = y[i] + dy_dt * dt
        z[i + 1] = z[i] + dz_dt * dt

    return x, y, z


# 定义Lorenz系统的参数
sigma = 10.0
rho = 28.0
beta = 8.0 / 3.0

# 定义离散化的参数
dt = 0.01
num_steps = 10000

# 设置初始条件
x0 = 0.1
y0 = 0.0
z0 = 0.0

# 进行欧拉法离散化
x, y, z = euler_discretization(x0, y0, z0, sigma, rho, beta, dt, num_steps)

# 绘制Lorenz系统的轨迹
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot(x, y, z, lw=0.5)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
ax.set_title('Lorenz System - Euler Discretization')
plt.show()

2.改进的欧拉（Euler）方法python代码

import numpy as np
import matplotlib.pyplot as plt

def lorenz(t, x, sigma, rho, beta):
    """Lorenz系统的微分方程"""
    dx = sigma * (x[1] - x[0])
    dy = x[0] * (rho - x[2]) - x[1]
    dz = x[0] * x[1] - beta * x[2]
    return np.array([dx, dy, dz])

def improved_euler_lorenz(sigma, rho, beta, x0, dt, num_steps):
    """使用改进的欧拉方法离散化Lorenz系统"""
    t = np.zeros(num_steps + 1)
    x = np.zeros((num_steps + 1, 3))
    x[0] = x0

    for i in range(num_steps):
        t[i+1] = t[i] + dt
        k1 = lorenz(t[i], x[i], sigma, rho, beta)
        k2 = lorenz(t[i] + dt/2, x[i] + (dt/2) * k1, sigma, rho, beta)
        x[i+1] = x[i] + dt * k2

    return t, x

# 参数设置
sigma = 10
rho = 28
beta = 8/3
x0 = np.array([0.1, 0, 0])
dt = 0.01
num_steps = 10000

# 使用改进的欧拉方法离散化Lorenz系统
t, x = improved_euler_lorenz(sigma, rho, beta, x0, dt, num_steps)

# 绘制轨迹
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot(x[:, 0], x[:, 1], x[:, 2], lw=0.5)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
ax.set_title('Lorenz System - Improved Euler Discretization')
plt.show()

3.四阶龙格-库塔法（Fourth-Order Runge-Kutta method）python代码

import numpy as np
import matplotlib.pyplot as plt

def lorenz(t, x, sigma, rho, beta):
    """Lorenz 系统的微分方程"""
    dx = sigma * (x[1] - x[0])
    dy = x[0] * (rho - x[2]) - x[1]
    dz = x[0] * x[1] - beta * x[2]
    return np.array([dx, dy, dz])

def runge_kutta_lorenz(sigma, rho, beta, x0, dt, num_steps):
    """使用四阶龙格-库塔法离散化 Lorenz 系统"""
    t = np.zeros(num_steps + 1)
    x = np.zeros((num_steps + 1, 3))
    x[0] = x0

    for i in range(num_steps):
        t[i+1] = t[i] + dt

        k1 = lorenz(t[i], x[i], sigma, rho, beta)
        k2 = lorenz(t[i] + dt/2, x[i] + (dt/2) * k1, sigma, rho, beta)
        k3 = lorenz(t[i] + dt/2, x[i] + (dt/2) * k2, sigma, rho, beta)
        k4 = lorenz(t[i] + dt, x[i] + dt * k3, sigma, rho, beta)

        x[i+1] = x[i] + (dt/6) * (k1 + 2*k2 + 2*k3 + k4)

    return t, x

# 参数设置
sigma = 10
rho = 28
beta = 8/3
x0 = np.array([0.1, 0, 0])
dt = 0.01
num_steps = 10000

# 使用四阶龙格-库塔法离散化 Lorenz 系统
t, x = runge_kutta_lorenz(sigma, rho, beta, x0, dt, num_steps)

# 绘制轨迹
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot(x[:, 0], x[:, 1], x[:, 2], lw=0.5)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
ax.set_title('Lorenz System - runge_kutta_lorenz')
plt.show()