基于拉丁方的混沌图像加密算法设计

Design of a Chaotic Image Encryption Algorithm Based on Latin Squares

摘要: 本文研究了一类二维时变离散时空系统的混沌性，并分析了它的多种伪随机性能。在此基础上，结合拉丁方构造的基本密码系统，设计了一种新的多元流密码算法。将该算法应用于数字图像加密之中，并与二元流密码系统的加密效果进行了对比。实验结果表明新算法具有足够大的密钥空间，并具有良好的加密性能和安全性。

Abstract: This paper studies the chaos of a class of two-dimensional time-varying discrete spatiotemporal systems and analyzes its several pseudorandom properties. On this basis, combined with a basic cipher system constructed by a Latin square, a new multi-bit stream cipher algorithm is designed out and is applied in digital image encryption. Compared with the encryption effect of the single-bit stream cipher algorithm, simulations show that the new algorithm has a large enough key space, good encryption performance and security.

1. 引言

随着信息技术的发展，密码学已经成为防止数据泄露、保护数据安全的重要手段。流密码是密码学中常见的一种密码算法，在信息安全领域中有着广泛应用 [1] [2] [3]。根据现代密码学可知，常见的流密码算法都是以模加法运算所建立的完善保密通信模型为理论基础，这导致了近几十年以来模2加法或2元加法或单比特流密码算法得到广泛研究和应用。但是，常见的二元流密码算法的基本加解密变换或运算过于简单而导致其性能存在一些不足 [4] [5]。近几年中，文献 [2] 将模加法流密码的理论模型推广为更一般的基于拉丁方变换的理论模型，将会极大促进流密码算法的理论及其应用的研究。文献 [2] 明确指出流密码算法设计分为两个部分：基本系统和密钥流生成器设计(或主密钥空间设计)，且任意阶拉丁方都可用于设计基本系统，并有可能实现理想安全的完善保密通信。不难发现，二元加法流密码的基本加解密变换本质上是由2阶拉丁方决定的。由于拉丁方的阶数决定了基本系统的复杂度和加解密难度，因而寻求高于2阶或1比特的拉丁方基本密码系统将有利于设计出安全性更高和性能更好的密码算法。然而，当前利用高阶或多比特拉丁方来设计基本密码系统和流密码算法的研究还未充分发展，有许多问题都值得进一步研究。

另一方面，密钥流序列也会极大地影响流密码算法的安全性和性能，其中的一个关键问题是如何设计出随机性优良的密钥流生成器。一个良好的密钥流生成器应满足输出的密钥流具有较强的不可预测性以及抵抗统计分析的能力。当前，利用离散混沌系统来构造密钥流序列及其相应的流密码算法受到了广泛关注和研究。在混沌序列密码算法研究中，根据各类文献发现，离散混沌系统的构造类型多样且类随机性能优良，但对多维时变离散混沌时空系统的研究相对较少。因此，本文拟对一类二维时变离散时空系统的混沌性进行分析或研究其在流密码算法中的应用，可参见图1。

Figure 1. Schematic diagram of key flow generator based on discrete spatiotemporal system

图1. 基于二维离散时空系统的密钥流生成器结构示意图

设Z为全体整数集， $t\in Z$， $N_t=\left\{t,t+1,\cdots \right\}$ 为单边整数集，I为有界实数集，且

$\{\begin{array}{l}{x}_{m+1,n}=f\left(m,y_{m,n},x_{m,n+1}\right)\\ y_{m+1,n}=g\left(m,x_{m,n},y_{m,n+1}\right)\end{array}$ (1-1)

其中， $m,n\in N_0$， $f:N_0\times I^3\to I$ 和 $g:N_0\times I^3\to I$ 是两个三元函数，称 $\left(f,g\right)$ 为系统(1-1)的系统函数。不难验证：给定两个序列 $\varphi ={\left\{{\varphi }_{0,n}\right\}}_{n=0}^{\infty }$ 和 $\phi ={\left\{{\phi }_{0,n}\right\}}_{n=0}^{\infty }$，存在一个二维离散时空序列 $z={\left\{\left(x_{m,n},y_{m,n}\right)\right\}}_{m,n=0}^{\infty }$ 满足(1-1)，且 $x_{0,n}={\varphi }_{0,n}$ 和 $y_{0,n}={\phi }_{0,n}$， $n=0,1,2,\cdots$。称z为系统(1-1)初值为 $\left(\varphi ,\phi \right)$ 的一个解。为了方便，下面将列向量 ${\left(a_0,a_1,\cdots \right)}^{\text{T}}$ 与行向量 $\left(a_0,a_1,\cdots \right)$ 不加区别。

记

$I_2^{\infty }=\left\{{\left\{{\left(a_n,b_n\right)}^{\text{T}}\right\}}_{n=0}^{\infty }=\left(\begin{array}{c}{a}_0{a}_1\cdots a_n\cdots \\ b_0{b}_1\cdots b_n\cdots \end{array}\right)|a_i,b_i\in I,i=0,1,2,\cdots \right\}$ (1-2)

不难验证，在 $I_2^{\infty }$ 上定义一个映射 $d_1:I_2^{\infty }\times I_2^{\infty }\to I_2^{\infty }$ 可得到一个度量空间 $\left(I_2^{\infty },d_1\right)$ ：

$d_1\left(z_1,z_2\right)={\displaystyle \underset{n=0}{\overset{\infty }{\sum }}\frac{\left|x_{1,n}-x_{2,n}\right|+\left|y_{1,n}-y_{2,n}\right|}{2^n}}$, $z_i={\left\{\left(x_{i,n},y_{i,n}\right)\right\}}_{n=0}^{\infty }$, $i=1,2$. (1-3)

设 $z={\left\{\left(x_{m,n},y_{m,n}\right)\right\}}_{m,n=0}^{\infty }$ 是系统(1-1)的一个解， $x_{0,n},y_{0,n}\in I$， $n\in N_0$，且

$z_m={\left\{\left(x_{m,n},y_{m,n}\right)\right\}}_{n=0}^{\infty }$, $m=0,1,2,\cdots$ (1-4)

则系统(1-1)等价于下式(1-5)中的无穷维离散系统：

$z_{m+1}={\left\{\left(x_{m+1,n},y_{m+1,n}\right)=\left(f\left(m,y_{m,n},x_{m,n+1}\right),g\left(m,x_{m,n},y_{m,n+1}\right)\right)\right\}}_{n=0}^{\infty }=G_{m+1}\left(z_m\right)$,(1-5)

其中， $m=0,1,2,\cdots$， $G_1,G_2,\cdots$ 是由 $\left(f,g\right)$ 所决定的 $I_2^{\infty }$ 上的一列映射。称(1-5)或 $G={\left\{G_m\right\}}_{m=0}^{\infty }$ 是由系统(1-1)或 $\left(f,g\right)$ 所导出的离散系统或映射列。

2. 拉丁方的构造方法

下面先给出一种拉丁方的构造方法，以便为基本系统设计做好准备。

定理2.1设n是一个自然数。对任意 $m=m_1{m}_2\cdots m_{2n}\in Z_2^{2n}$ 和 $k=k_1{k}_2\cdots k_{2n}\in Z_2^{2n}$，记

$\begin{array}{c}{c}_{mk}=c_1{c}_2{c}_3{c}_4\cdots c_{2n-1}{c}_{2n}\\ =E\left(m,k\right)=E\left(m_1{m}_2{m}_3{m}_4\cdots m_{2n-1}{m}_{2n},k_1{k}_2{k}_3{k}_4\cdots k_{2n-1}{k}_{2n}\right)\\ =[m_1\left({\bar{m}}_1\oplus m_2\right)m_3\left({\bar{m}}_3\oplus m_4\right)\cdots m_{2n-1}\left({\bar{m}}_{2n-1}\oplus m_{2n}\right)\\ +k_1\left({\bar{k}}_1\oplus k_2\right)k_3\left({\bar{k}}_3\oplus k_4\right)\cdots k_{2n-1}\left({\bar{k}}_{2n-1}\oplus k_{2n}\right)]\mathrm{mod}{2}^{2n}\\ =\left[W\left(m\right)+W\left(k\right)\right]\mathrm{mod}{2}^{2n}\end{array}$ (2-1)

则 $L={\left(c_{mk}\right)}_{2n\times 2n}$ 是一个2n阶拉丁方，其中， $\oplus$ 是异或运算， $\bar{a}=1-a$ 表示比特a的补运算，并将 $Z_{2^{2n}}$ 与 $Z_2^{2n}$ 对应的数不加区别，且对任意 $u=u_1{u}_2\cdots u_{2n-1}{u}_{2n}\in Z_2^{2n}$，映射 $W:Z_2^{2n}\to Z_2^{2n}$ 定义为

$W\left(u\right)=u_1\left({\bar{u}}_1\oplus u_2\right)u_3\left({\bar{u}}_3\oplus u_4\right)\cdots u_{2n-1}\left({\bar{u}}_{2n-1}\oplus u_{2n}\right)\in Z_2^{2n}$.

证明：对任一固定的 $m\in Z_{2^{2n}}$，下面证明当k依次取遍 $Z_{2^{2n}}=\left\{0,1,2,\cdots ,2^{2n}-1\right\}$ 中的每一个整数时，L第m行的所有元素 $c_{m1},c_{m2},\cdots ,c_{m,2n-1},c_{m,2n}$ 将是互不相同的整数，即为 $\left\{0,1,2,\cdots ,2^{2n}-1\right\}$。否则，将存在两个不同的整数 $k=k_1{k}_2\cdots k_{2n}\in Z_2^{2n}$ 和 $\tilde{k}={\tilde{k}}_1{\tilde{k}}_2\cdots {\tilde{k}}_{2n}\in Z_2^{2n}$，使得

$c_{mk}=E\left(m,k\right)=E\left(m,\tilde{k}\right)=c_{m\tilde{k}}$,

即 $m+k_1\left({\bar{k}}_1\oplus k_2\right)\cdots k_{2n-1}\left({\bar{k}}_{2n-1}\oplus k_{2n}\right)=m+{\tilde{k}}_1\left({\overline{\stackrel{˜}{k}}}_1\oplus {\tilde{k}}_2\right)\cdots {\tilde{k}}_{2n-1}\left({\overline{\stackrel{˜}{k}}}_{2n-1}\oplus {\tilde{k}}_{2n}\right)\mathrm{mod}{2}^{2n}$。因此，存在一个整数r，使得 $k_1\left({\bar{k}}_1\oplus k_2\right)\cdots k_{2n-1}\left({\bar{k}}_{2n-1}\oplus k_{2n}\right)={\tilde{k}}_1\left({\overline{\stackrel{˜}{k}}}_1\oplus {\tilde{k}}_2\right)\cdots {\tilde{k}}_{2n-1}\left({\overline{\stackrel{˜}{k}}}_{2n-1}\oplus {\tilde{k}}_{2n}\right)+r{2}^{2n}$。由于

$0\le k_1\left({\bar{k}}_1\oplus k_2\right)\cdots k_{2n-1}\left({\bar{k}}_{2n-1}\oplus k_{2n}\right)<2^{2n}$ 和 $0\le {\tilde{k}}_1\left({\overline{\stackrel{˜}{k}}}_1\oplus {\tilde{k}}_2\right)\cdots {\tilde{k}}_{2n-1}\left({\overline{\stackrel{˜}{k}}}_{2n-1}\oplus {\tilde{k}}_{2n}\right)<2^{2n}$,

因而 $r=0$ 和 $k_1\left({\bar{k}}_1\oplus k_2\right)\cdots k_{2n-1}\left({\bar{k}}_{2n-1}\oplus k_{2n}\right)={\tilde{k}}_1\left({\overline{\stackrel{˜}{k}}}_1\oplus {\tilde{k}}_2\right)\cdots {\tilde{k}}_{2n-1}\left({\overline{\stackrel{˜}{k}}}_{2n-1}\oplus {\tilde{k}}_{2n}\right)$。于是，依次可推出

$k_1={\tilde{k}}_1$, $\left({\bar{k}}_1\oplus k_2\right)=\left({\overline{\stackrel{˜}{k}}}_1\oplus {\tilde{k}}_2\right)$, $\cdots$, $k_{2n-1}={\tilde{k}}_{2n-1}$, $\left({\bar{k}}_{2n-1}\oplus k_{2n}\right)=\left({\overline{\stackrel{˜}{k}}}_{2n-1}\oplus {\tilde{k}}_{2n}\right)$.

这样就有 $k_1={\tilde{k}}_1$， $k_2={\tilde{k}}_2$， $\cdots$， $k_{2n-1}={\tilde{k}}_{2n-1}$， $k_{2n}={\tilde{k}}_{2n}$，即 $k=\tilde{k}$，矛盾！由此证明了方阵L的每一行取遍 $Z_{2^{2n}}=\left\{0,1,2,\cdots ,2^{2n}-1\right\}$ 中的每一个整数。类似可证明L的每一列也取遍 $Z_{2^{2n}}$ 中的每一个整数。因

此，L是一个2n阶拉丁方。证毕！

由定理2.1可知，当 $n=4$ 时可得到如下2⁸阶拉丁方L：

$L=\left[\begin{array}{ccc}\begin{array}{ccc}84& 83& 82\\ 83& 82& 81\\ 82& 81& 80\end{array}& \cdots & \begin{array}{ccc}167& 166& 169\\ 166& 165& 168\\ 165& 164& 167\end{array}\\ ⋮& \ddots & ⋮\\ \begin{array}{ccc}167& 166& 165\\ 166& 165& 164\\ 169& 168& 167\end{array}& \cdots & \begin{array}{ccc}250& 249& 252\\ 249& 248& 251\\ 252& 251& 254\end{array}\end{array}\right]$, $H=\left[\begin{array}{ccc}\begin{array}{ccc}0& 128& 64\\ 1& 129& 65\\ 2& 130& 66\end{array}& \cdots & \begin{array}{ccc}191& 127& 255\\ 190& 126& 254\\ 189& 125& 253\end{array}\\ ⋮& \ddots & ⋮\\ \begin{array}{ccc}253& 125& 189\\ 254& 126& 190\\ 255& 127& 191\end{array}& \cdots & \begin{array}{ccc}66& 130& 2\\ 65& 129& 1\\ 64& 128& 0\end{array}\end{array}\right]$

其中，方阵是文献 [6] 中所构造的2⁸阶拉丁方，且L比H 的构造方法要更复杂一些。本文将基于L 所设计的基本密码系统设为 $\left(M,K,C,E,D\right)$，其中，L的第k行是将基本明文空间 $M=\left\{0,1,\cdots ,255\right\}$ 依次作加密变换所得到的基本密文单元，即对任一 $m\in M$，按照式(2-1)中当 $n=4$ 时的公式加密可得到密文 $c=E\left(m,k\right)$，且基本解密变换如下：对任意 $u=u_1{u}_2\cdots u_7{u}_8\in Z_2^8$，

$m_{mk}=D\left(m,k\right)=W\left[\left(c-W\left(k\right)\right)\mathrm{mod}256\right]$, $W\left(u\right)=u_1\left({\bar{u}}_1\oplus u_2\right)\cdots u_7\left({\bar{u}}_7\oplus u_8\right)\in Z_2^8$.

3. 新混沌离散系统及其伪随机性分析

设 $I=\left[0,1\right)$，并定义两个函数 $f:N_0\times I^3\to I$ 和 $g:N_0\times I^3\to I$ ：

$f\left(m,y,x\right)=\langle s_m{y}^2+a_{m}x\rangle$ 和 $g\left(m,x,y\right)=\langle t_m{x}^2+b_{m}y\rangle$， (3-1)

其中， ${\left\{a_m\right\}}_{m=0}^{\infty }$， ${\left\{b_m\right\}}_{m=0}^{\infty }$， ${\left\{s_m\right\}}_{m=0}^{\infty }$ 和 ${\left\{t_m\right\}}_{m=0}^{\infty }$ 都是周期数列，即存在正整数p，使得 $a_m=a_{m+p}$， $b_m=b_{m+p}$， $s_m=s_{m+p}$ 和 $t_m=t_{m+p}$，对一切 $m\in \left\{0,1,2,\cdots \right\}$，且 $\langle u\rangle$ 表示实数u的小数。

不难发现，上述映射 $\left(f,g\right)$ 可产生如下的二维时变离散时空系统

$\{\begin{array}{l}{x}_{m+1,n}=f\left(m,y_{m,n},x_{m,n+1}\right)=\langle s_m{y}_{m,n}^2+a_m{x}_{m,n+1}\rangle \\ y_{m+1,n}=g\left(m,x_{m,n},y_{m,n+1}\right)=\langle t_m{x}_{m,n}^2+b_m{y}_{m,n+1}\rangle \end{array}$ (3-2)

其中， $x_{m,n},y_{m,n}\in I=\left[0,1\right)$，对任一 $m,n=0,1,2,\cdots$。参照现有文献可以发现，该二维时变离散时空系统(3-2)

的Devaney混沌性还没有被研究过。参照系统(1-1)和系统(1-5)之间的关系可知，系统(3-2)会等价于如下

离散系统：对 $z={\left\{\left(x_{m,n},y_{m,n}\right)\right\}}_{m,n=0}^{\infty }\in I_2^{\infty }$，

$z_{m+1}={\left\{\left(x_{m+1,n},y_{m+1,n}\right)=\left(\langle s_m{y}_{m,n}^2+a_m{x}_{m,n+1}\rangle ,\langle t_m{x}_{m,n}^2+b_m{y}_{m,n+1}\rangle \right)\right\}}_{n=0}^{\infty }=G_{m+1}\left(z_m\right)$,(3-3)

其中， $G_{m+1}:I_2^{\infty }\to I_2^{\infty }$ 是由f和g决定的映射，即对任一 $\beta ={\left\{\left(u_m,v_m\right)\right\}}_{m=0}^{\infty }\in I_2^{\infty }$，都有

$G_{m+1}\left(\beta \right)={\left\{\langle s_m{v}_n^2+a_m{u}_{n+1}\rangle ,\langle t_m{u}_n^2+b_m{v}_{n+1}\rangle \right\}}_{n=0}^{\infty }$, $m=1,2,\cdots$. (3-4)

参照文献 [7] [8] [9] [10]，易得如下引理及推论。

引理3.1式(3-3)中映射列 $G={\left\{G_m\right\}}_{m=1}^{\infty }$ 是周期为p的周期序列，即 $G_m=G_{m+p}$， $m=1,2,\cdots$。而且，对一切 $m,s,t\in N_1$，都有

$G_{s+mp-1}\circ G_{s+mp-2}\circ \cdots \circ G_s=G_{t+mp-1}\circ G_{t+mp-2}\circ \cdots \circ G_t$.

记式(3-3)定义的映射列 $G={\left\{G_m\right\}}_{m=1}^{\infty }$ 所确定的复合映射为

$H_m=G_m\circ G_{m-1}\circ \cdots \circ G_1$, $m=1,2,\cdots$. (3-5)

引理3.2对任一给定的常数 $c_m\ge 1$， $h,r\in I$ 和 $m=1,2,\cdots$，存在 $u\in I$，使得 $\langle h+c_{m}u\rangle =r$。

定义3.1若式(1-1)的系统函数 $\left(f,g\right)$ 所确定的映射列 $G={\left\{G_m\right\}}_{m=1}^{\infty }$ 或系统(1-5)在度量空间 $\left(I_2^{\infty },d_1\right)$ 上具有传递性、周期点的稠密性和初值敏感依赖性，则称 $G={\left\{G_m\right\}}_{m=1}^{\infty }$ 或系统(1-5)是Devaney混沌的，也称与系统(1-5)相应的等价系统(1-1)在 $\left(I_2^{\infty },d_1\right)$ 上是Devaney混沌的。

定理3.1在上述条件下，离散系统(3-2)是 $\left(I_2^{\infty },d_1\right)$ 上的Devaney混沌系统。

证明：由定义3.1可知，只需证明系统(3-3)在 $\left(I_2^{\infty },d_1\right)$ 上是Devaney混沌的，即该系统具有周期点的稠密性、传递性和初值的敏感依赖性。首先将证明该系统具有周期点的稠密性。

对于给定的任一点 $\alpha ={\left\{\left(u_n,v_n\right)\right\}}_{n=0}^{\infty }\in I_2^{\infty }$ 和 $\alpha$ 的任一邻域 $U\subseteq I_2^{\infty }$，都存在 $\delta >0$，使得

$B_{\delta }\left(\alpha \right)={\left\{x={\left\{\left(x_n,y_n\right)\right\}}_{n=0}^{\infty }\in I_2^{\infty }|d_1\left(x,\alpha \right)<\delta \right\}}_{n=0}^{\infty }\subseteq U$.

根据式(1-3)中 $d_1$ 的定义知，存在整数 $M\in \left\{1,2,\cdots \right\}$，使得

${\left\{x={\left\{\left(x_n,y_n\right)\right\}}_{n=0}^{\infty }\in I_2^{\infty }|x_i=u_i,y_i=v_i;i=0,1,\cdots ,M-1;\forall x_j,y_j\in I,j\ge M\right\}}_{n=0}^{\infty }\subseteq B_{\delta }\left(\alpha \right)$ (3-6)

对任意一点 $x={\left\{\left(x_n,y_n\right)\right\}}_{n=0}^{\infty }\in I_2^{\infty }$，记

$H_1\left(x\right)=G_1\left(x\right)={\left\{\left(x_n^{\left(1\right)},y_n^{\left(1\right)}\right)=\left(\langle s_0{y}_n^2+a_0{x}_{n+1}\rangle ,\langle t_0{x}_n^2+b_0{y}_{n+1}\rangle \right)\right\}}_{n=0}^{\infty }=x^{\left(1\right)}$,

$H_2\left(x\right)=G_2\left(x^{\left(1\right)}\right)={\left\{\left(x_n^{\left(2\right)},y_n^{\left(2\right)}\right)=\left(\langle s_1{\left(y_n^{\left(1\right)}\right)}^2+a_1{x}_{n+1}^{\left(1\right)}\rangle ,\langle t_1{\left(x_n^{\left(1\right)}\right)}^2+b_1{y}_{n+1}^{\left(1\right)}\rangle \right)\right\}}_{n=0}^{\infty }=x^{\left(2\right)}$,

其中， $x_n^{\left(2\right)}=\langle s_1{\left(t_0{x}_n^2+b_0{y}_{n+1}\right)}^2+a_1{s}_0{y}_{n+1}^2+a_1{a}_0{x}_{n+2}\rangle =\langle f_1\left(x_n,x_{n+1},y_n,y_{n+1}\right)+a_1{a}_0{x}_{n+2}\rangle$，且

$f_1\left(x_n,x_{n+1},y_n,y_{n+1}\right)=\langle s_1{\left(t_0{x}_n^2+b_0{y}_{n+1}\right)}^2+a_1{s}_0{y}_{n+1}^2\rangle$,

以及 $y_n^{\left(2\right)}=\langle g_1\left(y_n,y_{n+1},x_n,x_{n+1}\right)+b_1{b}_0{y}_{n+2}\rangle$ 和 $g_1\left(y_n,y_{n+1},x_n,x_{n+1}\right)=\langle t_1{\left(s_0{y}_n^2+a_0{x}_{n+1}\right)}^2+b_1{t}_0{x}_{n+1}^2\rangle$。

由递推法可知，对任意 $m=1,2,\cdots$ 和 $x={\left\{\left(x_n,y_n\right)\right\}}_{n=0}^{\infty }\in I_2^{\infty }$，有

$H_m\left(x\right)=G_m\left(x^{\left(m-1\right)}\right)=G_m\circ \cdots \circ G_1\left(x\right)={\left\{\left(x_n^{\left(m\right)},y_n^{\left(m\right)}\right)\right\}}_{n=0}^{\infty }=x^{\left(m\right)}$,

其中， $f_m\left(x_n,\cdots ,x_{n+m},y_n,\cdots ,y_{n+m}\right)=f_m\left({\lambda }_n\right)$ 和 $g_m\left(y_n,\cdots ,y_{n+m},x_n,\cdots ,x_{n+m}\right)=g_m\left({\rho }_n\right)$ 是两个多元函数， ${\lambda }_n=\left(x_n,\cdots ,x_{n+m},y_n,\cdots ,y_{n+m}\right)$ 和 ${\rho }_n=\left(y_n,\cdots ,y_{n+m},x_n,\cdots ,x_{n+m}\right)$，且

$x^{\left(m+1\right)}={\left\{\left(\langle f_m\left({\lambda }_n\right)+a_m\cdots a_0{x}_{n+m+1}\rangle ,\langle g_m\left({\rho }_n\right)+b_m\cdots b_0{y}_{n+m+1}\rangle \right)\right\}}_{n=0}^{\infty }$. (3-7)

由已知条件可知，对任一 $m=1,2,\cdots$，都有 $a_m\cdots a_0\ge 1$ 和 $b_m\cdots b_0\ge 1$。由引理3.2，对任一固定的 $\eta \in I$， $\lambda ,\rho \in I^{2\left(m+1\right)}$ 和 $m=1,2,\cdots$，都存在实数 $c=c_m,e=e_m\in I$，使得

$\langle f_m\left(\lambda \right)+a_m\cdots a_{0}c\rangle =\eta$ 和 $\langle g_m\left(\rho \right)+b_m\cdots b_{0}e\rangle =\eta$. (3-8)

对于 $\alpha ={\left\{\left(u_n,v_n\right)\right\}}_{n=0}^{\infty }\in U$，由式(3-6)、(3-7)、(3-8)，存在一点 $\beta ={\left\{\left({\sigma }_n,{\tau }_n\right)\right\}}_{n=0}^{\infty }\in U$，使得 ${\sigma }_n=u_n,{\tau }_n=v_n$， $n=0,1,\cdots ,M-1$，且 ${\sigma }_M,{\tau }_M,{\sigma }_{M+1},{\tau }_{M+1},\cdots$ 依次满足：

$\langle f_{M-1}\left({\lambda }_n\right)+a_M\cdots a_0{\sigma }_{n+M+1}\rangle ={\sigma }_n$, $\langle g_{M-1}\left({\rho }_n\right)+b_M\cdots b_0{\tau }_{n+M+1}\rangle ={\tau }_n$, $n=0,1,2,\cdots$.

因此， $H_M\left(\beta \right)=G_M\circ \cdots \circ G_1\left(\beta \right)=\beta \in V$ 和 $\beta \in U$，即 $\beta$ 是 $G={\left\{G_m\right\}}_{m=1}^{\infty }$ 的一个周期点。于是，系统(3-3)在 $\left(I_2^{\infty },d_1\right)$ 上具有周期点的稠密性。同样地，可利用与上述类似的方法证明系统(3-3)还具有传递性和初值敏感依赖性。故系统(3-3)在 $\left(I_2^{\infty },d_1\right)$ 上是Devaney混沌的。证毕！

例1：设二维时变离散时空系统为

$\{\begin{array}{l}{x}_{m+1,n}=f\left(m,y_{m,n},x_{m,n+1}\right)=\langle 4y_{m,n}^2+a_m{x}_{m,n+1}\rangle \\ y_{m+1,n}=g\left(m,x_{m,n},y_{m,n+1}\right)=\langle 4x_{m,n}^2+b_m{y}_{m,n+1}\rangle \end{array}$ (3-9)

其中， $a_m=3+{\left(-1\right)}^m$ 和 $b_m=3+{\left(-1\right)}^{m+1}$，对任意 $m=0,1,2,\cdots$。

基于现有文献无法断定上述系统是否具有Devaney混沌性。但由于该系统是(3-2)的特殊情形，结合定理3.1可知，该系统是 $\left(I_2^{\infty },d_1\right)$ 上的Devaney混沌系统。

下面通过SP800-22Rvlla标准系统检验例1系统混沌解序列的伪随机性能，其中，SP800-22Rvlla标准由美国国家标准技术研究所提出的专门针对随机数或伪随机数发生器产生的二进制随机序列的统计测试法。对长度为1×10⁶的混沌序列进行伪随机性检验，设置比特序列长度为10⁶，显著性水平为0.01，当测试p值大于显著性水平时，检验通过，测试结果见表1，且混乱性可见图2。

Table 1. Test results of SP800-22Rvlla

表1. SP800-22Rvlla测试结果

Figure 2. Confusion and autocorrelation of solutions

图2. 解的混乱性和自相关性

4. 基于拉丁方的多比特流密码算法

4.1. 算法描述

先将基本系统设计为上述8比特的拉丁方基本系统 $\left(M=C=K=Z_2^8,E,D\right)$，再结合定理3.1，将完整的加解密变换具体设计如下：

1) 取一幅数字图像作为明文m (如Lena图像)，并将m转换为二元序列 $m=m_1{m}_2\cdots$， $m_j\in Z_2$，将该二元序列按每8比特进行分组，将分组后得到的明文单元序列设为 $\tilde{m}={\tilde{m}}_1{\tilde{m}}_2\cdots$，其中 ${\tilde{m}}_1=m_{11}{m}_{12}\cdots m_{18}\in Z_2^8$，等等。必要时可对m的最后一个分组填充而组成一个8比特分组。

2) 加密变换E： ${\tilde{c}}_j=E\left({\tilde{k}}_j,{\tilde{m}}_j\right)$， $j=1,2,\cdots$，其中，将系统(3-9)的一个解序列作为256元密钥流序列 $\tilde{k}={\tilde{k}}_1{\tilde{k}}_2\cdots$，可得到256元密文单元序列 $\tilde{c}={\tilde{c}}_1{\tilde{c}}_2\cdots$。

3) 解密变换D： ${\tilde{m}}_j=D\left({\tilde{k}}_j,{\tilde{c}}_j\right)$， $j=1,2,\cdots$，可得到明文256元序列 $\tilde{m}={\tilde{m}}_1{\tilde{m}}_2\cdots$。然后将每个明文单元 ${\tilde{m}}_j$ 表示为8比特明文分组就得到解密后的原始二元明文序列 $m=m_1{m}_2\cdots$。最后可恢复出原数字图像m。

4.2. 实验结果及分析

将上述流密码算法用于数字图像加解密，从直方图、相邻像素相关性、密钥敏感性以及信息熵四个方面将它与单比特流密码系统进行比较，仿真效果一一介绍如下。

4.2.1. 直方图分析

直方图可以直观地反映一幅图像的像素值分布特征。如果密文图像的像素值分布越均匀，那么其抵抗统计攻击的能力越强，算法的安全性越高。在Matlab中，两种流密码算法仿真的加解密效果图及直方图如下图3所示：

Figure 3. Simulation renderings of two algorithms

图3. 算法仿真效果图

模2加法或单比特流密码算法密文图像的像素分布区间为[0, 255]，分布近似呈现周期性或不具有均匀分布，且周期内分布杂乱；密文图像(d)的像素在区间[0, 255]接近均匀分布。因此，与单比特算法密文图像相比，本加密算法加密后的统计特性更好，因而抵抗统计攻击的能力会更强。

4.2.2. 相邻像素相关性分析

像素相关性表示一幅图像中相邻像素之间的关联程度。当数字图像的原图或加密图像相邻像素点之间有较高相关性时，可利用相邻像素预测该点的像素值，进而有可能获取图像信息。可靠的加密算法应保证加密之后的密文图像具有足够弱的像素相关性。衡量图像相邻像素间的相关性主要有定性和定量两种方法：定性的方法是指从水平、垂直以及对角三个方向上直观地显示像素相关性分布图；定量的方法可通过计算相邻像素的相关系数，可利用该量化指标进行衡量。

从Lena原图像和两种加密算法得到的密文图像中随机选取4000相邻像素对，并从水平、垂直以及对角三个方向绘制像素分布图，结果见图4。

Figure 4. Pixel correlation distribution

图4. 像素相关性分布图

从像素间相关性分布图来看，原图像的像素点在水平、垂直和对角方向上都分布在对角线附近，表明这些方向上相邻像素间的相关性较强；模2算法的密文图像像素点分布在250 × 250整个平面区域内，但像素点在45˚方向上表现出较强的线性关系；本文算法的密文图像像素点在整个区域内均匀分布，相邻像素点的相关性大大减弱，可较好地掩盖原始图像的相关特征。

为定量测试密文图像的相邻像素相关性，在图像各个方向上选取4000个相邻像素点，根据下面的式(4-1)和(4-2)计算100组相邻像素的相关系数，其均值结果如表2所示。

$r_{xy}=\frac{cov\left(x,y\right)}{\sqrt{D\left(x\right)}\times \sqrt{D\left(y\right)}}$ (4-1)

$\{\begin{array}{l}cov\left(x,y\right)=\frac{1}{N}{\displaystyle {\sum }_{i=1}^N\left(x_i-E\left(x\right)\right)\left(y_i-E\left(y\right)\right)}\\ E\left(x\right)=\frac{1}{N}{\displaystyle {\sum }_{i=1}^N{x}_i}\\ D\left(x\right)=\frac{1}{N}{\displaystyle {\sum }_{i=1}^N{\left(x_i-E\left(x\right)\right)}^2}\end{array}$ (4-2)

Table 2. Correlation simulation data

表2. 相关性仿真数据

对比单比特流密码算法得到的密文图像，本文加密算法得到的密文图像在各个方向上的相邻像素相关系数都更接近0，几乎不存在相关性。

4.2.3. 密钥敏感性

一个好的加密算法应对密钥的变化极其敏感。当密钥仅发生微小变化时，也无法还原出明文图像。为了评估本文所提出的加密算法的密钥敏感度，利用单一变量法，每次仅改变密钥 $x_0$ 、 $y_0$ 、 $u_1$ 和 $u_2$ 其中一个密钥值，改动幅度为1 × 10⁻¹⁴，得到key1、key2、key3和key4四个新密钥。当攻击方获取原密钥加密的原加密图像后，分别利用四个微小改变的新密钥对原加密图像进行解密也无法恢复出明文图像，且解密图像中无有效信息。因此，本文加密系统具有良好的密钥敏感性，可保证密钥空间具有足够多的密钥量(图5)。

(a) 原始图像 (b) 密文图像 (c) 原密钥的解密图像(d) key1的解密图像 (e) key2的解密图像 (f) key3的解密图像 (g) key4的解密图像

Figure 5. Key sensitivity test

图5. 密钥敏感性测试

4.2.4. 信息熵分析

信息熵为图像所含信息量，可以反映图像的随机性或不确定性。一幅图像的不确定性越强，信息量越大，且信息熵越大。根据图像信息熵的计算公式(4-3)和最大熵原理 [11], 当图像灰度值范围为[0, 255]时，理论最大信息熵达到8。明文图像和密文图像的信息熵结果如表3所示。

$H\left(x\right)=-{\displaystyle {\sum }_{i=1}^{256}P\left(x_i\right){\log }_{2}P\left(x_i\right)}$ (4-3)

其中， $x_i=i$， $P\left(x_i\right)$ 为像素 $x_i$ 出现的概率。

Table 3. Comparison of information entropy

表3. 信息熵对比结果

明文图像和单比特流加密算法得到的密文图像信息熵均为7.4442，而本文加密算法得到的密文图像信息熵为7.9971，接近理想值8。因此，本文加密算法可以更好地抵御信息熵的攻击，具有更好的加密性能。

5. 小结

本文构建了一类新的时变离散时空混沌系统，基于该系统和高阶拉丁方构造了与常见模2加法或1比特流密码算法明显不同的一种多比特流密码算法。通过对数字图像加解密的多种效果分析与对比，验证了本文提出的算法具有可行性和良好的性能与安全性，具有一定的实际应用参考价值。不过，也要指出：复杂度适中的拉丁方及其变换组的流密码系统性构造方法还值得今后进一步研究。