写在前面

众所周知，代打密码学实验只有一次和无数次（逃）。一开始其实蛮期待下一个密码学实验的，大抵是由于上一次做实验 0-3 大数运算实在太痛苦，花了大量时间调试，最优的那个算法还是没有调试出来，博客也没法写，所以导致这次看到 AES 和马上就要到的提交代码 DDL 就有点犹豫。

可惜，因为有同学考公申请迟交，提交代码的期限延长了。说干就干，结果花了一天一夜就完成了所有 AES 实验（包括选做部分）。一眼丁真，评价为比大数运算简单。

实验描述

完整的实验描述可以在实验网站上看到。完整实验代码放到了GitHub上。

AES 是一种分组密码，意味着它将数据分成固定大小的块，然后对每个块应用加密算法。AES 支持的分组大小为 128 位，而密钥长度可以是 128 位、192 位或 256 位。我们这次要实现的是 AES-128，也就是分组大小和密钥长度都是 128 位。为了简便，下面将 32 位称为一个字，128 位称为一个块。

抛开各种工作模式的差异，AES 的核心就是用某个密钥对某个（明文或者处理过的）分组进行加密/解密。加密过程可以描述如下：

  1. 对原始密钥进行拓展，最终生成 11 个轮密钥，每个轮密钥 128 位，依次用于接下来的每轮加密过程。

  2. 第 0 轮，进行轮密钥加，也即将分组与轮密钥逐位异或。

  3. 第 1-9 轮，每轮用上一轮过后的分组，依次进行字节替代、行移位、列混淆和轮密钥加四个步骤。

  4. 第 10 轮，进行字节替代、行移位和轮密钥加，最后得到密文。

解密过程可以从后往前逆着回去，操作都用相应的逆操作。也有另一种等价的解密过程，操作顺序与加密过程一致，只是每一步的操作都用逆操作，且部分轮密钥要额外经过一次列混淆。

以上加密解密过程可以参考来自 NIST 的 AES 标准文档FIPS 197。

AES-128 CBC 模式的实现

以 CBC 模式为例，实现普通的 AES-128，这对应着实验 3-1，以及仓库中的aes.c。

输入输出的处理

与之前的实验相同，本次实验的数据也是二进制给出的。由于明文、密钥等都可以理解为字节流，所以并无端序的区别。但要注意的是，AES 是可以用状态机描述的，其中的每个状态都是一个 4*4 的字节矩阵，且是列优先的。

例如，若输入给出明文：

1
00 01 02 03 04 05 06 07 08 09 0A 0B 0C 0D 0E 0F

则初始状态矩阵应该为：

1
2
3
4
00 04 08 0C
01 05 09 0D
02 06 0A 0E
03 07 0B 0F

同理，轮密钥矩阵的排列方式也是相似的。所以注意一开始读入的时候要处理好，直接用uint8_t[4][4]来存储，并注意是列优先的。（一开始我用的是uint8_t[16]来读入，中间用uint32_t[4]存储状态，还弄成行优先了，每次都要先拆成uint8_t[4][4]再运算，被自己蠢哭，而且直到写这篇文章的时候才意识到这点。鉴定为没上过理论课，并且不认真学前置知识导致的）

笔者写到这里的时候，又花了两个小时去调整了所有的接口/函数参数类型，再一次被自己蠢哭。

总之，最终主函数中处理输入输出的代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
int main(void)
{
#ifndef ONLINE_JUDGE
    FILE *in = fopen("aesi.bin", "rb");
    FILE *out = fopen("aeso.bin", "wb");
#else
#define in stdin
#define out stdout
#endif

    uint8_t mode;
    uint8_t key[4][4];
    uint8_t IV[4][4];
    uint32_t len;
    uint8_t res[4][4];

    fread(&mode, sizeof(uint8_t), 1, in);
    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            fread(&key[j][i], sizeof(uint8_t), 1, in);
    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            fread(&IV[j][i], sizeof(uint8_t), 1, in);
    fread(&len, sizeof(uint32_t), 1, in);

    uint8_t last[4][4];
    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            last[i][j] = IV[i][j];

    init_tables();

    if (mode == 0x01)
    {
        uint32_t i = 0;
        uint8_t data[16];
        uint8_t state[4][4];
        while (i < len)
        {
            if (i + 16 > len)
            {
                fread(data, sizeof(uint8_t), len - i, in);
                padding(data, 16 - (len - i));
            }
            else
                fread(data, sizeof(uint8_t), 16, in);

            for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                    state[j][i] = data[i * 4 + j];

            encrypt(state, key, last, res);
            fwrite(res, sizeof(uint8_t), 16, out);
            i += 16;
        }

        if (i == len)
        {
            padding(data, 16);

            for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                    state[j][i] = data[i * 4 + j];

            encrypt(state, key, last, res);
            fwrite(res, sizeof(uint8_t), 16, out);
        }
    }
    else
    {
        uint32_t i = 0;
        uint8_t data[16];
        uint8_t state[4][4];

        while (i < len)
        {
            fread(data, sizeof(uint8_t), 16, in);
            for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                    state[j][i] = data[i * 4 + j];

            decrypt(state, key, last, res);

            int len_ = 16;
            if (i + 16 == len)
                len_ = 16 - (int)res[3][3];

            fwrite(res, sizeof(uint8_t), len_, out);

            i += 16;
            for (int i = 0; i < 4; i++)
                for (int j = 0; j < 4; j++)
                    last[j][i] = data[i * 4 + j];
        }
    }

    return 0;
}

其中有一些奇怪的变量和函数，例如IV、init_tables()、padding()等，这些我们后面会提到，现在先不用管。

上述代码中，key、IV我们都特别处理成了列优先的形式，整个加密过程中，也是以列优先方式进行的。因为多次调用fread和fwrite的开销较大，所以最开始读入data是用uint8_t[16]，然后再转成uint8_t[4][4]。最后输出res的时候，也是提前转成了行优先的形式，这样就可以只调用一次fwrite。

密钥拓展

伪代码可以参考《密码学原理与实践（第三版）》一书的算法 3.6：

  1. 一共有 11 个轮密钥，也即 44 个字。前四个字W[0]W[0]到W[3]W[3]就是原始密钥。

  2. 计算W[i]W[i]。首先令W[i]=W[i−1]W[i]=W[i-1]，如果i≡0mod  4i\equiv 0 \mod 4，则依次进行字循环、字节替代、轮常量异或，否则不进行任何操作；最后再与W[i−4]W[i-4]异或。

字循环

也即一个字内以字节为单位，左循环移位一位：

1
2
3
4
5
6
7
8
9
void rot_word(uint8_t word[4])
{
    uint8_t temp[4];
    for (int i = 0; i < 4; i++)
        temp[i] = word[(i + 1) % 4];

    for (int i = 0; i < 4; i++)
        word[i] = temp[i];
}

字节替代

也对输入的每个字节，高四位作为行号，低四位作为列号，查表替换：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
const uint8_t SBOX[16][16] = {
    {0x63, 0x7c, 0x77, 0x7b, 0xf2, 0x6b, 0x6f, 0xc5, 0x30, 0x01, 0x67, 0x2b, 0xfe, 0xd7, 0xab, 0x76},
    {0xca, 0x82, 0xc9, 0x7d, 0xfa, 0x59, 0x47, 0xf0, 0xad, 0xd4, 0xa2, 0xaf, 0x9c, 0xa4, 0x72, 0xc0},
    {0xb7, 0xfd, 0x93, 0x26, 0x36, 0x3f, 0xf7, 0xcc, 0x34, 0xa5, 0xe5, 0xf1, 0x71, 0xd8, 0x31, 0x15},
    {0x04, 0xc7, 0x23, 0xc3, 0x18, 0x96, 0x05, 0x9a, 0x07, 0x12, 0x80, 0xe2, 0xeb, 0x27, 0xb2, 0x75},
    {0x09, 0x83, 0x2c, 0x1a, 0x1b, 0x6e, 0x5a, 0xa0, 0x52, 0x3b, 0xd6, 0xb3, 0x29, 0xe3, 0x2f, 0x84},
    {0x53, 0xd1, 0x00, 0xed, 0x20, 0xfc, 0xb1, 0x5b, 0x6a, 0xcb, 0xbe, 0x39, 0x4a, 0x4c, 0x58, 0xcf},
    {0xd0, 0xef, 0xaa, 0xfb, 0x43, 0x4d, 0x33, 0x85, 0x45, 0xf9, 0x02, 0x7f, 0x50, 0x3c, 0x9f, 0xa8},
    {0x51, 0xa3, 0x40, 0x8f, 0x92, 0x9d, 0x38, 0xf5, 0xbc, 0xb6, 0xda, 0x21, 0x10, 0xff, 0xf3, 0xd2},
    {0xcd, 0x0c, 0x13, 0xec, 0x5f, 0x97, 0x44, 0x17, 0xc4, 0xa7, 0x7e, 0x3d, 0x64, 0x5d, 0x19, 0x73},
    {0x60, 0x81, 0x4f, 0xdc, 0x22, 0x2a, 0x90, 0x88, 0x46, 0xee, 0xb8, 0x14, 0xde, 0x5e, 0x0b, 0xdb},
    {0xe0, 0x32, 0x3a, 0x0a, 0x49, 0x06, 0x24, 0x5c, 0xc2, 0xd3, 0xac, 0x62, 0x91, 0x95, 0xe4, 0x79},
    {0xe7, 0xc8, 0x37, 0x6d, 0x8d, 0xd5, 0x4e, 0xa9, 0x6c, 0x56, 0xf4, 0xea, 0x65, 0x7a, 0xae, 0x08},
    {0xba, 0x78, 0x25, 0x2e, 0x1c, 0xa6, 0xb4, 0xc6, 0xe8, 0xdd, 0x74, 0x1f, 0x4b, 0xbd, 0x8b, 0x8a},
    {0x70, 0x3e, 0xb5, 0x66, 0x48, 0x03, 0xf6, 0x0e, 0x61, 0x35, 0x57, 0xb9, 0x86, 0xc1, 0x1d, 0x9e},
    {0xe1, 0xf8, 0x98, 0x11, 0x69, 0xd9, 0x8e, 0x94, 0x9b, 0x1e, 0x87, 0xe9, 0xce, 0x55, 0x28, 0xdf},
    {0x8c, 0xa1, 0x89, 0x0d, 0xbf, 0xe6, 0x42, 0x68, 0x41, 0x99, 0x2d, 0x0f, 0xb0, 0x54, 0xbb, 0x16}};

void sub_byte(uint8_t *byte)
{
    int row = *byte >> 4;
    int col = *byte & 0x0f;
    *byte = SBOX[row][col];
}

void sub_word(uint8_t word[4])
{
    for (int i = 0; i < 4; i++)
        sub_byte(&word[i]);
}

这里的字节替代和稍后加密过程中的字节替代，用的是同一张表SBOX，称为 S 盒。S 盒的构造涉及到有限域F28=Z2[x]/(x8+x4+x3+x+1)\mathbb{F}_{2^8}=\mathbb{Z}_2[x]/(x^8+x^4+x^3+x+1)中的运算。

如果记输入的字为a=(a7,a6,a5,a4,a3,a2,a1,a0)T,ai∈{0,1}\mathbf{a}=(a_7, a_6, a_5, a_4, a_3, a_2, a_1, a_0)^T, a_i\in\{0, 1\}，则经过 利用 S 盒对输入做字节代换等价一个仿射变换：

S(a)=La⋅a−1+μS(\mathbf{a}) = \mathbf{L_a} \cdot \mathbf{a}^{-1} + \mathbf{\mu}

其中：

  • a−1\mathbf{a}^{-1}表示a\mathbf{a}在有限域中的逆元

  • La\mathbf{L_a}是仿射变换矩阵，它就是把下面的过程：

    bi=ai+a(i+4)%8+a(i+5)%8+a(i+6)%8+a(i+7)%8,i=0,1,…,7b_i = a_i + a_{(i+4)\%8} + a_{(i+5)\%8} + a_{(i+6)\%8} + a_{(i+7)\%8}, i=0, 1, \dots, 7

    写成矩阵运算形式的结果：

    La=[1000111111001011111001111111000111110000011110000011110000011111]\mathbf{L_a} = \begin{bmatrix} 1 & 0 & 0 & 0 & 1 & 1 & 1 & 1 \\ 1 & 1 & 0 & 0 & 1 & 0 & 1 & 1 \\ 1 & 1 & 1 & 0 & 0 & 1 & 1 & 1 \\ 1 & 1 & 1 & 1 & 0 & 0 & 0 & 1 \\ 1 & 1 & 1 & 1 & 0 & 0 & 0 & 0 \\ 0 & 1 & 1 & 1 & 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 1 & 1 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 1 & 1 & 1 & 1 \\ \end{bmatrix}

  • μ\mathbf{\mu}是仿射变换的平移项，它是一个常数向量：

    μ=(0,1,1,0,0,0,1,1)T\mathbf{\mu} = (0, 1, 1, 0, 0, 0, 1, 1)^T

所以如果将一个矩阵初始化为：

[0x00000x0001⋯0x000f0x00100x0011⋯0x001f⋮⋮⋱⋮0x00f00x00f1⋯0x00ff]\begin{bmatrix} \text{0x0000} & \text{0x0001} & \cdots & \text{0x000f} \\ \text{0x0010} & \text{0x0011} & \cdots & \text{0x001f} \\ \vdots & \vdots & \ddots & \vdots \\ \text{0x00f0} & \text{0x00f1} & \cdots & \text{0x00ff} \\ \end{bmatrix}

然后：

  1. 将矩阵中的每个元素映射为其在有限域中的逆（0x00 映射为自身）
  2. 对矩阵中的每个元素应用仿射变换

就得到了我们代码中的常量表SBOX。我们这里直接用现成的，没有自己构造。

以上关于 S 盒构造的内容参考了这篇文章。

轮常量异或

有一个轮常量数组，是一系列 32 位的常量。与对应的轮常量逐位异或即可：

1
2
3
4
5
6
7
8
9
10
11
12
13
const uint32_t RCON[11][4] = {{0x00, 0x00, 0x00, 0x00},
    {0x01, 0x00, 0x00, 0x00}, {0x02, 0x00, 0x00, 0x00},
    {0x04, 0x00, 0x00, 0x00}, {0x08, 0x00, 0x00, 0x00},
    {0x10, 0x00, 0x00, 0x00}, {0x20, 0x00, 0x00, 0x00},
    {0x40, 0x00, 0x00, 0x00}, {0x80, 0x00, 0x00, 0x00},
    {0x1b, 0x00, 0x00, 0x00}, {0x36, 0x00, 0x00, 0x00}};


void rcon(uint8_t word[4], int round)
{
    for (int i = 0; i < 4; i++)
        word[i] ^= RCON[round][i];
}

实际用到的只有RCON[1]到RCON[10]。

完整密钥拓展代码

写完以上几个函数后，可以写出完整的拓展代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
void key_expand(uint8_t key[4][4], uint8_t expand_key[11][4][4])
{
    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            expand_key[0][i][j] = key[i][j];

    for (int i = 1; i <= 10; i++)
    {
        for (int j = 0; j < 4; j++)
        {
            uint8_t temp[4];
            int p = i * 4 + j;
            for (int k = 0; k < 4; k++)
                temp[k] = expand_key[(p - 1) / 4][k][(p - 1) % 4];

            if ((i * 4 + j) % 4 == 0)
            {
                rot_word(temp);
                sub_word(temp);
                rcon(temp, i);
            }

            for (int k = 0; k < 4; k++)
                expand_key[i][k][j] = temp[k] ^ expand_key[(p - 4) / 4][k][(p - 4) % 4];
        }
    }
}

本质上还是按照伪代码一个字一个字处理的。务必注意，我们存储的时候是列优先的，所以这里给temp赋值包括将结果写回expand_key数组时，要注意下标。

加密过程

有了轮密钥过程就可以进行密钥拓展了。

填充

由于 CBC 模式要求明文/密文的长度都必须是 16 字节的倍数（与密钥长度一致），所以需要有填充机制。实验要求采用 PKCS#7 Padding：如果需要填充 N 个字节长度才是 16 bytes 的倍数，则需要填充 N 个 0x0N。

例如：明文长度为 10 字节，需要填充 6 个 0x06；明文长度为 20 字节，就需要填充 12 个 0x0C。特别地，若明文刚好是 16 字节的倍数，也要填充 16 个 0x10。

代码如下：

1
2
3
4
5
6
void padding(uint8_t *data, uint32_t padding_len)
{
    uint8_t padding_content = (uint8_t)padding_len;
    for (uint32_t i = 0; i <= padding_len; i++)
        data[16 - i] = padding_content;
}

轮密钥加

有限域内的加法就是异或：

1
2
3
4
5
6
void add_round_key(uint8_t state[4][4], uint8_t round_key[4][4])
{
    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] ^= round_key[i][j];
}

字节替代

具体原理就不再赘述了，只需再写一个函数方便处理一个块就行：

1
2
3
4
5
void sub_bytes(uint8_t state[4][4])
{
    for (int i = 0; i < 4; i++)
        sub_word(state[i]);
}

行移位

将状态矩阵的第 i 行循环左移 i 个字节：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
void shift_row(uint8_t state[4][4])
{
    uint32_t temp[4][4];
    for (int i = 0; i < 4; i++)
    {
        temp[0][i] = state[0][i];
        temp[1][i] = state[1][(i + 1) % 4];
        temp[2][i] = state[2][(i + 2) % 4];
        temp[3][i] = state[3][(i + 3) % 4];
    }

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] = temp[i][j];
}

列混淆

列混淆的本质是做有限域内的矩阵乘法。列混淆矩阵是固定的：

[0x020x030x010x010x010x020x030x010x010x010x020x030x030x010x010x02]\begin{bmatrix} \text{0x02} & \text{0x03} & \text{0x01} & \text{0x01} \\ \text{0x01} & \text{0x02} & \text{0x03} & \text{0x01} \\ \text{0x01} & \text{0x01} & \text{0x02} & \text{0x03} \\ \text{0x03} & \text{0x01} & \text{0x01} & \text{0x02} \\ \end{bmatrix}

用其左乘状态矩阵即可：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
void mix_column(uint8_t state[4][4])
{
    uint32_t temp[4][4];
    for (int i = 0; i < 4; i++)
    {
        temp[0][i] = gf_mul(0x02, state[0][i]) ^ gf_mul(0x03, state[1][i]) ^ state[2][i] ^ state[3][i];
        temp[1][i] = state[0][i] ^ gf_mul(0x02, state[1][i]) ^ gf_mul(0x03, state[2][i]) ^ state[3][i];
        temp[2][i] = state[0][i] ^ state[1][i] ^ gf_mul(0x02, state[2][i]) ^ gf_mul(0x03, state[3][i]);
        temp[3][i] = gf_mul(0x03, state[0][i]) ^ state[1][i] ^ state[2][i] ^ gf_mul(0x02, state[3][i]);
    }

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] = temp[i][j];
}

其中gf_mul()是有限域内的乘法运算。虽然之前写过一个有限域运算的视线，但相比那个，这里的有限域太小了，可以直接利用生成元构造正表和反表，将有限域乘法转化为查表。

生成元是指有限域中的一个元素gg，它的幂次可以生成整个有限域，或者说有限域中的每个元素都能够被表示为生成元的幂次（0 除外）。对于有限域中的两个元素aa和bb，它们的乘法可以表示为：

a⋅b=gm⋅gn=gm+na \cdot b = g^m\cdot g^n = g^{m+n}

所以只需构造一张正表，能通过指数mm得到gmg^m，再构造一张反表，能通过gmg^m得到mm，就可以实现有限域的乘法运算。

我们选取生成元x+1x+1。有限域中乘以x+1x+1，相当于先左移一位，再加上原数。如果有进位，就减去 x8+x4+x3+x+1x^8+x^4+x^3+x+1，也即与0x11b异或。所以可以写出生成正表和反表的代码：

1
2
3
4
5
6
7
8
9
10
11
12
void init_tables()
{
    table[0] = 1; // g^0
    for (int i = 1; i < 255; i++)
    {
        table[i] = (table[i - 1] << 1) ^ table[i - 1];
        if (table[i] & 0x100)
            table[i] ^= 0x11b;
    }
    for (int i = 0; i < 255; i++)
        inverse_table[table[i]] = i;
}

当然，你也可以生成一次之后，将其作为常量数组写进代码里。这样效率并无明显提高，所以最终还是选择了动态生成。

最终，可以写出有限域乘法的代码：

1
2
3
4
5
6
7
uint32_t gf_mul(uint32_t x, uint32_t y)
{
    if (x == 0 || y == 0)
        return 0;

    return table[(inverse_table[x] + inverse_table[y]) % 255];
}

完整加密代码

由于 CBC 模式下，每一个明文块被加密前，要先与上一个密文块异或。第一个明文块与初始随机向量IV疑惑。

我们用last这个变量来存储，并将其初始化为IV。每次加密完一个块后，将其赋值给last，以便下一次异或。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
void encrypt(uint8_t state[4][4], uint8_t key[4][4], uint8_t last[4][4], uint8_t res[4][4])
{
    uint8_t expand_key[11][4][4];
    key_expand(key, expand_key);

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] ^= last[i][j];

    // Round 0
    add_round_key(state, expand_key[0]);

    // Round 1-10
    for (int j = 1; j <= 10; j++)
    {
        sub_bytes(state);
        shift_row(state);
        if (j != 10)
            mix_column(state);
        add_round_key(state, expand_key[j]);
    }

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            res[i][j] = state[j][i];

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            last[i][j] = state[i][j];
}

解密过程

解密过程我实现的是普通的“逆着回去”，沿着加密过程的反方向来。首先要实现一系列的逆操作

逆轮密钥加

根据异或的性质，再异或一个相同的数就可以还原，直接用原来的add_round_key()函数即可

逆字节替代

有 S 盒自然也有逆 S 盒，也是直接查表即可

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
const uint8_t INV_SBOX[16][16] = {
    {0x52, 0x09, 0x6a, 0xd5, 0x30, 0x36, 0xa5, 0x38, 0xbf, 0x40, 0xa3, 0x9e, 0x81, 0xf3, 0xd7, 0xfb},
    {0x7c, 0xe3, 0x39, 0x82, 0x9b, 0x2f, 0xff, 0x87, 0x34, 0x8e, 0x43, 0x44, 0xc4, 0xde, 0xe9, 0xcb},
    {0x54, 0x7b, 0x94, 0x32, 0xa6, 0xc2, 0x23, 0x3d, 0xee, 0x4c, 0x95, 0x0b, 0x42, 0xfa, 0xc3, 0x4e},
    {0x08, 0x2e, 0xa1, 0x66, 0x28, 0xd9, 0x24, 0xb2, 0x76, 0x5b, 0xa2, 0x49, 0x6d, 0x8b, 0xd1, 0x25},
    {0x72, 0xf8, 0xf6, 0x64, 0x86, 0x68, 0x98, 0x16, 0xd4, 0xa4, 0x5c, 0xcc, 0x5d, 0x65, 0xb6, 0x92},
    {0x6c, 0x70, 0x48, 0x50, 0xfd, 0xed, 0xb9, 0xda, 0x5e, 0x15, 0x46, 0x57, 0xa7, 0x8d, 0x9d, 0x84},
    {0x90, 0xd8, 0xab, 0x00, 0x8c, 0xbc, 0xd3, 0x0a, 0xf7, 0xe4, 0x58, 0x05, 0xb8, 0xb3, 0x45, 0x06},
    {0xd0, 0x2c, 0x1e, 0x8f, 0xca, 0x3f, 0x0f, 0x02, 0xc1, 0xaf, 0xbd, 0x03, 0x01, 0x13, 0x8a, 0x6b},
    {0x3a, 0x91, 0x11, 0x41, 0x4f, 0x67, 0xdc, 0xea, 0x97, 0xf2, 0xcf, 0xce, 0xf0, 0xb4, 0xe6, 0x73},
    {0x96, 0xac, 0x74, 0x22, 0xe7, 0xad, 0x35, 0x85, 0xe2, 0xf9, 0x37, 0xe8, 0x1c, 0x75, 0xdf, 0x6e},
    {0x47, 0xf1, 0x1a, 0x71, 0x1d, 0x29, 0xc5, 0x89, 0x6f, 0xb7, 0x62, 0x0e, 0xaa, 0x18, 0xbe, 0x1b},
    {0xfc, 0x56, 0x3e, 0x4b, 0xc6, 0xd2, 0x79, 0x20, 0x9a, 0xdb, 0xc0, 0xfe, 0x78, 0xcd, 0x5a, 0xf4},
    {0x1f, 0xdd, 0xa8, 0x33, 0x88, 0x07, 0xc7, 0x31, 0xb1, 0x12, 0x10, 0x59, 0x27, 0x80, 0xec, 0x5f},
    {0x60, 0x51, 0x7f, 0xa9, 0x19, 0xb5, 0x4a, 0x0d, 0x2d, 0xe5, 0x7a, 0x9f, 0x93, 0xc9, 0x9c, 0xef},
    {0xa0, 0xe0, 0x3b, 0x4d, 0xae, 0x2a, 0xf5, 0xb0, 0xc8, 0xeb, 0xbb, 0x3c, 0x83, 0x53, 0x99, 0x61},
    {0x17, 0x2b, 0x04, 0x7e, 0xba, 0x77, 0xd6, 0x26, 0xe1, 0x69, 0x14, 0x63, 0x55, 0x21, 0x0c, 0x7d}};

void inv_sub_byte(uint8_t *byte)
{
    int row = *byte >> 4;
    int col = *byte & 0x0f;
    *byte = INV_SBOX[row][col];
}

void inv_sub_word(uint8_t word[4])
{
    for (int i = 0; i < 4; i++)
        inv_sub_byte(&word[i]);
}

void inv_sub_bytes(uint8_t state[4][4])
{
    for (int i = 0; i < 4; i++)
        inv_sub_word(state[i]);
}

逆行移位

原来是左移，逆操作就右移回去即可：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
void inv_shift_row(uint8_t state[4][4])
{
    uint32_t temp[4][4];
    for (int i = 0; i < 4; i++)
    {
        temp[0][i] = state[0][i];
        temp[1][i] = state[1][(i + 3) % 4];
        temp[2][i] = state[2][(i + 2) % 4];
        temp[3][i] = state[3][(i + 1) % 4];
    }

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] = temp[i][j];
}

逆列混淆

列混淆是乘一个举证，逆操作就是左乘逆矩阵。列混淆矩阵的逆矩阵如下：

[0x0e0x0b0x0d0x090x090x0e0x0b0x0d0x0d0x090x0e0x0b0x0b0x0d0x090x0e]\begin{bmatrix} \text{0x0e} & \text{0x0b} & \text{0x0d} & \text{0x09} \\ \text{0x09} & \text{0x0e} & \text{0x0b} & \text{0x0d} \\ \text{0x0d} & \text{0x09} & \text{0x0e} & \text{0x0b} \\ \text{0x0b} & \text{0x0d} & \text{0x09} & \text{0x0e} \\ \end{bmatrix}

所以可写出代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
void inv_mix_column(uint8_t state[4][4])
{
    uint32_t temp[4][4];
    for (int i = 0; i < 4; i++)
    {
        temp[0][i] = gf_mul(0x0e, state[0][i]) ^ gf_mul(0x0b, state[1][i]) ^ gf_mul(0x0d, state[2][i]) ^ gf_mul(0x09, state[3][i]);
        temp[1][i] = gf_mul(0x09, state[0][i]) ^ gf_mul(0x0e, state[1][i]) ^ gf_mul(0x0b, state[2][i]) ^ gf_mul(0x0d, state[3][i]);
        temp[2][i] = gf_mul(0x0d, state[0][i]) ^ gf_mul(0x09, state[1][i]) ^ gf_mul(0x0e, state[2][i]) ^ gf_mul(0x0b, state[3][i]);
        temp[3][i] = gf_mul(0x0b, state[0][i]) ^ gf_mul(0x0d, state[1][i]) ^ gf_mul(0x09, state[2][i]) ^ gf_mul(0x0e, state[3][i]);
    }

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] = temp[i][j];
}

完整解密代码

是 CBC 模式，记得与上一个密文块异或才得到原来的明文块：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
void decrypt(uint8_t state[4][4], uint8_t key[4][4], uint8_t last[4][4], uint8_t res[4][4])
{
    uint8_t expand_key[11][4][4];
    key_expand(key, expand_key);

    // Round 10~1
    for (int j = 10; j >= 1; j--)
    {
        add_round_key(state, expand_key[j]);
        if (j != 10)
            inv_mix_column(state);
        inv_shift_row(state);
        inv_sub_bytes(state);
    }

    // Round 0
    add_round_key(state, expand_key[0]);

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            state[i][j] ^= last[i][j];

    for (int i = 0; i < 4; i++)
        for (int j = 0; j < 4; j++)
            res[i][j] = state[j][i];
}

SIMD 加速

可以利用 AES-NI 中的专用指令来加速，这份参考资料中有详细的说明以及 C++示例代码。我将其改写为了 C 语言版本。

需要以下头文件：

1
2
#include <emmintrin.h>
#include <immintrin.h>

拓展密钥

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
void key_expand_simd(uint8_t *key, __m128i *expand_key)
{
#define EXPAND(n, Rcon)                           \
    t = _mm_slli_si128(work, 4);                  \
    t = _mm_xor_si128(work, t);                   \
    t2 = _mm_slli_si128(t, 8);                    \
    t = _mm_xor_si128(t, t2);                     \
    work = _mm_aeskeygenassist_si128(work, Rcon); \
    work = _mm_shuffle_epi32(work, 0xFF);         \
    work = _mm_xor_si128(t, work);                \
    expand_key[n] = work;

    __m128i work = _mm_loadu_si128((const __m128i *)(key));
    expand_key[0] = work;
    __m128i t, t2;

    EXPAND(1, 0x01);
    EXPAND(2, 0x02);
    EXPAND(3, 0x04);
    EXPAND(4, 0x08);
    EXPAND(5, 0x10);
    EXPAND(6, 0x20);
    EXPAND(7, 0x40);
    EXPAND(8, 0x80);
    EXPAND(9, 0x1B);
    EXPAND(10, 0x36);
}

其中_mm_aeskeygenassist_si128()是生成轮密钥的辅助函数。总之按EXPAND宏的定义，经过一系列操作就可以得到下一个轮密钥。

加密

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
void encrypt_simd(uint8_t *data, uint8_t *key, __m128i *last, uint8_t *res)
{
    __m128i expand_key[11];
    key_expand_simd(key, expand_key);

    __m128i state = _mm_loadu_si128((const __m128i *)(data));
    state = _mm_xor_si128(state, expand_key[0]);
    state = _mm_xor_si128(state, *last);

    for (int i = 1; i <= 9; i++)
        state = _mm_aesenc_si128(state, expand_key[i]);

    state = _mm_aesenclast_si128(state, expand_key[10]);

    _mm_storeu_si128((__m128i *)(res), state);
    _mm_storeu_si128(last, state);
}

可以发现，key和data无需调整端序/行优先，直接用按字节顺序读入的uint8_t数组就可以正确地初始化所需的__m128i，非常方便（当然，要在主函数调整传入的参数）。

中间轮用一个函数_mm_aesenc_si128()就可以完成，最后一轮没有列混淆，所以用_mm_aesenclast_si128()。

还要注意，这里的last变为了__m128i类型，主函数在做初始化的时候，也要做相应的调整。

解密

AES-NI 仅支持等价的解密过程：操作顺序与加密过程中相同，只是操作要换为相应的逆操作，且用到的轮密钥不同。

先要生成解密过程中用的轮密钥：

1
2
3
4
5
6
7
void inv_key_expand_simd(uint8_t *key, __m128i *expand_key)
{
    key_expand_simd(key, expand_key);

    for (int i = 1; i <= 9; i++)
        expand_key[i] = _mm_aesimc_si128(expand_key[i]);
}

其中_mm_aesimc_si128()是进行逆列混淆。

具体的解密过程如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
void eq_decrypt_simd(uint8_t *data, uint8_t *key, __m128i *last, uint8_t *res)
{
    __m128i expand_key[11];
    inv_key_expand_simd(key, expand_key);

    __m128i state = _mm_loadu_si128((const __m128i *)(data));
    state = _mm_xor_si128(state, expand_key[10]);

    for (int i = 9; i >= 1; i--)
        state = _mm_aesdec_si128(state, expand_key[i]);

    state = _mm_aesdeclast_si128(state, expand_key[0]);
    state = _mm_xor_si128(state, *last);
    _mm_storeu_si128((__m128i *)(res), state);
}

可见与加密的代码几乎完全一样。其中_mm_aesdec_si128()和_mm_aesdeclast_si128()分别是解密过程中的中间轮和最后一轮。

性能比较

至此 AES-128 CBC 模式的普通版本和 SIMD 加速版本就都实现完了。前者在实验平台的评测机上，通过所有样例的耗时是 1751ms；后者的耗时则是 143ms，提升了一个数量级。

另外，这个版本中只需要实现一个模式，因此一些预处理过程与核心加密解密过程耦合了，直接写到了主函数中，这样是不太好的。接下来实现所有模式，会有所改进。

AES-128 所有模式实现

对应着实验 3-2，以及仓库中的aes_all.c。

解耦

首先我们将核心加密/解密过程解耦出来，让加密函数只管对一个块进行加密，输出密文块；解密函数只管对一个块进行解密，输出明文块。

加密和解密函数修改如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
void encrypt_simd(__m128i *data, __m128i *expand_key, __m128i *res)
{
    __m128i state = _mm_xor_si128(*data, expand_key[0]);

    for (int i = 1; i <= 9; i++)
        state = _mm_aesenc_si128(state, expand_key[i]);

    state = _mm_aesenclast_si128(state, expand_key[10]);

    _mm_storeu_si128(res, state);
}

void eq_decrypt_simd(__m128i *data, __m128i *expand_key, __m128i *res)
{
    __m128i state = _mm_xor_si128(*data, expand_key[10]);

    for (int i = 9; i >= 1; i--)
        state = _mm_aesdec_si128(state, expand_key[i]);

    state = _mm_aesdeclast_si128(state, expand_key[0]);

    _mm_storeu_si128(res, state);
}

这就非常舒服！

主函数

解耦后，在主函数中，只根据类型进行判断，调用相应的加密解密函数：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
int main(void)
{
#ifndef ONLINE_JUDGE
    FILE *in = fopen("aes_ctr_i2.bin", "rb");
    FILE *out = fopen("aeso.bin", "wb");
#else
#define in stdin
#define out stdout
#endif

    uint8_t mode;
    uint8_t key[16];
    uint8_t IV[16];
    uint32_t len;

    fread(&mode, sizeof(uint8_t), 1, in);
    fread(key, sizeof(uint8_t), 16, in);
    fread(IV, sizeof(uint8_t), 16, in);
    fread(&len, sizeof(uint32_t), 1, in);

    __m128i expand_key[11];
    if (mode == 0x80 || mode == 0x81)
        inv_key_expand_simd(key, expand_key);
    else
        key_expand_simd(key, expand_key);

    switch (mode)
    {
    case 0x00:
        ecb_encrypt(expand_key, len, in, out);
        break;
    case 0x01:
        cbc_encrypt(expand_key, IV, len, in, out);
        break;
    case 0x02:
        cfb_encrypt(expand_key, IV, len, in, out);
        break;
    case 0x03:
        ofb_encrypt(expand_key, IV, len, in, out);
        break;
    case 0x04:
        ctr_encrypt(expand_key, IV, len, in, out);
        break;
    case 0x80:
        ecb_decrypt(expand_key, len, in, out);
        break;
    case 0x81:
        cbc_decrypt(expand_key, IV, len, in, out);
        break;
    case 0x82:
        cfb_decrypt(expand_key, IV, len, in, out);
        break;
    case 0x83:
        ofb_decrypt(expand_key, IV, len, in, out);
        break;
    case 0x84:
        ctr_decrypt(expand_key, IV, len, in, out);
        break;
    }

    return 0;
}

判断调用那种密钥拓展函数，后面会介绍到。

CBC 模式

CBC 模式要求，每一个明文块在加密前，要先与上一个密文块异或。第一个明文块与初始向量IV异或。其他特点可以自行搜索，这里只讲述与实现相关的部分。

CBC 模式已经实现过了，我们只要简单修改一下参数类型即可。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
void cbc_encrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16];
    __m128i last = _mm_loadu_si128((const __m128i *)(IV));
    __m128i res;

    uint32_t i = 0;
    while (i < len)
    {
        if (i + 16 > len)
        {
            fread(data, sizeof(uint8_t), len - i, in);
            padding(data, 16 - (len - i));
        }
        else
            fread(data, sizeof(uint8_t), 16, in);

        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));
        data_ = _mm_xor_si128(data_, last);
        encrypt_simd(&data_, expand_key, &res);
        last = res;

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), 16, out);
        i += 16;
    }

    if (i == len)
    {
        padding(data, 16);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));
        data_ = _mm_xor_si128(data_, last);
        encrypt_simd(&data_, expand_key, &res);
        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), 16, out);
    }
}

void cbc_decrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16];
    __m128i last = _mm_loadu_si128((const __m128i *)(IV));
    __m128i res;

    uint32_t i = 0;
    while (i < len)
    {
        fread(data, sizeof(uint8_t), 16, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        eq_decrypt_simd(&data_, expand_key, &res);
        res = _mm_xor_si128(res, last);

        int len_ = 16;
        if (i + 16 == len)
        {
            uint8_t *res_ = (uint8_t *)&res;
            len_ = 16 - (int)res_[15];
        }

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);

        i += 16;
        last = data_;
    }
}

ECB 模式

ECB 模式是最简单的模式，每一个明文块都独立加密即可。与 CBC 模式一样，明文也需要填充，使明文长度是 16 字节的倍数。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
void ecb_encrypt(__m128i *expand_key, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16];
    __m128i res;

    uint32_t i = 0;
    while (i < len)
    {
        if (i + 16 > len)
        {
            fread(data, sizeof(uint8_t), len - i, in);
            padding(data, 16 - (len - i));
        }
        else
            fread(data, sizeof(uint8_t), 16, in);

        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));
        encrypt_simd(&data_, expand_key, &res);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), 16, out);
        i += 16;
    }

    if (i == len)
    {
        padding(data, 16);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));
        encrypt_simd(&data_, expand_key, &res);
        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), 16, out);
    }
}

void ecb_decrypt(__m128i *expand_key, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16];
    __m128i res;

    uint32_t i = 0;
    while (i < len)
    {
        fread(data, sizeof(uint8_t), 16, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        eq_decrypt_simd(&data_, expand_key, &res);

        int len_ = 16;
        if (i + 16 == len)
        {
            uint8_t *res_ = (uint8_t *)&res;
            len_ = 16 - (int)res_[15];
        }

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);

        i += 16;
    }
}

CFB 模式

CFB 模式是一种流模式。记明文块为xix_i，密文块为yiy_i，加密函数为EE，y0=IVy_0=IV。第一个明文块是x1x_1。

通过加密以前的密文分组得到密钥流元素：

zi=E(yi−1)z_i = E(y_{i-1})

然后让密钥流元素与明文异或得到密文：

yi=xi⊕ziy_i = x_i \oplus z_i

解密过程也是先计算出密钥流元素，然后与密文异或得到明文。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
void cfb_encrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16] = {0};
    __m128i z = _mm_loadu_si128((const __m128i *)(IV));
    __m128i res;

    for (uint32_t i = 0; i < len; i += 16)
    {
        encrypt_simd(&z, expand_key, &res);

        int len_ = 16;
        if (i + 16 > len)
            len_ = len - i;

        fread(data, sizeof(uint8_t), len_, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        res = _mm_xor_si128(res, data_);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);

        z = res;
    }
}

void cfb_decrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16] = {0};
    __m128i z = _mm_loadu_si128((const __m128i *)(IV));
    __m128i res;

    for (uint32_t i = 0; i < len; i += 16)
    {
        encrypt_simd(&z, expand_key, &res);

        int len_ = 16;
        if (i + 16 > len)
            len_ = len - i;
        fread(data, sizeof(uint8_t), len_, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        res = _mm_xor_si128(res, data_);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);

        z = data_;
    }
}

CFB 模式不要求明文长度是 16 字节的倍数，所以不需要填充。由于初始与 IV 异或，得到的密文块是 128 位的，之后的密文块会与 128 位的上一个密文块异或。所以密文块总是 128 位的。

解密时，要记得根据长度截断。

OFB 模式

OFB 模式也是一种流模式但是是通过反复加密 IV 得到密钥流元素。

记z0=IVz_0=IV，则有：

zi=E(zi−1)z_i = E(z_{i-1})

然后让密钥流元素与明文异或得到密文：

yi=xi⊕ziy_i = x_i \oplus z_i 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
void ofb_encrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16] = {0};
    __m128i z = _mm_loadu_si128((const __m128i *)(IV));
    __m128i res;

    for (uint32_t i = 0; i < len; i += 16)
    {
        encrypt_simd(&z, expand_key, &res);

        z = res;

        int len_ = 16;
        if (i + 16 > len)
            len_ = len - i;

        fread(data, sizeof(uint8_t), len_, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        res = _mm_xor_si128(res, data_);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);
    }
}

void ofb_decrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16] = {0};
    __m128i z = _mm_loadu_si128((const __m128i *)(IV));
    __m128i res;

    for (uint32_t i = 0; i < len; i += 16)
    {
        encrypt_simd(&z, expand_key, &res);

        z = res;

        int len_ = 16;
        if (i + 16 > len)
            len_ = len - i;
        fread(data, sizeof(uint8_t), len_, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        res = _mm_xor_si128(res, data_);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);
    }
}

同样，解密的时候要注意截断。

CTR 模式

CTR 模式也是一个流模式，反复加密一个不断递增的计时器得到密钥流元素。

实验中给的 CTR 模式下的样例是来源于 NIST 的文档 Recommendation for Block Cipher Modes of Operation附录。在样例中，计时器的初始值就是 IV，且是以大端序给出。计时器递增就是直接加 1。

网络上也有资料，说计时器由长度相同的两部分组成。前半部分是一个随机向量nonce，后半部分是一个计数器counter。计数器从 0 开始，每次加 1。

本次实验中，我们采用前一种方式。后面一种方式只要做简单的修改就能实现了。

具体而言，记初始计数器ctr=IVctr=IV，则有：

Ti=ctr+i−1mod  2128T_i=ctr+i-1 \mod 2^{128}

密钥流元素：

zi=E(Ti)z_i=E(T_i)

密文：

yi=xi⊕ziy_i=x_i \oplus z_i 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
void ctr_encrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16] = {0};
    __m128i ctr = _mm_loadu_si128((const __m128i *)IV);
    __m128i mask = _mm_set_epi8(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15);
    ctr = _mm_shuffle_epi8(ctr, mask);

    __m128i res;

    for (uint32_t i = 0; i < len; i += 16)
    {
        __m128i T = _mm_shuffle_epi8(ctr, mask);

        encrypt_simd(&T, expand_key, &res);

        int len_ = 16;
        if (i + 16 > len)
            len_ = len - i;

        fread(data, sizeof(uint8_t), len_, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        res = _mm_xor_si128(res, data_);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);

        ctr = _mm_add_epi64(ctr, _mm_set_epi64x(0, 1));
    }
}

void ctr_decrypt(__m128i *expand_key, uint8_t *IV, uint32_t len, FILE *in, FILE *out)
{
    uint8_t data[16] = {0};
    __m128i ctr = _mm_loadu_si128((const __m128i *)IV);
    __m128i mask = _mm_set_epi8(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15);
    ctr = _mm_shuffle_epi8(ctr, mask);

    __m128i res;

    for (uint32_t i = 0; i < len; i += 16)
    {
        __m128i T = _mm_shuffle_epi8(ctr, mask);

        encrypt_simd(&T, expand_key, &res);

        int len_ = 16;
        if (i + 16 > len)
            len_ = len - i;

        fread(data, sizeof(uint8_t), len_, in);
        __m128i data_ = _mm_loadu_si128((const __m128i *)(data));

        res = _mm_xor_si128(res, data_);

        uint8_t *res_ = (uint8_t *)&res;
        fwrite(res_, sizeof(uint8_t), len_, out);

        ctr = _mm_add_epi64(ctr, _mm_set_epi64x(0, 1));
    }
}

解密函数中，计时器是大端读入的，所以要用_mm_shuffle_epi8()调整一下。然后利用_mm_add_epi64()实现计数器加 1。

至此所有模式都实现完了。在实验平台上成为第二个通过所有样例的程序。