Chapter 15: The Data Encryption Standard (SageMath)#

“For 20 years, the DES was the most widely used encryption algorithm in the world.” — Bruce Schneier, Applied Cryptography, 1996

In this chapter we implement the complete Data Encryption Standard (DES) from scratch using SageMath, study its Feistel structure, and analyze its security properties — including weak keys discovered via Gröbner bases and the avalanche effect.

Python Version Available

This notebook contains the SageMath implementation. A pure Python/NumPy version is available in Chapter 15: The Data Encryption Standard.

15.1 Historical Context#

In 1972, the U.S. National Bureau of Standards (NBS, now NIST) issued a public call for a standard encryption algorithm. IBM submitted a modified version of its Lucifer cipher, designed by Horst Feistel and refined by a team including Don Coppersmith. After review by the NSA — which controversially reduced the key length from 128 to 56 bits and modified the S-boxes — the algorithm was published as FIPS 46 in January 1977.

DES became the dominant symmetric cipher for over two decades. Its eventual demise came not from cryptanalytic weakness but from its 56-bit key length:

Year

Event

1977

DES published as FIPS 46

1993

Matsui’s linear cryptanalysis (\(2^{43}\) known plaintexts)

1997

RSA DES Challenge: broken in 96 days by distributed computing

1998

EFF Deep Crack: broken in 56 hours ($250,000 hardware)

1999

Deep Crack + distributed.net: broken in 22 hours

2001

AES selected as replacement

Tip

Although DES is now obsolete, its Feistel structure influenced nearly every subsequent block cipher. Understanding DES is essential for understanding the attacks (linear and differential cryptanalysis) that shaped modern cipher design.

15.2 The Feistel Structure#

Definition 15.1 — Feistel Network

A Feistel network of \(r\) rounds splits the input block into two halves \((L_0, R_0)\) and applies:

\[ L_i = R_{i-1}, \quad R_i = L_{i-1} \oplus F(R_{i-1}, K_i)\]

where \(F\) is the round function and \(K_i\) is the \(i\)-th round key.

Theorem 15.1 — Feistel Invertibility

A Feistel cipher is always invertible regardless of whether \(F\) is invertible. Decryption uses the same structure with round keys in reverse order.

Definition 15.2 — DES Round Function

The DES round function \(F(R, K)\) consists of:

  1. Expansion \(E\): 32 bits → 48 bits

  2. Key mixing: XOR with 48-bit round key

  3. S-box substitution: eight 6→4-bit S-boxes

  4. Permutation \(P\): 32-bit straight permutation

\[ F(R, K) = P\bigl(S_1(B_1) \| S_2(B_2) \| \cdots \| S_8(B_8)\bigr)\]

where \(B_1 \| B_2 \| \cdots \| B_8 = E(R) \oplus K\).

15.3 Implementation#

We implement the complete DES algorithm using SageMath. The initial permutation is computed using SymmetricGroup and PermutationGroupElement, and the XOR operation uses SageMath’s native ^^ operator over \(\mathrm{GF}(2)\).

All standard tables (IP, IP\(^{-1}\), E, P, S-boxes, PC-1, PC-2) are included.

Take a look at 0x10001/des to see a more optimized Python implementation.

Permutation Tables#

The initial permutation IP and its inverse are standard DES tables. We compute IP\(^{-1}\) algebraically from the group inverse using SageMath’s SymmetricGroup.

#Setup
INITIAL_PERMUTATION = [
    57, 49, 41, 33, 25, 17, 9,  1,
    59, 51, 43, 35, 27, 19, 11, 3,
    61, 53, 45, 37, 29, 21, 13, 5,
    63, 55, 47, 39, 31, 23, 15, 7,
    56, 48, 40, 32, 24, 16, 8,  0,
    58, 50, 42, 34, 26, 18, 10, 2,
    60, 52, 44, 36, 28, 20, 12, 4,
    62, 54, 46, 38, 30, 22, 14, 6,
]
G = SymmetricGroup(range(0,64))
IPperm=PermutationGroupElement(INITIAL_PERMUTATION,G)

INVERSE_PERMUTATION=[(IPperm**(-1))(i) for i in range(0,64)]

def IP(data):
    return [data[INITIAL_PERMUTATION[k]] for k in range(0,64)]
def IPinv(data):
    return [data[INVERSE_PERMUTATION[k]] for k in range(0,64)]

EXPANSION =[
    31, 0,  1,  2,  3,  4,
    3,  4,  5,  6,  7,  8,
    7,  8,  9,  10, 11, 12,
    11, 12, 13, 14, 15, 16,
    15, 16, 17, 18, 19, 20,
    19, 20, 21, 22, 23, 24,
    23, 24, 25, 26, 27, 28,
    27, 28, 29, 30, 31, 0,
    ]
def Expansion(data):
    return [data[EXPANSION[k]] for k in range(0,48)]
    
PERMUTATION_P = [
    15, 6,  19, 20, 28, 11, 27, 16,
    0,  14, 22, 25, 4,  17, 30, 9,
    1,  7,  23, 13, 31, 26, 2,  8,
    18, 12, 29, 5,  21, 10, 3,  24]

def PermuteP(R):
    return [R[PERMUTATION_P[k]] for k in range(0,32)]

print("DES permutation tables loaded.")
print(f"IP as group element: {IPperm}")
DES permutation tables loaded.
IP as group element: (0,57,54,12,27,39)(1,49,52,28,31,7)(2,41,50,44,26,47)(3,33,48,60,30,15)(4,25,55)(5,17,53,20,29,23)(6,9,51,36,24,63)(8,59,38)(10,43,34,40,58,46)(11,35,32,56,62,14)(13,19,37,16,61,22)(18,45)
IPperm
(0,57,54,12,27,39)(1,49,52,28,31,7)(2,41,50,44,26,47)(3,33,48,60,30,15)(4,25,55)(5,17,53,20,29,23)(6,9,51,36,24,63)(8,59,38)(10,43,34,40,58,46)(11,35,32,56,62,14)(13,19,37,16,61,22)(18,45)

Helper Functions#

Utility functions for splitting blocks and performing bitwise XOR. Note that SageMath uses ^^ for XOR (the ^ operator is exponentiation in Sage).

def Split(data):
    return data[:32],data[32:]

def Unsplit(L,R):
    return L+R

def BSplit(data):
    return [data[6*i:6*i+6] for i in range(0,8)]
    
def ListXOR(li1,li2):
    return [li1[i]^^li2[i] for i in range(0,len(li1))]

The S-Box Function#

The S-box function is the nonlinear core of DES. Each of the eight S-boxes maps a 6-bit input to a 4-bit output. The outer bits (bits 0 and 5) select the row, and the inner bits (bits 1–4) select the column.

In SageMath, we use the .digits(base=2, padto=4) method on integers to convert S-box output values to 4-bit lists.

SUBSTITUTION_BOX = [
    [
        [14, 4,  13, 1,  2,  15, 11, 8,  3,  10, 6,  12, 5,  9,  0,  7],
        [0,  15, 7,  4,  14, 2,  13, 1,  10, 6,  12, 11, 9,  5,  3,  8],
        [4,  1,  14, 8,  13, 6,  2,  11, 15, 12, 9,  7,  3,  10, 5,  0],
        [15, 12, 8,  2,  4,  9,  1,  7,  5,  11, 3,  14, 10, 0,  6,  13],
    ],
    [
        [15, 1,  8,  14, 6,  11, 3,  4,  9,  7,  2,  13, 12, 0,  5,  10],
        [3,  13, 4,  7,  15, 2,  8,  14, 12, 0,  1,  10, 6,  9,  11, 5],
        [0,  14, 7,  11, 10, 4,  13, 1,  5,  8,  12, 6,  9,  3,  2,  15],
        [13, 8,  10, 1,  3,  15, 4,  2,  11, 6,  7,  12, 0,  5,  14, 9],
    ],
    [
        [10, 0,  9,  14, 6,  3,  15, 5,  1,  13, 12, 7,  11, 4,  2,  8],
        [13, 7,  0,  9,  3,  4,  6,  10, 2,  8,  5,  14, 12, 11, 15, 1],
        [13, 6,  4,  9,  8,  15, 3,  0,  11, 1,  2,  12, 5,  10, 14, 7],
        [1,  10, 13, 0,  6,  9,  8,  7,  4,  15, 14, 3,  11, 5,  2,  12],
    ],
    [
        [7,  13, 14, 3,  0,  6,  9,  10, 1,  2,  8,  5,  11, 12, 4,  15],
        [13, 8,  11, 5,  6,  15, 0,  3,  4,  7,  2,  12, 1,  10, 14, 9],
        [10, 6,  9,  0,  12, 11, 7,  13, 15, 1,  3,  14, 5,  2,  8,  4],
        [3,  15, 0,  6,  10, 1,  13, 8,  9,  4,  5,  11, 12, 7,  2,  14],
    ],
    [
        [2,  12, 4,  1,  7,  10, 11, 6,  8,  5,  3,  15, 13, 0,  14, 9],
        [14, 11, 2,  12, 4,  7,  13, 1,  5,  0,  15, 10, 3,  9,  8,  6],
        [4,  2,  1,  11, 10, 13, 7,  8,  15, 9,  12, 5,  6,  3,  0,  14],
        [11, 8,  12, 7,  1,  14, 2,  13, 6,  15, 0,  9,  10, 4,  5,  3],
    ],
    [
        [12, 1,  10, 15, 9,  2,  6,  8,  0,  13, 3,  4,  14, 7,  5,  11],
        [10, 15, 4,  2,  7,  12, 9,  5,  6,  1,  13, 14, 0,  11, 3,  8],
        [9,  14, 15, 5,  2,  8,  12, 3,  7,  0,  4,  10, 1,  13, 11, 6],
        [4,  3,  2,  12, 9,  5,  15, 10, 11, 14, 1,  7,  6,  0,  8,  13],
    ],
    [
        [4,  11,  2, 14, 15, 0,  8,  13, 3,  12, 9,  7,  5,  10, 6,  1],
        [13, 0,  11, 7,  4,  9,  1,  10, 14, 3,  5,  12, 2,  15, 8,  6],
        [1,  4,  11, 13, 12, 3,  7,  14, 10, 15, 6,  8,  0,  5,  9,  2],
        [6,  11, 13, 8,  1,  4,  10, 7,  9,  5,  0,  15, 14, 2,  3,  12],
    ],
    [
        [13, 2,  8,  4,  6,  15, 11, 1,  10, 9,  3,  14, 5,  0,  12, 7],
        [1,  15, 13, 8,  10, 3,  7,  4,  12, 5,  6,  11, 0,  14, 9,  2],
        [7,  11, 4,  1,  9,  12, 14, 2,  0,  6,  10, 13, 15, 3,  5,  8],
        [2,  1,  14, 7,  4,  10, 8,  13, 15, 12, 9,  0,  3,  5,  6,  11],
    ],
]

def SBox(data48):
    Blist=BSplit(data48)
    Sout=[]
    for i in range(0,8):
        B=Blist[i]
        S=SUBSTITUTION_BOX[i] 
        b0,b1,b2,b3,b4,b5=B
        r=2*b0+b5 #row selection
        c=8*b1+4*b2+2*b3+b4 #column selection
        Sout+=S[r][c].digits(base=2,padto=4)[::-1] #pad to 4 bits
    return Sout

The Feistel Function#

The Feistel function is a composition of the S-box substitution with XOR. It applies to the right half (\(R\)) combined with the round key (\(K\)), followed by the \(P\)-permutation.

def FeistelFunction(R,K):
    return PermuteP(SBox(ListXOR(Expansion(R),K)))

Key Schedule#

The key schedule expands the 64-bit key (with 8 parity bits) into sixteen 48-bit round keys via Permuted Choice 1 (PC-1), cyclic left rotations, and Permuted Choice 2 (PC-2).

#in rounds zeroth, first, eight we shift by 1
ROTATES = [1, 1, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 1]

PERMUTED_CHOICE1 = [
    56, 48, 40, 32, 24, 16, 8,
    0,  57, 49, 41, 33, 25, 17,
    9,  1,  58, 50, 42, 34, 26,
    18, 10, 2,  59, 51, 43, 35,
    62, 54, 46, 38, 30, 22, 14,
    6,  61, 53, 45, 37, 29, 21,
    13, 5,  60, 52, 44, 36, 28,
    20, 12, 4,  27, 19, 11, 3,
]

PERMUTED_CHOICE2 = [
    13, 16, 10, 23, 0,  4,
    2,  27, 14, 5,  20, 9,
    22, 18, 11, 3,  25, 7,
    15, 6,  26, 19, 12, 1,
    40, 51, 30, 36, 46, 54,
    29, 39, 50, 44, 32, 47,
    43, 48, 38, 55, 33, 52,
    45, 41, 49, 35, 28, 31,
]

#cyclic rotation by n positions to the left
def CyclicRot(li,n):
    l=len(li)
    return li[n%l:]+li[:-((l-n)%l)]

def KeySchedule(key):
    '''key has 64 bits; output: 16 keys Ki, each of size 48 bits
    '''
    PC1=[key[el] for el in PERMUTED_CHOICE1]
    C0=PC1[:28]
    D0=PC1[28:]
    Cim1=C0
    Dim1=D0
    Kli=[]
    for i in range(1,17):
        v=ROTATES[i-1]
        Ci=CyclicRot(Cim1,v)
        Di=CyclicRot(Dim1,v)
        CDi=Ci+Di
        PC2=[CDi[el] for el in PERMUTED_CHOICE2]
        Ki=PC2
        Kli.append(Ki)
        Cim1=Ci
        Dim1=Di
    return Kli

DES Encryption#

The first and final instruction of DES is a permutation IP and its inverse. They do not affect the security of the protocol and are kept only for historical and compatibility reasons (they are essential from the point of view of a hardware implementation; in the software version it is better to ignore them, however the results will be incompatible with standard DES).

In 16 rounds we add the round keys with the help of a Feistel function.

def DESEncryption(plaintext,key):
    expkey=KeySchedule(key)
    L=[]
    R=[]
    L0,R0=Split(IP(plaintext))
    L.append(L0)
    R.append(R0)
    Lim1=L0
    Rim1=R0
    for i in range(1,17):
        Li=Rim1
        Ki=expkey[i-1]
        Ri=ListXOR(Lim1,FeistelFunction(Rim1,Ki))
        L.append(Li)
        R.append(Ri)
        Lim1=Li
        Rim1=Ri
    return IPinv(Unsplit(Ri,Li))

DES Decryption#

In decryption we only reverse the order of round keys, exactly as predicted by the Feistel invertibility theorem (Theorem 15.1).

def DESDecryption(cipher,key):
    expkey=KeySchedule(key)
    L=[]
    R=[]
    L0,R0=Split(IP(cipher))
    L.append(L0)
    R.append(R0)
    Lim1=L0
    Rim1=R0
    for i in range(16,0,-1):
        Li=Rim1
        Ki=expkey[i-1]
        Ri=ListXOR(Lim1,FeistelFunction(Rim1,Ki))
        L.append(Li)
        R.append(Ri)
        Lim1=Li
        Rim1=Ri
    return IPinv(Unsplit(Ri,Li))

Key Generation#

Every 8th bit in the key is an odd parity bit. To generate a full 64-bit key we only need 56 bits of input; the other 8 are generated from parity.

def GenerateKey(key56input):
    key=[]
    for i in range(0,8):
        genpart=key56input[7*i:7*i+7]
        parity=(sum(genpart)+1)%2 #sum of every 7 bits should be odd
        key+=genpart+[parity]
    
    return key

Correctness Test#

A composition of encryption with decryption should recover the original plaintext in each case.

#testing
for _ in range(0,1000):
    data=[choice([0,1]) for _ in range(0,64)]
    key=GenerateKey([choice([0,1]) for _ in range(0,56)])
    cipher=DESEncryption(data,key)
    plain=DESDecryption(cipher,key)
    assert plain==data

print("Tested 1000 random encrypt/decrypt cycles: ALL PASSED")
Tested 1000 random encrypt/decrypt cycles: ALL PASSED

15.4 Key Space and Security#

Key Space

DES uses a 56-bit effective key, giving \(2^{56} \approx 7.2 \times 10^{16}\) possible keys. By modern standards this is far too small — an exhaustive search is feasible with specialized hardware or distributed computing.

Danger

DES should never be used for new applications. Triple-DES (3DES) provides 112-bit effective security but is also being phased out in favor of AES.

15.5 Experiments#

Experiment 1: Avalanche Effect#

A single bit flip in the plaintext should affect approximately half the ciphertext bits. We measure the avalanche across rounds using a parameterized DES encryption and SageMath’s list_plot.

We first define a version of DES that accepts a variable number of rounds.

def DESEncryptionParam(plaintext,key,n=17):
    expkey=KeySchedule(key)
    L=[]
    R=[]
    L0,R0=Split(IP(plaintext))
    L.append(L0)
    R.append(R0)
    Lim1=L0
    Rim1=R0
    for i in range(1,n):
        Li=Rim1
        Ki=expkey[i-1]
        Ri=ListXOR(Lim1,FeistelFunction(Rim1,Ki))
        L.append(Li)
        R.append(Ri)
        Lim1=Li
        Rim1=Ri
    return IPinv(Unsplit(Ri,Li))

def DetectDiff(plain,key,bitchange,n=17):
    cipher1=DESEncryptionParam(plain,key,n)
    plain[bitchange]=1-plain[bitchange]
    cipher2=DESEncryptionParam(plain,key,n)
    return cipher1,cipher2
def ChangeInfluence(pos,r):
    c1,c2=DetectDiff([choice([0,1]) for _ in range(0,64)],GenerateKey([choice([0,1]) for _ in range(0,56)]),pos,r)
    c1c2diff=ListXOR(c1,c2)
    return c1c2diff.count(1)/64.0
# Avalanche effect: fraction of differing bits vs. number of rounds
n_samples = 50
avalanche_data = [(i, sum([ChangeInfluence(0, i) for _ in range(n_samples)]) / n_samples)
                  for i in range(2, 18)]

P = list_plot(avalanche_data,
              plotjoined=True,
              marker='o',
              color='red',
              axes_labels=['Number of rounds', 'Fraction of differing bits'],
              title='DES Avalanche Effect: 1-Bit Plaintext Change')
P += line([(2, 0.5), (17, 0.5)], color='green', linestyle='--',
          legend_label='Ideal (50%)')
P.show(figsize=(8, 4))
../_images/c3e5e5d4e869390d18c62515601aa47feee8a9a00f8cb44a2d000b6a3392b50f.png

Observation

Full diffusion is achieved by approximately round 5–6. By round 16, a single bit change in the plaintext produces an average of ~32 differing ciphertext bits — exactly the ideal 50% threshold.

Experiment 2: Weak Keys via Gröbner Bases#

DES has 4 weak keys where all 16 round keys are identical. Encrypting twice with a weak key returns the original plaintext.

We discover these algebraically using SageMath: model the 56 key bits as variables in PolynomialRing(GF(2), 56, 'k'), require all round keys to be equal, compute a Gröbner basis, and find the kernel of the resulting linear system.

The answer is:

  • 01010101 01010101

  • E0E0E0E0 F1F1F1F1

  • 1F1F1F1F 0E0E0E0E

  • FEFEFEFE FEFEFEFE

Below we explain where this answer comes from (a linear algebra task!).

#Find a system of equations which the 56 free variables must satisfy.
KR=PolynomialRing(GF(2),56,'k')
key=list(KR.gens())
for i in range(1,8): #every 8th bit is the parity (standard)
    key.insert(8*i-1,sum([key[8*i-(j+2)] for j in range(0,7)])+1)
expkey=KeySchedule(key)
vli=[vector(ek) for ek in expkey]
eqs=[]
for i in range(1,16):
    eqs+=list(vli[0]-vli[i])
    
I=KR.ideal(eqs)
bas=I.groebner_basis()
kgens=list(KR.gens())
m=matrix(GF(2),[[b.coefficient(k) for k in kgens] for b in bas])
ker=m.transpose().kernel()
assert ker.dimension()==2 #so we have four keys

print(f"Kernel dimension: {ker.dimension()}")
print(f"Number of weak keys: 2^{ker.dimension()} = {2**ker.dimension()}")
Kernel dimension: 2
Number of weak keys: 2^2 = 4

Weak Key Solutions#

kerbas=ker.basis()
k1=0*kerbas[0]
k2=kerbas[0]
k3=kerbas[1]
k4=k2+k3
def KeyToHex(key):
    size=len(key)//4
    hlist=[]
    for i in range(0,size):
        kpart=key[4*i:4*i+4]
        hlist.append(hex(int(''.join(map(str,kpart)),base=2))[2:].upper())
    return ''.join(hlist)
def HexToKey(h):
    key=[]
    for hel in h:
        key+=map(int,bin(int(h,base=16))[2:])
    return key
k=k1
key=GenerateKey(list(k))
print("Initial key")
print(KeyToHex(key))
roundkeys1=KeySchedule(key)
print("\nRound key")
for el in roundkeys1:
    #print(''.join(map(str,el)))
    print(KeyToHex(el))
Initial key
0101010101010101

Round key
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
000000000000
k=k2
key=GenerateKey(list(k))
print("Initial key")
print(KeyToHex(key))
roundkeys1=KeySchedule(key)
print("\nRound key")
for el in roundkeys1:
    #print(''.join(map(str,el)))
    print(KeyToHex(el))
Initial key
E0E0E0E0F1F1F1F1

Round key
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
FFFFFF000000
k=k3
key=GenerateKey(list(k))
print("Initial key")
print(KeyToHex(key))
roundkeys1=KeySchedule(key)
print("\nRound key")
for el in roundkeys1:
    #print(''.join(map(str,el)))
    print(KeyToHex(el))
Initial key
1F1F1F1F0E0E0E0E

Round key
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
000000FFFFFF
k=k4
key=GenerateKey(list(k))
print("Initial key")
print(KeyToHex(key))
roundkeys1=KeySchedule(key)
print("\nRound key")
for el in roundkeys1:
    #print(''.join(map(str,el)))
    print(KeyToHex(el))
Initial key
FEFEFEFEFEFEFEFE

Round key
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF
FFFFFFFFFFFF

Weak Key Self-Inversion#

A property of a weak key is that the encryption procedure is self-invertible: \(E_{K_{\text{weak}}}(E_{K_{\text{weak}}}(P)) = P\).

weak_key=list(map(int,GenerateKey(list(k1))))

for _ in range(0,1000):
    data=[choice([0,1]) for _ in range(0,64)]
    cipher=DESEncryption(data,weak_key)
    plain=DESEncryption(cipher,weak_key)
    assert plain==data
    
print("All assertions are satisfied - a weak key provides a self-invertible procedure")
All assertions are satisfied - a weak key provides a self-invertible procedure

Experiment 3: Complementation Property#

DES satisfies \(E_K(P) = \overline{E_{\overline{K}}(\overline{P})}\), which halves the effective search space for brute force.

This property allows a small reduction in the algorithm for finding the encryption key. Knowing two pairs \((P_1,C_1)\) and \((\overline{P}_{1},C_2)\), it follows that

\[ \overline{C}_{2} = E_{\overline{K}}(P_1)\]

so checking whether a given key \(K\) with message \(P_1\) generates \(C_1\) or \(\overline{C}_{2}\) eliminates two keys at once: \(K\) and \(\overline{K}\).

n_tests = 50
all_pass = True

for _ in range(n_tests):
    key = GenerateKey([choice([0,1]) for _ in range(56)])
    pt = [choice([0,1]) for _ in range(64)]
    
    key_comp = [1 - b for b in key]   # complement the key
    pt_comp = [1 - b for b in pt]     # complement the plaintext
    
    ct1 = DESEncryption(pt, key)
    ct2 = DESEncryption(pt_comp, key_comp)
    ct2_comp = [1 - b for b in ct2]   # complement of ct2
    
    if ct1 != ct2_comp:
        all_pass = False
        break

print(f"Complementation property E_K(P) = complement(E_comp(K)(comp(P)))")
print(f"Tested {n_tests} random cases: {'ALL PASSED' if all_pass else 'FAILED'}")
print()
print("Implication: brute force search needs only 2^55 encryptions (not 2^56)")
print("because for each key K tried, we get the result for complement(K) for free.")
Complementation property E_K(P) = complement(E_comp(K)(comp(P)))
Tested 50 random cases: ALL PASSED

Implication: brute force search needs only 2^55 encryptions (not 2^56)
because for each key K tried, we get the result for complement(K) for free.

Experiment 4: Diffusion Across Rounds#

We visualize the fraction of differing bits as a function of the number of rounds using the original SageMath list_plot. This is the same avalanche measurement from Experiment 1, shown with the user’s original code.

list_plot([(i,sum([ChangeInfluence(0,i) 
for _ in range(0,50)])/50) 
for i in range(2,18)])
../_images/78fc8bdc6d2e9da2ac4fb863a8f0f234844c9c4de88fcb9e76a23e9553b3f105.png

Observation

After round 1, the Feistel structure ensures only half the bits are affected (the right half influences the left half via \(F\)). By round 5, the plot shows near-uniform influence — every input bit affects every output bit with probability close to 0.5. This is the diffusion property that Shannon identified as essential for security.

15.6 Exercises#

Exercise 15.1 (Warm-up) Verify the Feistel invertibility theorem by implementing a generic 4-round Feistel cipher with a random round function \(F\) and confirming that decryption recovers the plaintext.

Exercise 15.2 (Applied — Algebraic) DES has 4 weak keys and 6 pairs of semi-weak keys (12 keys total). A semi-weak pair \((K_1, K_2)\) satisfies \(E_{K_1}(E_{K_2}(P)) = P\).

Use the Gröbner basis approach from Experiment 2: instead of requiring all round keys to be equal, require that the round key sequence of \(K_1\) is the reverse of that of \(K_2\). Find all 6 semi-weak key pairs.

Exercise 15.3 (Analysis) Show that 2DES (double DES) with two independent 56-bit keys provides only \(2^{57}\) effective security, not \(2^{112}\), due to a meet-in-the-middle attack.

Exercise 15.4 (Theory) Prove the complementation property: \(E_K(P) = \overline{E_{\bar{K}}(\bar{P})}\).

Exercise 15.5 (Challenge) Implement a reduced DES with only 6 rounds using DESEncryptionParam. Generate 100 random plaintext-ciphertext pairs and verify that the avalanche is incomplete (some bit correlations remain). What is the minimum number of rounds for full avalanche?

15.7 Summary#

Concept

Key Point

Structure

16-round Feistel network on 64-bit blocks

Key length

56 bits effective (plus 8 parity bits)

Round function

Expansion → XOR → 8 S-boxes → P-permutation

Feistel property

Always invertible, regardless of round function

Weak keys

4 keys where all round keys are identical (found via Gröbner bases)

Complementation

\(E_K(P) = \overline{E_{\bar{K}}(\bar{P})}\) halves search space

Avalanche

Full diffusion by round 5–6

SageMath tools

SymmetricGroup, PolynomialRing(GF(2)), ideal(), groebner_basis(), matrix(), kernel()

Tip

DES is now the standard target for both linear and differential cryptanalysis. In Chapters 16–18, we will use linear cryptanalysis (Matsui 1993) to attack an SPN cipher modeled on DES. In Chapters 19–21, differential cryptanalysis (Biham & Shamir 1991) will provide an even more powerful attack.

References#

  1. Feistel, H. (1973). Cryptography and Computer Privacy. Scientific American, 228(5), 15–23.

  2. National Bureau of Standards (1977). Data Encryption Standard. FIPS Publication 46.

  3. Biham, E. and Shamir, A. (1991). Differential Cryptanalysis of DES-like Cryptosystems. Journal of Cryptology, 4(1), 3–72. [9]

  4. Matsui, M. (1994). Linear Cryptanalysis Method for DES Cipher. EUROCRYPT ‘93, LNCS 765, 386–397. [10]

  5. Heys, H. M. (2002). A Tutorial on Linear and Differential Cryptanalysis. Cryptologia, 26(3), 189–221. [11]