HMAC-SHA-256 implementation in Python 3

A hash-based message authentication code (HMAC) is an algorithm for generating a message authentication code (MAC), which can be used to verify both the integrity and the authentication of a given message. Although both constructs, HMAC and MAC, are based on a cryptographic hash function (such as SHA-1, Whirlpool or RIPEMD-160), the former requires a key (shared between the sender and the receiver of the message) while the latter doesn’t. The HMAC concept was proposed by Bellare, Canetti, and Krawczyk in 1996 and is described in RFC 2104.

As seen from its name, HMAC-SHA-256 uses as its engine the SHA-256 cryptographic hash function, which produces message digests of 256 bits in length. Like the other members of the SHA-2 family (and also MD-5 and SHA-1), SHA-256 is an iterative hash function (based on the Merkle–Damgård scheme) that works by breaking up the input message into blocks of a fixed size (512 bits for SHA-256) and iterating over them with a compression function.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: HMAC-SHA256 implementation in Python 3
#
# Date: 2013-06-10
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#================================================================

from functools import reduce
from math import log,ceil

def intToList2(number,length):
    """Convert a number into a byte list
       with specified length"""
    return [(number >> i) & 0xff
            for i in reversed(range(0,length*8,8))]

def intToList(number):
    """Converts an integer of any length into an integer list"""
    L1 = log(number,256)
    L2 = ceil(L1)
    if L1 == L2:
        L2 += 1
    return [(number&(0xff<<8*i))>>8*i for i in reversed(range(L2))] 

def listToInt(lst):
    """Convert a byte list into a number"""
    return reduce(lambda x,y:(x<<8)+y,lst)

def bitList32ToList4(lst):
    """Convert a 32-bit list into a 4-byte list"""
    def bitListToInt(lst):
        return reduce(lambda x,y:(x<<1)+y,lst)    
    
    lst2 = []
    for i in range(0,len(lst),8):
        lst2.append(bitListToInt(lst[i:i+8]))
    return list([0]*(4-len(lst2)))+lst2

def list4ToBitList32(lst):
    """Convert a 4-byte list into a 32-bit list"""
    def intToBitList2(number,length):
        """Convert an integer into a bit list
        with specified length"""
        return [(number>>n) & 1
                for n in reversed(range(length))]    
    
    lst2 = []
    for e in lst:
        lst2 += intToBitList2(e,8)
    return list([0]*(32-len(lst2)))+lst2

def add32(p,q,r=None,s=None,t=None):
    """Add up to five 32-bit numbers"""
    mask32 = (1<<32)-1
    p2,q2 = listToInt(p), listToInt(q)
    if t is None:
        if s is None:
            if r is None:
                return intToList2((p2+q2)&mask32,4)
            else:
                r2 = listToInt(r)
                return intToList2((p2+q2+r2)&mask32,4)
        else:
            r2,s2 = listToInt(r),listToInt(s)
            return intToList2((p2+q2+r2+s2)&mask32,4)
    else:
        r2,s2,t2 = listToInt(r),listToInt(s),listToInt(t)
        return intToList2((p2+q2+r2+s2+t2)&mask32,4)

def xor(x,y,z=None):
    """Evaluate the XOR on two or three operands"""
    if z is None:
        return list(i^j for i,j in zip(x,y))
    else:
        return list(i^j^k for i,j,k in zip(x,y,z))
   
def sha256(m):
    """Return the SHA-256 digest of input"""
    def padding(m):
        """Pad message according to SHA-256 rules"""
        def bitListToList(lst):
            """Convert a bit list into a byte list"""
            lst2 = [0]*((8-len(lst)%8)%8)+lst
            return [reduce(lambda x,y:(x<<1)+y,lst2[i*8:i*8+8])
                    for i in range(len(lst2)//8)]
    
        def intToBitList(number):
            """Convert an integer into a bit list"""
            return list(map(int,list(bin(number)[2:])))
    
        if type(m) is int:
            m1 = intToBitList(m)
            L = len(m1)
            k = (447-L)%512
            return bitListToList(m1+[1]+list([0]*k))+intToList2(L,8)
        else:
            m1 = m
            if type(m) is str:
                m1 = list(map(ord,m))
            if not(type(m) is list):
                raise TypeError
            L = len(m1)*8
            k = (447-L)%512
            return m1+bitListToList([1]+list([0]*k))+intToList2(L,8)    
    
    def compress(m):
        """Evaluates SHA-256 compression function to input"""
        def Ch(x,y,z):
            return list([(i&j)^((i^0xff)&k) for i,j,k in zip(x,y,z)])
        
        def Maj(x,y,z):
            return list([(i&j)^(i&k)^(j&k) for i,j,k in zip(x,y,z)])              
        
        def rotRight(p,n):
            """Rotate 32-bit word right by n bits"""
            p2 = list4ToBitList32(p)
            return bitList32ToList4(p2[-n:]+p2[:-n])
        
        def shiftRight(p,n):
            """Shift 32-bit right by n bits"""
            p2 = list4ToBitList32(p)
            return bitList32ToList4(list(bytes(n))+p2[:-n])
        
        def Sigma0(p):
            """SHA-256 function"""
            return xor(rotRight(p,2),rotRight(p,13),rotRight(p,22))
        
        def Sigma1(p):
            """SHA-256 function"""
            return xor(rotRight(p,6),rotRight(p,11),rotRight(p,25))
        
        def sigma0(p):
            """SHA-256 function"""
            return xor(rotRight(p,7),rotRight(p,18),shiftRight(p,3))
        
        def sigma1(p):
            """SHA-256 function"""
            return xor(rotRight(p,17),rotRight(p,19),shiftRight(p,10))            
        
        nonlocal H
        [a,b,c,d,e,f,g,h] = H
        K = [0x428a2f98, 0x71374491, 0xb5c0fbcf, 0xe9b5dba5,
             0x3956c25b, 0x59f111f1, 0x923f82a4, 0xab1c5ed5,
             0xd807aa98, 0x12835b01, 0x243185be, 0x550c7dc3,
             0x72be5d74, 0x80deb1fe, 0x9bdc06a7, 0xc19bf174,
             0xe49b69c1, 0xefbe4786, 0x0fc19dc6, 0x240ca1cc,
             0x2de92c6f, 0x4a7484aa, 0x5cb0a9dc, 0x76f988da,
             0x983e5152, 0xa831c66d, 0xb00327c8, 0xbf597fc7,
             0xc6e00bf3, 0xd5a79147, 0x06ca6351, 0x14292967,
             0x27b70a85, 0x2e1b2138, 0x4d2c6dfc, 0x53380d13,
             0x650a7354, 0x766a0abb, 0x81c2c92e, 0x92722c85,
             0xa2bfe8a1, 0xa81a664b, 0xc24b8b70, 0xc76c51a3,
             0xd192e819, 0xd6990624, 0xf40e3585, 0x106aa070,
             0x19a4c116, 0x1e376c08, 0x2748774c, 0x34b0bcb5,
             0x391c0cb3, 0x4ed8aa4a, 0x5b9cca4f, 0x682e6ff3,
             0x748f82ee, 0x78a5636f, 0x84c87814, 0x8cc70208,
             0x90befffa, 0xa4506ceb, 0xbef9a3f7, 0xc67178f2]    
        W = [None]*64
        for t in range(16):
            W[t] = m[t*4:t*4+4]
        for t in range(16,64):
            W[t] = add32(sigma1(W[t-2]),W[t-7],sigma0(W[t-15]),W[t-16])
        for t in range(64):
            T1 = add32(h,Sigma1(e),Ch(e,f,g),intToList2(K[t],4),W[t])
            T2 = add32(Sigma0(a),Maj(a,b,c))
            h = g; g = f; f = e; e = add32(d,T1)
            d = c; c = b; b = a; a = add32(T1,T2)
        H = [add32(x,y) for x,y in zip([a,b,c,d,e,f,g,h],H)]
        
    H0 = [0x6a09e667, 0xbb67ae85, 0x3c6ef372, 0xa54ff53a,
         0x510e527f, 0x9b05688c, 0x1f83d9ab, 0x5be0cd19]
    H = list(map(lambda x:intToList2(x,4),H0))
    mp = padding(m)
    for i in range(0,len(mp),64):
        compress(mp[i:i+64])
    return listToInt([s2 for s1 in H for s2 in s1])

def hmac_sha256(k,m):
    """Return the HMAC-SHA-256 of the input
       HMAC(k,m)=SHA-256((k⊕opad)∥SHA-256((k⊕ipad)∥m))"""
    opad = list([0x5c]*64); ipad = list([0x36]*64)
    if type(k) is int:
        k1 = intToList(k)
        L = len(k1)
        if L > 64:
            K = intToList2(sha256(k),32)+list([0]*32)
        else:
            K = k1+list([0]*(64-L))  
    else:
        k1 = list(map(ord,k))      
        L = len(k1)
        if L > 64:
            K = intToList(sha256(k1))
        else:
            K = k1+list([0]*(64-L))
    if type(m) is int:
        M = intToList(m)
    else:
        M = list(map(ord,m))
    arg1 = xor(K,opad)
    arg2 = xor(K,ipad)
    return sha256(arg1+intToList(sha256(arg2+M)))

if __name__ == '__main__':

    # Wikipedia's test case #1
    assert hmac_sha256("","") == 0xb613679a0814d9ec772f95d778c35fc5ff1697c493715653c6c712144292c5ad
    
    # Wikipedia's test case #2
    assert hmac_sha256("key", "The quick brown fox jumps over the lazy dog") == \
           0xf7bc83f430538424b13298e6aa6fb143ef4d59a14946175997479dbc2d1a3cd8
    
    # RFC 4231 - Identifiers and Test Vectors for HMAC-SHA-224, HMAC-SHA-256,
    # HMAC-SHA-384, and HMAC-SHA-512
    
    # RFC 4231 Test case 1
    Key1 = 0x0b0b0b0b0b0b0b0b0b0b0b0b0b0b0b0b0b0b0b0b
    Data1 = 0x4869205468657265
    HMAC1 = 0xb0344c61d8db38535ca8afceaf0bf12b881dc200c9833da726e9376c2e32cff7    
    assert hmac_sha256(Key1,Data1) == HMAC1      

    # RFC 4231 Test case 2
    Key2 = 0x4a656665
    Data2 = 0x7768617420646f2079612077616e7420666f72206e6f7468696e673f
    HMAC2 = 0x5bdcc146bf60754e6a042426089575c75a003f089d2739839dec58b964ec3843
    assert hmac_sha256(Key2,Data2) == HMAC2
    
    # RFC 4231 Test case 3    
    Key3 = 0xaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa
    Data3 = 0xdddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddddd
    HMAC3 = 0x773ea91e36800e46854db8ebd09181a72959098b3ef8c122d9635514ced565fe
    assert hmac_sha256(Key3,Data3) == HMAC3

    # RFC 4231 Test case 4
    Key4 = 0x0102030405060708090a0b0c0d0e0f10111213141516171819
    Data4 = 0xcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcdcd
    HMAC4 = 0x82558a389a443c0ea4cc819899f2083a85f0faa3e578f8077a2e3ff46729665b
    assert hmac_sha256(Key4,Data4) == HMAC4

    # RFC 4231 Test case 5
    Key5 = 0x0c0c0c0c0c0c0c0c0c0c0c0c0c0c0c0c0c0c0c0c
    Data5 = 0x546573742057697468205472756e636174696f6e 
    HMAC5 = 0xa3b6167473100ee06e0c796c2955552b
    assert hmac_sha256(Key5,Data5)>>128 == HMAC5

    # RFC 4231 Test case 6
    Key6 = 0xaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa
    Data6 = 0x54657374205573696e67204c6172676572205468616e20426c6f636b2d53697a65204b6579202d2048617368204b6579204669727374
    HMAC6 = 0x60e431591ee0b67f0d8a26aacbf5b77f8e0bc6213728c5140546040f0ee37f54
    assert hmac_sha256(Key6,Data6) == HMAC6

    # RFC 4231 Test case 7    
    Key7 = 0xaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa
    Data7 = 0x5468697320697320612074657374207573696e672061206c6172676572207468616e20626c6f636b2d73697a65206b657920616e642061206c6172676572207468616e20626c6f636b2d73697a6520646174612e20546865206b6579206e6565647320746f20626520686173686564206265666f7265206265696e6720757365642062792074686520484d414320616c676f726974686d2e
    HMAC7 = 0x9b09ffa71b942fcb27635fbcd5b0e944bfdc63644f0713938a7f51535c3a35e2
    assert hmac_sha256(Key7,Data7) == HMAC7

    print("Ok!")

Validation of an AES implementation in Python 3

The Cryptographic Algorithm Validation Program (CAVP) defines validation testing for cryptographic algorithms approved by the US National Institute of Standards and Technology (NIST). All of the tests under CAVP, established by NIST and its Canadian counterpart (CSEC) in 1995, are handled by accredited third-party laboratories.

To assist prospective vendors in checking their implementations, NIST provides electronic versions of the vectors for the Known Answer Test (KAT) for the three NIST-approved symmetric cryptographic algorithms: AES, Triple-DES, and Skipjack. Also available are sample values for the Monte Carlo (MCT) test and the Multiblock Message (MMT) test for the same algorithms, thus completing the set of tests a cryptographic implementation (dubbed “Implementation Under Test”) will face during a formal validation. A detailed account of the procedures involved in validating AES implementations can be found in the NIST document The Advanced Encryption Standard Algorithm Validation Suite (AESAVS).

The NIST KAT validation suite for AES contains 72 files describing test vectors for different AES modes of operation: ECB (Electronic Codebook), CBC (Cipher Block Chaining), CFB (Cipher Feedback) and OFB (Output Feedback). Besides this partition, there are also separated tests for AES encryption and decryption, exemplified by the abridged contents of the file ECBKeySbox128.rsp (see below).

# CAVS 11.1
# Config info for aes_values
# AESVS KeySbox test data for ECB
# State : Encrypt and Decrypt
# Key Length : 128
# Generated on Fri Apr 22 15:11:26 2011

[ENCRYPT]

COUNT = 0
KEY = 10a58869d74be5a374cf867cfb473859
PLAINTEXT = 00000000000000000000000000000000
CIPHERTEXT = 6d251e6944b051e04eaa6fb4dbf78465

COUNT = 1
KEY = caea65cdbb75e9169ecd22ebe6e54675
PLAINTEXT = 00000000000000000000000000000000
CIPHERTEXT = 6e29201190152df4ee058139def610bb

... (18 test cases omitted)

[DECRYPT]

COUNT = 0
KEY = 10a58869d74be5a374cf867cfb473859
CIPHERTEXT = 6d251e6944b051e04eaa6fb4dbf78465
PLAINTEXT = 00000000000000000000000000000000

COUNT = 1
KEY = caea65cdbb75e9169ecd22ebe6e54675
CIPHERTEXT = 6e29201190152df4ee058139def610bb
PLAINTEXT = 00000000000000000000000000000000

... (18 test cases omitted)

These vectors can be used to informally verify the correctness of an AES implementation, such as the one presented below, which successfully passed all 4,156 KAT tests involving ECB and CBC modes. This Python code differs from the one already presented in a previous blog only with respect to the input/output types accepted: if the input plaintext (ciphertext) is of integer type, so will be the correspondent ciphertext (plaintext) output. This feature simplifies validation testing, since the integer plaintexts and ciphertexts usually employed in those scenarios can be deal with directly without any type conversion.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: AES implementation in Python 3
#              (sundAES)
#
# Date: 2013-06-02 (version 1.1)
#       2012-01-16 (version 1.0)
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================
import sys
from itertools import repeat
from functools import reduce
from copy import copy

__all__ = ["setKey","encrypt","decrypt"]

def memoize(func):
    """Memoization function"""
    memo = {}
    def helper(x):
        if x not in memo:
            memo[x] = func(x)
        return memo[x]
    return helper

def mult(p1,p2):
    """Multiply two polynomials in GF(2^8)/x^8+x^4+x^3+x+1"""
    p = 0
    while p2:
        if p2&0x01:
            p ^= p1
        p1 <<= 1
        if p1&0x100:
            p1 ^= 0x1b
        p2 >>= 1
    return p&0xff

# Auxiliary one-parameter functions defined for memoization
# (to speed up multiplication in GF(2^8))

@memoize
def x2(y):
    """Multiplication by 2"""
    return mult(2,y)

@memoize
def x3(y):
    """Multiplication by 3"""
    return mult(3,y)

@memoize
def x9(y):
    """Multiplication by 9"""
    return mult(9,y)

@memoize
def x11(y):
    """Multiplication by 11"""
    return mult(11,y)

@memoize
def x13(y):
    """Multiplication by 13"""
    return mult(13,y)

@memoize
def x14(y):
    """Multiplication by 14"""
    return mult(14,y)

class AES:
    """Class definition for AES objects"""
    keySizeTable = {"SIZE_128":16,
                    "SIZE_192":24,
                    "SIZE_256":32}
    wordSizeTable = {"SIZE_128":44,
                     "SIZE_192":52,
                     "SIZE_256":60}
    numberOfRoundsTable = {"SIZE_128":10,
                           "SIZE_192":12,
                           "SIZE_256":14}
    cipherModeTable = {"MODE_ECB":1,
                       "MODE_CBC":2}
    paddingTable = {"NoPadding":0,
                    "PKCS7Padding":1}
    # S-Box
    sBox = (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)
    # Inverse S-Box
    invSBox = (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)

    # Instance variables
    wordSize = None
    w = [None]*60 # Round subkeys list
    keyDefined = None # Key definition flag
    numberOfRounds = None
    cipherMode = None
    padding = None # Padding scheme
    ivEncrypt = None # Initialization
    ivDecrypt = None #  vectors

    def __init__(self,mode,padding = "NoPadding"):
        """Create a new instance of an AES object"""
        try:
            assert mode in AES.cipherModeTable
        except AssertionError:
            print("Cipher mode not supported:",mode)
            sys.exit("ValueError")
        self.cipherMode = mode
        try:
            assert padding in AES.paddingTable
        except AssertionError:
            print("Padding scheme not supported:",padding)
            sys.exit(ValueError)
        self.padding = padding
        self.keyDefined = False

    def intToList(self,number):
        """Convert an 16-byte number into a 16-element list"""
        return [(number>>i)&0xff for i in reversed(range(0,128,8))]

    def intToList2(self,number):
        """Converts an integer into one (or more) 16-element list"""
        lst = []
        while number:
            lst.append(number&0xff)
            number >>= 8
        m = len(lst)%16
        if m == 0 and len(lst) != 0:
            return lst[::-1]
        else:
            return list(bytes(16-m)) + lst[::-1]

    def listToInt(self,lst):
        """Convert a list into a number"""
        return reduce(lambda x,y:(x<<8)+y,lst)

    def wordToState(self,wordList):
        """Convert list of 4 words into a 16-element state list"""
        return [(wordList[i]>>j)&0xff
                for j in reversed(range(0,32,8)) for i in range(4)]

    def listToState(self,list):
        """Convert a 16-element list into a 16-element state list"""
        return [list[i+j] for j in range(4) for i in range(0,16,4)]

    stateToList = listToState # this function is an involution

    def subBytes(self,state):
        """SubBytes transformation"""
        return [AES.sBox[e] for e in state]

    def invSubBytes(self,state):
        """Inverse SubBytes transformation"""
        return [AES.invSBox[e] for e in state]

    def shiftRows(self,s):
        """ShiftRows transformation"""
        return s[:4]+s[5:8]+s[4:5]+s[10:12]+s[8:10]+s[15:]+s[12:15]

    def invShiftRows(self,s):
        """Inverse ShiftRows transformation"""
        return s[:4]+s[7:8]+s[4:7]+s[10:12]+s[8:10]+s[13:]+s[12:13]

    def mixColumns(self,s):
        """MixColumns transformation"""
        return [x2(s[i])^x3(s[i+4])^   s[i+8] ^   s[i+12]  for i in range(4)]+ \
               [   s[i] ^x2(s[i+4])^x3(s[i+8])^   s[i+12]  for i in range(4)]+ \
               [   s[i] ^   s[i+4] ^x2(s[i+8])^x3(s[i+12]) for i in range(4)]+ \
               [x3(s[i])^   s[i+4] ^   s[i+8] ^x2(s[i+12]) for i in range(4)]

    def invMixColumns(self,s):
        """Inverse MixColumns transformation"""
        return [x14(s[i])^x11(s[i+4])^x13(s[i+8])^ x9(s[i+12]) for i in range(4)]+ \
               [ x9(s[i])^x14(s[i+4])^x11(s[i+8])^x13(s[i+12]) for i in range(4)]+ \
               [x13(s[i])^ x9(s[i+4])^x14(s[i+8])^x11(s[i+12]) for i in range(4)]+ \
               [x11(s[i])^x13(s[i+4])^ x9(s[i+8])^x14(s[i+12]) for i in range(4)]

    def addRoundKey (self,subkey,state):
        """AddRoundKey transformation"""
        return [i^j for i,j in zip(subkey,state)]

    xorLists = addRoundKey

    def rotWord(self,number):
        """Rotate subkey left"""
        return (((number&0xff000000)>>24) +
                ((number&0xff0000)<<8) +
                ((number&0xff00)<<8) +
                ((number&0xff)<<8))

    def subWord(self,key):
        """Substitute subkeys bytes using S-box"""
        return ((AES.sBox[(key>>24)&0xff]<<24) +
                (AES.sBox[(key>>16)&0xff]<<16) +
                (AES.sBox[(key>>8)&0xff]<<8) +
                 AES.sBox[key&0xff])

    def setKey(self,keySize,key,iv = None):
        """KeyExpansion transformation"""
        rcon = (0x00,0x01,0x02,0x04,0x08,0x10,0x20,0x40,0x80,0x1B,0x36)
        try:
            assert keySize in AES.keySizeTable
        except AssertionError:
            print("Key size identifier not valid")
            sys.exit("ValueError")
        try:
            assert isinstance(key,int)
        except AssertionError:
            print("Invalid key")
            sys.exit("ValueError")
        klen = len("{:02x}".format(key))//2
        try:
            assert klen <= AES.keySizeTable[keySize]
        except AssertionError:
            print("Key size mismatch")
            sys.exit("ValueError")
        try:
            assert ((self.cipherMode == "MODE_CBC" and isinstance(iv,int)) or
                     self.cipherMode == "MODE_ECB")
        except AssertionError:
                print("IV is mandatory for CBC mode")
                sys.exit(ValueError)

        if self.cipherMode == "MODE_CBC":
            temp = self.intToList(iv)
            self.ivEncrypt = copy(temp)
            self.ivDecrypt = copy(temp)
        nr = AES.numberOfRoundsTable[keySize]
        self.numberOfRounds = nr
        self.wordSize = AES.wordSizeTable[keySize]
        if nr == 10:
            nk = 4
            keyList = self.intToList(key)
        elif nr == 12:
            nk = 6
            keyList =  self.intToList(key>>64) + \
                      (self.intToList(key&int("ff"*32,16)))[8:]
        else:
            nk = 8
            keyList =  self.intToList(key>>128) + \
                       self.intToList(key&int("ff"*64,16))
        for index in range(nk):
            self.w[index] =  (keyList[4*index]<<24) + \
                             (keyList[4*index+1]<<16) + \
                             (keyList[4*index+2]<<8) +\
                              keyList[4*index+3]
        for index in range(nk,self.wordSize):
            temp = self.w[index - 1]
            if index % nk == 0:
                temp = (self.subWord(self.rotWord(temp)) ^
                        rcon[index//nk]<<24)
            elif self.numberOfRounds == 14 and index%nk == 4:
                temp = self.subWord(temp)
            self.w[index] = self.w[index-nk]^temp
        self.keyDefined = True
        return

    def getKey(self,operation):
        """Return next round subkey for encryption or decryption"""
        if operation == "encryption":
            for i in range(0,self.wordSize,4):
                yield self.wordToState(self.w[i:i+4])
        else: # operation = "decryption":
            for i in reversed(range(0,self.wordSize,4)):
                yield self.wordToState(self.w[i:i+4])

    def encryptBlock(self,plaintextBlock):
        """Encrypt a 16-byte block with key already defined"""
        key = self.getKey("encryption")
        state = self.listToState(plaintextBlock)
        state = self.addRoundKey(next(key),state)
        for _ in repeat(None,self.numberOfRounds - 1):
            state = self.subBytes(state)
            state = self.shiftRows(state)
            state = self.mixColumns(state)
            state = self.addRoundKey(next(key),state)
        state = self.subBytes(state)
        state = self.shiftRows(state)
        state = self.addRoundKey(next(key),state)
        return self.stateToList(state)

    def decryptBlock(self,ciphertextBlock):
        """Decrypt a 16-byte block with key already defined"""
        key = self.getKey("decryption")
        state = self.listToState(ciphertextBlock)
        state = self.addRoundKey(next(key),state)
        for _ in repeat(None,self.numberOfRounds - 1):
            state = self.invShiftRows(state)
            state = self.invSubBytes(state)
            state = self.addRoundKey(next(key),state)
            state = self.invMixColumns(state)
        state = self.invShiftRows(state)
        state = self.invSubBytes(state)
        state = self.addRoundKey(next(key),state)
        return self.stateToList(state)

    def padData(self,data):
        """Add PKCS7 padding to plaintext (or just add bytes to fill a block)"""
        paddingLength = 16-(len(data)%16)
        if self.padding == "NoPadding":
            paddingLength %= 16
        if type(data) is bytes:
            return data+bytes(list([paddingLength]*paddingLength))
        else:
            return [ord(s) for s in data]+[paddingLength]*paddingLength

    def unpadData(self,byteList):
        """Remove PKCS7 padding (if present) from plaintext"""
        if self.padding == "PKCS7Padding":
            return "".join(chr(e) for e in byteList[:-byteList[-1]])
        else:
            return "".join(chr(e) for e in byteList)

    def encrypt(self,input):
        """Encrypt plaintext passed as a string or as an integer"""
        try:
            assert self.keyDefined
        except AssertionError:
            print("Key not defined")
            sys.exit("ValueError")

        if type(input) is int:
            inList = self.intToList2(input)
        else:
            inList = self.padData(input)
        outList = []
        if self.cipherMode == "MODE_CBC":
            outBlock = self.ivEncrypt
            for i in range(0,len(inList),16):
                auxList = self.xorLists(outBlock,inList[i:i+16])
                outBlock = self.encryptBlock(auxList)
                outList += outBlock
            self.ivEncrypt = outBlock
        else:
            for i in range(0,len(inList),16):
                outList += self.encryptBlock(inList[i:i+16])
        if type(input) is int:
            return self.listToInt(outList)
        else:
            return outList

    def decrypt(self,input):
        """Decrypt ciphertext passed as a string or as an integer"""
        try:
            assert self.keyDefined
        except AssertionError:
            print("Key not defined")
            sys.exit("ValueError")
        if type(input) is int:
            inList = self.intToList2(input)
        else:
            inList = input
        outList = []
        if self.cipherMode == "MODE_CBC":
            oldInBlock = self.ivDecrypt
            for i in range(0,len(inList),16):
                newInBlock = inList[i:i+16]
                auxList = self.decryptBlock(newInBlock)
                outList += self.xorLists(oldInBlock,auxList)
                oldInBlock = newInBlock
            self.ivDecrypt = oldInBlock
        else:
            for i in range(0,len(inList),16):
                outList += self.decryptBlock(inList[i:i+16])
        if type(input) is int:
            return self.listToInt(outList)
        else:
            return self.unpadData(outList)

The Python program below extracts data from the files contained in the KAT_AES directory, then builds the test cases and finally executes them. Only the files applicable to ECB and CBC operation modes (the modes supported by this AES implementation) are considered.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: Validation of an AES implementation in Python
#
# Date: 2013-06-02
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================

import os,sys,re
from functools import reduce
from glob import glob
import sundAES

# Global counters
noFilesTested = noFilesSkipped = 0
counterOK = counterNOK = 0

class AEStester:
    """"""
    def buildTestCases(self,filename):
        """Build test cases described in a given file"""
        global noFilesTested,noFilesSkipped

        self.basename = os.path.basename(filename)
        if self.basename.startswith('ECB'):
            self.mode = "MODE_ECB"
            noFilesTested += 1
        elif self.basename.startswith('CBC'):
            self.mode = "MODE_CBC"
            noFilesTested += 1
        else:
            noFilesSkipped += 1
            return

        digits = re.search("\d{3}",self.basename)
        self.keysize = 'SIZE_' + digits.group()
        result = re.search("CFB\d*(\D{6,})\d{3}",self.basename)
        if result != None:
            self.typeTest = result.group(1)
        self.typeTest = re.search("\w{3}(\D{6,})\d{3}",self.basename).group(1)
        self.iv = None
        for line in open(filename):
            line = line.strip()
            if (line == "") or line.startswith('#'):
                continue
            elif line == '[ENCRYPT]':
                self.operation = 'encrypt'
                continue
            elif line == '[DECRYPT]':
                self.operation = 'decrypt'
                continue
            param,_,value = line.split(' ',2)
            if param == "COUNT":
                self.count = int(value)
                continue
            else:
                self.__setattr__(param.lower(),int(value,16))
                if (self.operation == 'encrypt') and (param == "CIPHERTEXT") or \
                   (self.operation == 'decrypt' and param == "PLAINTEXT"):
                    self.runTestCase()

    def runTestCase(self):
        """Execute test case and report result"""
        global counterOK,counterNOK

        def printTestCase(result):
            print("Type={0:s} Mode={1:s} Keysize={2:s} Function={3:s} Count={4:03d} {5:s}"\
                  .format(self.typeTest,self.mode[5:],\
                          self.keysize[5:],self.operation.upper(),\
                          self.count,result))

        obj = sundAES3.AES(self.mode)
        obj.setKey(self.keysize,self.key,self.iv)
        if self.operation == 'encrypt':
            CIPHERTEXT = obj.encrypt(self.plaintext)
            try:
                assert self.ciphertext == CIPHERTEXT
                counterOK += 1
                printTestCase("OK")
            except AssertionError:
                counterNOK +=1
                print(self.basename)
                printTestCase("failed")
                print("Expected ciphertext={0:0x}".format(self.ciphertext))
                print("Returned ciphertext={0:0x}".format(CIPHERTEXT))
        else:
            PLAINTEXT = obj.decrypt(self.ciphertext)
            try:
                assert self.plaintext == PLAINTEXT
                counterOK += 1
                printTestCase("OK")
            except AssertionError:
                counterNOK +=1
                print(self.basename)
                printTestCase("failed")
                print("Expected plaintext={0:0x}".format(self.plaintext))
                print("Returned plaintext={0:0x}".format(PLAINTEXT))

# Main program
path = os.path.dirname(__file__)
files = sys.argv[1:]
if not files:
    files = glob(os.path.join(path,'KAT_AES','*.rsp'))
    files.sort()
for file in files:
    AEStester().buildTestCases(file)
print("Files tested={0:d}".format(noFilesTested))
print("Files skipped={0:d}".format(noFilesSkipped))
print("Test cases OK={0:d}".format(counterOK))
print("Test cases NOK={0:d}".format(counterNOK))

Executed without arguments, this program will apply all ECB and CBC tests described in the files against the AES implementation, producing the following output:

Type=GFSbox Mode=CBC Keysize=128 Function=ENCRYPT Count=000 OK
Type=GFSbox Mode=CBC Keysize=128 Function=ENCRYPT Count=001 OK
Type=GFSbox Mode=CBC Keysize=128 Function=ENCRYPT Count=002 OK
Type=GFSbox Mode=CBC Keysize=128 Function=ENCRYPT Count=003 OK
Type=GFSbox Mode=CBC Keysize=128 Function=ENCRYPT Count=004 OK

... (4,147 lines omitted)

Type=VarTxt Mode=ECB Keysize=256 Function=DECRYPT Count=124 OK
Type=VarTxt Mode=ECB Keysize=256 Function=DECRYPT Count=125 OK
Type=VarTxt Mode=ECB Keysize=256 Function=DECRYPT Count=126 OK
Type=VarTxt Mode=ECB Keysize=256 Function=DECRYPT Count=127 OK
Files tested=24
Files skipped=48
Test cases OK=4156
Test cases NOK=0

If, however, this program is run with one or more files in the KAT_AES directory as argument(s), it will instead build and execute the tests contained in the indicated file(s).

AES-GCM implementation in Python 3

Galois/Counter Mode (GCM) is a mode of operation for symmetric key cryptographic block ciphers that provides authenticated encryption.

Proposed by David McGrew and John Viega in 2005, GCM is suited for high-speed secure computing and communication. Acknowledging this fact, the US National Institute of Standards and Technology (NIST) standardized GCM and its companion algorithm, Galois Message Authentication Code (GMAC), in 2007. The latter is an authentication-only variant of GCM which can be used as an incremental message authentication code (MAC). AES-GCM (the Advanced Encryption Algorithm operating in Galois/Counter Mode) has also been included in NSA Suite B Cryptography.

GCM’s confidentiality service is based on a variation of the Counter mode (CTR) while its authenticity assurance relies on a universal hash function defined over the binary Galois field GF(2^{128})/x^{128}+x^7+x^2+x+1.

GCM is defined for block ciphers with block sizes of 128, 192, and 256 bits (AES uses 128-bit blocks). Both GCM and GMAC can accept initialization vectors (IVs) of arbitrary length (AES and other symmetric ciphers, on the other hand, require IVs to be of the same size as the cipher’s block size).

Its interesting to note that Intel high-end CPUs support a special instruction, PCLMULQDQ, that computes the 128-bit product (the carry-less multiplication) of two 64-bit operands. This instruction can be used as a building block to perform the carry-less multiplication of two 128-bit operands required by GCM. The other step in GCM is reduction modulo the irreducible polynomial x^{128}+x^7+x^2+x+1, an operation that takes advantage of the PSRLD, PSLLD and PSHUFD instructions. Together with Intel’s Advanced Encryption Standard Instructions (AES-NI), these instructions allow GCM to offer authenticated encryption services at very high rates without using hardware-based solutions (e.g. FPGAs or ASICs). According to Intel, “an AES-GCM implementation based on the AES-NI and PCLMULQDQ instructions delivered a 400% throughput performance gain when compared to a non-AES-NI enabled software solution on the same platform.”

The Python code below implements AES-GCM using the AES implementation already presented and supports the three key sizes used by AES (128, 192 and 256 bits). All eighteen test cases proposed by McGrew & Viega were used to validate this implementation. Please note that this code is not of production quality.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: AES-GCM (Galois Counter/Mode) implementation
#              in Python 3
#
# Date: 2013-05-30
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#================================================================

import pyAES
from functools import reduce

def xor(x,y):
    """Returns the exclusive or (xor) between two vectors"""
    return bytes(i^j for i,j in zip(x,y))

def intToList(number,listSize):
    """Convert a number into a byte list""" 
    return [(number >> i) & 0xff
            for i in reversed(range(0,listSize*8,8))]

def listToInt(list):
    """Convert a byte list into a number""" 
    return reduce(lambda x,y:(x<<8)+y,list)

def GHASH (hkey,aad,ctext):
    """GCM's GHASH function"""
    def xorMultH (p,q):
        """Multiply (p^q) by hash key"""
        def multGF2(x,y):
            """Multiply two polynomials in GF(2^m)/g(w)
               g(w) = w^128 + w^7 + w^2 + w + 1
               (operands and result bits reflected)"""      
            (x,y) = map(lambda z:listToInt(list(z)),(x,y))
            z = 0
            while y & ((1<<128)-1):
                if y & (1<<127):
                    z ^= x
                y <<= 1
                if x & 1:
                    x = (x>>1)^(0xe1<<120)
                else:
                    x >>= 1
            return bytes(intToList(z,16))

        return bytes(multGF2(hkey,xor(p,q)))
    
    def gLen(s):
        """Evaluate length of input in bits and returns
           it in the LSB bytes of a 64-bit string"""
        return bytes(intToList(len(s)*8,8))  
  
    x = bytes(16)    
    aadP = aad + bytes((16-len(aad)%16)%16)
    ctextP = ctext + bytes((16-len(ctext)%16)%16)
    for i in range(0,len(aadP),16):
        x = xorMultH(x,aadP[i:i+16])
    for i in range(0,len(ctextP),16):
        x = xorMultH(x,ctextP[i:i+16])
    return xorMultH(x,gLen(aad) + gLen(ctext))

def GCM_crypt(keysize,key,iv,input,aad):
    """GCM's Authenticated Encryption/Decryption Operations"""
    def incr(m):
        """Increment the LSB 32 bits of input counter"""
        n = list(m)
        n12 = bytes(n[:12])
        ctr = listToInt(n[12:])
        if ctr == (1<<32)-1:
            return n12 + bytes(4)
        else:
            return n12 + bytes(intToList(ctr+1,4))
       
    obj = pyAES.AES('MODE_ECB')
    obj.setKey(keysize,key)    
    h = bytes(obj.encrypt(bytes(16)))
    output = bytes()
    L = len(input)
    if len(iv) == 12:
        y0 = bytes(iv) + bytes(b'\x00\x00\x00\x01')
    else:
        y0 = bytes(GHASH(h,bytes(),iv))
    y = y0
    for i in range(0,len(input),16):
        y = incr(y)
        ctextBlock = xor(bytes(obj.encrypt(y)),
                         input[i:min(i+16,L)])
        output += bytes(ctextBlock)
    g = obj.encrypt(y0)
    tag = xor(GHASH(h,aad,output),g)
    return output,tag,g,h
 
def GCM_encrypt(keysize,key,iv,ptext,aad):
    """GCM's Authenticated Encryption Operation"""
    (ctext,tag,g,h) = GCM_crypt(keysize,key,iv,ptext,aad)
    return ctext,xor(GHASH(h,aad,ctext),g)

def GCM_decrypt(keysize,key,iv,ctext,aad,tag):
    """GCM's Authenticated Decryption Operation"""
    (ptext,_,g,h) = GCM_crypt(keysize,key,iv,ctext,aad)
    if tag == xor(GHASH(h,aad,ctext),g):
        return True,ptext
    else:
        return False,None

######################## Testing section ########################

if __name__ == '__main__':

    def printHex(s):
        """Prints a bytes string in hex format"""
        print('{0:#0x}'.format(bytesToInt(s)))
    
    def bytesToInt(b):
        """Converts a bytes string into an integer"""
        return listToInt(list(b))
    
    def listToInt(list):
        """Convert a byte list into a number""" 
        return reduce(lambda x,y:(x<<8)+y,list)
    
    def checkTestVector(id,keysize,key,ptext,aad,iv,ctext,tag):
        def intToBytes(n):
            """Converts an integer into a bytes string"""
            lst = []
            while n:
                lst.append(n&0xff)
                n >>= 8
            return bytes(reversed(lst))

        def convertType(v):
            """Convert input variable type to bytes"""
            if type(v) is int:
                return intToBytes(v)
            elif type(v) is bytes:
                return v
            else:
                return bytes(v,'ISO-8859-1')

        def printOutputs():
            """Prints expected/evaluated tags and
                      expected/evaluated ciphertexts"""
            print("Tag expected:  ", end="")
            printHex(tag)
            print("Tag evaluated: ", end="")
            printHex(TAG)
            print("Ciphertext expected:  ", end="")
            printHex(ctext)
            print("Ciphertext evaluated: ", end="")
            printHex(CTEXT)            
          
        (ptext,aad,iv,ctext,tag) = map(convertType,(ptext,aad,iv,ctext,tag))
        (CTEXT,TAG) = GCM_encrypt(keysize,key,iv,ptext,aad)
        (SUCCESS,PTEXT) = GCM_decrypt(keysize,key,iv,CTEXT,aad,TAG)
    
        try:
            assert SUCCESS & (CTEXT == ctext) & (TAG == tag)
            print("Test case {0:0d} succeeded".format(id))
        except AssertionError:
            print("Test case {0:0d} failed".format(id))
            printOutputs()           
    
    # Test cases extracted from
    # "The Galois/Counter Mode of Operation(GCM)", McGrew & Viega, 2005 (http://goo.gl/DWJPK)
    
    # Some useful constants
    emptyString = bytes() # zero length bit string
    nullBitString128 = bytes(16) # 128-bit null string 
    nullBitString96 = bytes(12) # 96-bit null IV
    
    # Test Case #1
    testcase1 = {'id':1,'keysize':"SIZE_128",'key':0x0,'ptext':emptyString,
                 'aad':emptyString,'iv':nullBitString96,'ctext':emptyString,
                 'tag':0x58e2fccefa7e3061367f1d57a4e7455a}
    
    # Test Case #2
    testcase2 = {'id':2,'keysize':"SIZE_128",'key':0x0,'ptext':nullBitString128,
                 'aad':emptyString,'iv':nullBitString96,'ctext':0x0388dace60b6a392f328c2b971b2fe78,
                 'tag':0xab6e47d42cec13bdf53a67b21257bddf}
    
    # Test Case #3
    testcase3 = {'id':3,'keysize':"SIZE_128",'key':0xfeffe9928665731c6d6a8f9467308308,
                 'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b391aafd255,
                 'aad':emptyString,'iv':0xcafebabefacedbaddecaf888,
                 'ctext':0x42831ec2217774244b7221b784d0d49ce3aa212f2c02a4e035c17e2329aca12e21d514b25466931c7d8f6a5aac84aa051ba30b396a0aac973d58e091473f5985,
                 'tag':0x4d5c2af327cd64a62cf35abd2ba6fab4}
    
    ## Test Case #4
    testcase4 = {'id':4,'keysize':"SIZE_128",'key':0xfeffe9928665731c6d6a8f9467308308,
                 'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                 'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,'iv':0xcafebabefacedbaddecaf888,
                 'ctext':0x42831ec2217774244b7221b784d0d49ce3aa212f2c02a4e035c17e2329aca12e21d514b25466931c7d8f6a5aac84aa051ba30b396a0aac973d58e091,
                 'tag':0x5bc94fbc3221a5db94fae95ae7121a47}

    ## Test Case #5
    testcase5 = {'id':5,'keysize':"SIZE_128",'key':0xfeffe9928665731c6d6a8f9467308308,
                 'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                 'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,'iv':0xcafebabefacedbad,
                 'ctext':0x61353b4c2806934a777ff51fa22a4755699b2a714fcdc6f83766e5f97b6c742373806900e49f24b22b097544d4896b424989b5e1ebac0f07c23f4598,
                 'tag':0x3612d2e79e3b0785561be14aaca2fccb}

    ## Test Case #6
    testcase6 = {'id':6,'keysize':"SIZE_128",'key':0xfeffe9928665731c6d6a8f9467308308,
                 'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                 'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,
                 'iv':0x09313225df88406e555909c5aff5269aa6a7a9538534f7da1e4c303d2a318a728c3c0c95156809539fcf0e2429a6b525416aedbf5a0de6a57a637b39b,
                 'ctext':0x8ce24998625615b603a033aca13fb894be9112a5c3a211a8ba262a3cca7e2ca701e4a9a4fba43c90ccdcb281d48c7c6fd62875d2aca417034c34aee5,
                 'tag':0x619cc5aefffe0bfa462af43c1699d050}
        
    ## Test Case #7
    testcase7 = {'id':7,'keysize':"SIZE_192",'key':0x0,'ptext':emptyString,
                 'aad':emptyString,'iv':nullBitString96,'ctext':emptyString,
                 'tag':0xcd33b28ac773f74ba00ed1f312572435}
   
    ## Test Case #8
    testcase8 = {'id':8,'keysize':"SIZE_192",'key':0x0,'ptext':nullBitString128,
                 'aad':emptyString,'iv':nullBitString96,
                 'ctext':0x98e7247c07f0fe411c267e4384b0f600,
                 'tag':0x2ff58d80033927ab8ef4d4587514f0fb}

    ## Test Case #9
    testcase9 = {'id':9,'keysize':"SIZE_192",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c,
                 'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b391aafd255,
                 'aad':emptyString,'iv':0xcafebabefacedbaddecaf888,
                 'ctext':0x3980ca0b3c00e841eb06fac4872a2757859e1ceaa6efd984628593b40ca1e19c7d773d00c144c525ac619d18c84a3f4718e2448b2fe324d9ccda2710acade256,
                 'tag':0x9924a7c8587336bfb118024db8674a14}
                   
    ## Test Case #10
    testcase10 = {'id':10,'keysize':"SIZE_192",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                  'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,'iv':0xcafebabefacedbaddecaf888,
                  'ctext':0x3980ca0b3c00e841eb06fac4872a2757859e1ceaa6efd984628593b40ca1e19c7d773d00c144c525ac619d18c84a3f4718e2448b2fe324d9ccda2710,
                  'tag':0x2519498e80f1478f37ba55bd6d27618c}
                                    
    ## Test Case #11
    testcase11 = {'id':11,'keysize':"SIZE_192",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                  'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,'iv':0xcafebabefacedbad,
                  'ctext':0x0f10f599ae14a154ed24b36e25324db8c566632ef2bbb34f8347280fc4507057fddc29df9a471f75c66541d4d4dad1c9e93a19a58e8b473fa0f062f7,
                  'tag':0x65dcc57fcf623a24094fcca40d3533f8}
                                                    
    ## Test Case #12
    testcase12 = {'id':12,'keysize':"SIZE_192",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                  'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,
                  'iv':0x9313225df88406e555909c5aff5269aa6a7a9538534f7da1e4c303d2a318a728c3c0c95156809539fcf0e2429a6b525416aedbf5a0de6a57a637b39b,
                  'ctext':0xd27e88681ce3243c4830165a8fdcf9ff1de9a1d8e6b447ef6ef7b79828666e4581e79012af34ddd9e2f037589b292db3e67c036745fa22e7e9b7373b,
                  'tag':0xdcf566ff291c25bbb8568fc3d376a6d9}

    ## Test Case #13
    testcase13 = {'id':13,'keysize':"SIZE_256",'key':0x0,'ptext':emptyString,
                  'aad':emptyString,'iv':nullBitString96,'ctext':emptyString,
                  'tag':0x530f8afbc74536b9a963b4f1c4cb738b}

    ## Test Case #14
    testcase14 = {'id':14,'keysize':"SIZE_256",'key':0x0,'ptext':nullBitString128,
                   'aad':emptyString,'iv':nullBitString96,'ctext':0xcea7403d4d606b6e074ec5d3baf39d18,
                   'tag':0xd0d1c8a799996bf0265b98b5d48ab919}             

    ## Test Case #15
    testcase15 = {'id':15,'keysize':"SIZE_256",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c6d6a8f9467308308,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b391aafd255,
                  'aad':emptyString,'iv':0xcafebabefacedbaddecaf888,
                  'ctext':0x522dc1f099567d07f47f37a32a84427d643a8cdcbfe5c0c97598a2bd2555d1aa8cb08e48590dbb3da7b08b1056828838c5f61e6393ba7a0abcc9f662898015ad,
                  'tag':0xb094dac5d93471bdec1a502270e3cc6c}

    ## Test Case #16
    testcase16 = {'id':16,'keysize':"SIZE_256",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c6d6a8f9467308308,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                  'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,'iv':0xcafebabefacedbaddecaf888,
                  'ctext':0x522dc1f099567d07f47f37a32a84427d643a8cdcbfe5c0c97598a2bd2555d1aa8cb08e48590dbb3da7b08b1056828838c5f61e6393ba7a0abcc9f662,
                  'tag':0x76fc6ece0f4e1768cddf8853bb2d551b}

    ## Test Case #17
    testcase17 = {'id':17,'keysize':"SIZE_256",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c6d6a8f9467308308,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                  'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,'iv':0xcafebabefacedbad,
                  'ctext':0xc3762df1ca787d32ae47c13bf19844cbaf1ae14d0b976afac52ff7d79bba9de0feb582d33934a4f0954cc2363bc73f7862ac430e64abe499f47c9b1f,
                  'tag':0x3a337dbf46a792c45e454913fe2ea8f2} 

    ## Test Case #18
    testcase18 = {'id':18,'keysize':"SIZE_256",'key':0xfeffe9928665731c6d6a8f9467308308feffe9928665731c6d6a8f9467308308,
                  'ptext':0xd9313225f88406e5a55909c5aff5269a86a7a9531534f7da2e4c303d8a318a721c3c0c95956809532fcf0e2449a6b525b16aedf5aa0de657ba637b39,
                  'aad':0xfeedfacedeadbeeffeedfacedeadbeefabaddad2,
                  'iv':0x9313225df88406e555909c5aff5269aa6a7a9538534f7da1e4c303d2a318a728c3c0c95156809539fcf0e2429a6b525416aedbf5a0de6a57a637b39b,
                  'ctext':0x5a8def2f0c9e53f1f75d7853659e2a20eeb2b22aafde6419a058ab4f6f746bf40fc0c3b780f244452da3ebf1c5d82cdea2418997200ef82e44ae7e3f,
                  'tag':0xa44a8266ee1c8eb0c8b5d4cf5ae9f19a}
                  
    for i in range(1,19):    
        checkTestVector(**(eval("testcase" + str(i))))

Multiplication over the binary finite field GF(2^m)

The multGF2() function shown in the Python script below implements the element (polynomial) multiplication over a binary finite field GF(2^{m}). The second function, setGF2(), sets the three constants needed for its colleague to perform its multiplication task: "mask1" and "mask2" (used in “and” operations) and "polyred", a polynomial constant used in the polynomial reduction of the product. These constants are defined from the parameters passed to setGF2(): "degree" (the extension degree of the binary field GF(2^{m})) and "irPoly" (the irreducible polynomial used for the product reduction). The version of multGF2() in the post AES implementation in Python, dubbed “mult()“, was configured specifically for the AES binary finite field used in the MixColumns operation, GF(2^{8})/x^{8} + x^{4} + x^{3} + x + 1.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: Binary finite field multiplication in Python 3
#
# Date: 2012-02-16
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================
from functools import reduce

# constants used in the multGF2 function
mask1 = mask2 = polyred = None

def setGF2(degree, irPoly):
    """Define parameters of binary finite field GF(2^m)/g(x)
       - degree: extension degree of binary field
       - irPoly: coefficients of irreducible polynomial g(x)
    """
    global mask1, mask2, polyred
    mask1 = mask2 = 1 << degree
    mask2 -= 1
    if sum(irPoly) <= len(irPoly):
        polyred = reduce(lambda x, y: (x << 1) + y, irPoly[1:])    
    else:
        polyred = poly2Int(irPoly[1:]) 
        
def multGF2(p1, p2):
    """Multiply two polynomials in GF(2^m)/g(x)"""
    p = 0
    while p2:
        if p2 & 1:
            p ^= p1
        p1 <<= 1
        if p1 & mask1:
            p1 ^= polyred
        p2 >>= 1
    return p & mask2

#=============================================================================
#                        Auxiliary formatting functions
#=============================================================================
def int2Poly(bInt):
    """Convert a "big" integer into a "high-degree" polynomial"""
    exp = 0
    poly = []
    while bInt:
        if bInt & 1:
            poly.append(exp)
        exp += 1
        bInt >>= 1
    return poly[::-1]

def poly2Int(hdPoly):
    """Convert a "high-degree" polynomial into a "big" integer"""
    bigInt = 0
    for exp in hdPoly:
        bigInt += 1 << exp
    return bigInt

def i2P(sInt):
    """Convert a "small" integer into a "low-degree" polynomial"""
    return [(sInt >> i) & 1 for i in reversed(range(sInt.bit_length()))]

def p2I(ldPoly):
    """Convert a "low-degree" polynomial into a "small" integer"""
    return reduce(lambda x, y: (x << 1) + y, ldPoly)

def ldMultGF2(p1, p2):
    """Multiply two "low-degree" polynomials in GF(2^n)/g(x)"""
    return multGF2(p2I(p1), p2I(p2))

def hdMultGF2(p1, p2):
    """Multiply two "high-degree" polynomials in GF(2^n)/g(x)"""
    return multGF2(poly2Int(p1), poly2Int(p2))

if __name__ == "__main__":
  
    # Define binary field GF(2^3)/x^3 + x + 1
    setGF2(3, [1, 0, 1, 1])
    
    # Alternative way to define GF(2^3)/x^3 + x + 1
    setGF2(3, i2P(0b1011))    
    
    # Check if (x + 1)(x^2 + 1) == x^2
    assert ldMultGF2([1, 1], [1, 0, 1]) == p2I([1, 0, 0])

    # Check if (x^2 + x + 1)(x^2 + 1) == x^2 + x
    assert ldMultGF2([1, 1, 1], [1, 0, 1]) == p2I([1, 1, 0])
    
    # Define binary field GF(2^8)/x^8 + x^4 + x^3 + x + 1
    setGF2(8, [1, 0, 0, 0, 1, 1, 0, 1, 1])
    
    # Alternative way to define GF(2^8)/x^8 + x^4 + x^3 + x + 1
    setGF2(8, i2P(0b100011011))
    
    # Check if (x)(x^7 + x^2 + x + 1) == x^4 + x^2 + 1
    assert ldMultGF2([1, 0], [1, 0, 0, 0, 0, 1, 1, 1]) == p2I([1, 0, 1, 0, 1])
    
    # Check if (x + 1)(x^6 + x^5 + x^3 + x^2 + x) == x^7 + x^5 + x^4 + x
    assert ldMultGF2([1, 1], [1, 1, 0, 1, 1, 1, 0]) == p2I([1, 0, 1, 1, 0, 0, 1, 0])

    # Define binary field GF(2^571)/x^571 + x^10 + x^5 + x^2 + x
    setGF2(571, [571, 10, 5, 2, 1])

    # Calculate the product of two polynomials in GF(2^571)/x^571 + x^10 + x^5 + x^2 + x,
    # x^518 + x^447 + x^320 + x^209 + x^119 + x + 1 and x^287 + x^145 + x^82 + + x^44
    print(int2Poly(hdMultGF2([518, 447, 320, 209, 119, 1, 0], [287, 145, 82, 44])))

The output of the above script is:

[562, 529, 496, 491, 465, 406, 402, 364, 354, 291, 288, 287, 264, 253, 244,
239, 235, 201, 173, 168, 165, 164, 163, 146, 145, 102, 97, 94, 93, 83, 82,
46, 45, 44, 41, 39, 38, 37, 34, 30, 26, 23, 22]

This means that the product of the polynomials x^{518} + x^{447} + x^{320} + x^{209} + x^{119} + x + 1 and x^{287} + x^{145} + x^{82} + x^{44} is:

x^562 + x^529 + x^496 + x^491 + x^465 + x^406 + x^402 + x^364 + x^354 +
x^291 + x^288 + x^287 + x^264 + x^253 + x^244 + x^239 + x^236 + x^235 +
x^201 + x^173 + x^168 + x^165 + x^164 + x^163 + x^146 + x^145 + x^102 +
x^97 + x^94 + x^93 + x^83 + x^82 + x^46 + x^45 + x^44 + x^41 + x^39 +
x^38 + x^37 + x^34 + x^30 + x^26 + x^23 + x^22

An elegant algorithm for calculating Hamming distance

Generally speaking, the Hamming distance between two strings of the same length is the number of positions in which the corresponding symbols are different. This metric, proposed by Richard W. Hamming in his seminal paper “Error Detecting and Error Correcting Codes”, Bell System Technical Journal 26(2):147-160 (1950), also expresses the least number of (intended or unintended) substitutions needed to change one string to another.

Besides being used in computer- and communications-related fields such as information theory, coding theory, and cryptography, the Hamming distance concept has also found its way into genomics for the comparison of genomic sequences. Its importance can also be judged by the fact that modern microprocessors have a specific instruction (POPCNT, “population count”) for counting the number of bits set to 1.

Assuming two bit strings of equal length, x and y, the “obvious” algorithm to calculate the Hamming distance between them is to count the number of “1” bits in the result of the expression “x xor y”, as shown in the following Python code:

def hamming1(x,y):
    """Calculate the Hamming distance between two bit strings"""
    assert len(x) == len(y)
    count,z = 0,int(x,2)^int(y,2)
    while z:
        if z&1:
            count += 1
        z >>= 1
    return count

While searching the web the other day, I have come across a very elegant algorithm for the same purpose, which can be coded in Python as follows:

def hamming2(x,y):
    """Calculate the Hamming distance between two bit strings"""
    assert len(x) == len(y)
    count,z = 0,int(x,2)^int(y,2)
    while z:
        count += 1
        z &= z-1 # magic!
    return count

Apart from being 50% faster than the first one, its simplicity is awesome! After being counted, the lowest-order nonzero bit is cleared like magic! When I started a new search later, this time looking for the algorithm’s authorship, I finally learned that it was proposed by Peter Wegner in 1960 (“A technique for counting ones in a binary computer“, Communications of the ACM 3 (5): 322).

Where did I find this information? The answer was, talking again about “obvious” things, in the very same Wikipedia’s entry for “Hamming distance“! This little gem has been there for at least a year (according to the entry’s history), just waiting to be discovered! Quoting Antoine de Saint Exupéry, “You know you have achieved perfection in design, not when you have nothing more to add, but when you have nothing more to take away.”

Simplified AES implementation in Python

Simplified AES, created by Edward Schaefer and two of his students at Santa Clara University in 2003, is described in the paper “A Simplified AES Algorithm and Its Linear and Differential Cryptoanalyses”, Cryptologia, Vol. XXVII (2), pages 148-177. Like Simplified DES, its purpose is educational, since its key and block size are very small (16-bit and 8-bit, respectively). Nonetheless, due to its simplicity, it is possible for students to encrypt or decrypt a block doing all operations by hand, making it easier for them to understand the structure and operation of its elder brother, AES.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: Simplified AES implementation in Python 3
#
# Date: 2012-02-11
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================
import sys

# S-Box
sBox  = [0x9, 0x4, 0xa, 0xb, 0xd, 0x1, 0x8, 0x5,
         0x6, 0x2, 0x0, 0x3, 0xc, 0xe, 0xf, 0x7]

# Inverse S-Box
sBoxI = [0xa, 0x5, 0x9, 0xb, 0x1, 0x7, 0x8, 0xf,
         0x6, 0x0, 0x2, 0x3, 0xc, 0x4, 0xd, 0xe]

# Round keys: K0 = w0 + w1; K1 = w2 + w3; K2 = w4 + w5
w = [None] * 6

def mult(p1, p2):
    """Multiply two polynomials in GF(2^4)/x^4 + x + 1"""
    p = 0
    while p2:
        if p2 & 0b1:
            p ^= p1
        p1 <<= 1
        if p1 & 0b10000:
            p1 ^= 0b11
        p2 >>= 1
    return p & 0b1111

def intToVec(n):
    """Convert a 2-byte integer into a 4-element vector"""
    return [n >> 12, (n >> 4) & 0xf, (n >> 8) & 0xf,  n & 0xf]            

def vecToInt(m):
    """Convert a 4-element vector into 2-byte integer"""
    return (m[0] << 12) + (m[2] << 8) + (m[1] << 4) + m[3]

def addKey(s1, s2):
    """Add two keys in GF(2^4)"""   
    return [i ^ j for i, j in zip(s1, s2)]
    
def sub4NibList(sbox, s):
    """Nibble substitution function"""
    return [sbox[e] for e in s]
    
def shiftRow(s):
    """ShiftRow function"""
    return [s[0], s[1], s[3], s[2]]

def keyExp(key):
    """Generate the three round keys"""
    def sub2Nib(b):
        """Swap each nibble and substitute it using sBox"""
        return sBox[b >> 4] + (sBox[b & 0x0f] << 4)

    Rcon1, Rcon2 = 0b10000000, 0b00110000
    w[0] = (key & 0xff00) >> 8
    w[1] = key & 0x00ff
    w[2] = w[0] ^ Rcon1 ^ sub2Nib(w[1])
    w[3] = w[2] ^ w[1]
    w[4] = w[2] ^ Rcon2 ^ sub2Nib(w[3])
    w[5] = w[4] ^ w[3]

def encrypt(ptext):
    """Encrypt plaintext block"""
    def mixCol(s):
        return [s[0] ^ mult(4, s[2]), s[1] ^ mult(4, s[3]),
                s[2] ^ mult(4, s[0]), s[3] ^ mult(4, s[1])]    
    
    state = intToVec(((w[0] << 8) + w[1]) ^ ptext)
    state = mixCol(shiftRow(sub4NibList(sBox, state)))
    state = addKey(intToVec((w[2] << 8) + w[3]), state)
    state = shiftRow(sub4NibList(sBox, state))
    return vecToInt(addKey(intToVec((w[4] << 8) + w[5]), state))
    
def decrypt(ctext):
    """Decrypt ciphertext block"""
    def iMixCol(s):
        return [mult(9, s[0]) ^ mult(2, s[2]), mult(9, s[1]) ^ mult(2, s[3]),
                mult(9, s[2]) ^ mult(2, s[0]), mult(9, s[3]) ^ mult(2, s[1])]
    
    state = intToVec(((w[4] << 8) + w[5]) ^ ctext)
    state = sub4NibList(sBoxI, shiftRow(state))
    state = iMixCol(addKey(intToVec((w[2] << 8) + w[3]), state))
    state = sub4NibList(sBoxI, shiftRow(state))
    return vecToInt(addKey(intToVec((w[0] << 8) + w[1]), state))

if __name__ == '__main__':
    # Test vectors from "Simplified AES" (Steven Gordon)
    # (http://hw.siit.net/files/001283.pdf)
    
    plaintext = 0b1101011100101000
    key = 0b0100101011110101
    ciphertext = 0b0010010011101100
    keyExp(key)
    try:
        assert encrypt(plaintext) == ciphertext
    except AssertionError:
        print("Encryption error")
        print(encrypt(plaintext), ciphertext)
        sys.exit(1)
    try:
        assert decrypt(ciphertext) == plaintext
    except AssertionError:
        print("Decryption error")
        print(decrypt(ciphertext), plaintext)
        sys.exit(1)
    print("Test ok!")
    sys.exit()

Generating Permutations From Factoradic Numbers

This is the title of an interesting tutorial by Prof. Joel Guilherme S. Filho. It describes how to generate permutations of a set of n elements from a given decimal index based on its representation as a “factoradic number”, i.e. a number represented in the mixed base [(n-1)!, (n-2)!, … , 2!, 1!, 0!].

Thanks to the list manipulation facilities and the support for arbitrary long integers provided by the Python programming language, it’s easy to write a short script (shown below) that generates permutations even for large indexes. Please read his tutorial for additional details.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: Generating Permutations From Factoradic Numbers
#
# Date: 2012-02-06
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#=============================================================
from sys import argv
from math import factorial

N = int(argv[1])

n = 0
while N > factorial(n):
    n += 1

S = [i for i in range(n)]
Factoradic, Permutation = [None] * n, [None] * n

v, s = N, S[:]
for i in range(n):
    r, v = divmod(v, factorial(n - i - 1))
    Factoradic[i], Permutation[i] = r, s.pop(r)

print("For N = {:d}\n".format(N))
print("n = {:d}\n".format(n))
print("Permutation[0] = {:s}\n".format(S))
print("Factoradic[N] = {:s}\n".format(Factoradic))
print("Permutation[N] = {:s}\n".format(Permutation))

Executing this script with 51 as the command-line argument, we get the output below:

For N = 51

n = 5

Perm[0] = [0, 1, 2, 3, 4]

Factoradic[N] = [2, 0, 1, 1, 0]

Permutation[N] = [2, 0, 3, 4, 1]

Trying now a big number like the Mersenne number \small 2^{607} - 1 (which happens to be the Mersenne prime # 14), we get the following results:

For N = 53113799281676709868958820655246862732959311772703192319944413 ... 29486246501015346579337652707239409519978766587351943831270835393219031728127

n = 113

Permutation[0] = [0, 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, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112]

Factoradic[N] = [2, 77, 30, 88, 71, 102, 97, 70, 4, 12, 75, 71, 27, 9, 91, 43, 75, 58, 30, 73, 10, 11, 68, 33, 77, 79, 50, 80, 41, 77, 48, 5, 32, 78, 25, 74, 0, 9, 49, 3, 43, 14, 6, 20, 2, 6, 61, 14, 29, 7, 37, 41, 5, 15, 30, 54, 26, 0, 41, 19, 29, 50, 6, 6, 2, 25, 6, 8, 44, 10, 26, 13, 16, 0, 14, 4, 5, 24, 24, 30, 22, 6, 6, 29, 10, 24, 12, 23, 3, 20, 15, 20, 16, 11, 12, 7, 2, 15, 10, 5, 5, 8, 6, 7, 2, 0, 0, 1, 0, 1, 0, 1, 0]

Permutation[N] = [2, 78, 31, 91, 73, 107, 102, 72, 5, 14, 82, 77, 30, 11, 104, 49, 87, 65, 36, 88, 13, 16, 84, 42, 98, 101, 61, 106, 52, 103, 60, 7, 43, 111, 34, 108, 0, 17, 69, 6, 63, 24, 12, 35, 4, 18, 99, 28, 53, 20, 67, 75, 15, 37, 58, 105, 54, 1, 85, 45, 64, 110, 22, 23, 9, 62, 26, 32, 112, 39, 74, 46, 51, 3, 50, 25, 29, 83, 86, 96, 80, 38, 40, 109, 55, 94, 59, 95, 21, 92, 76, 97, 81, 66, 70, 47, 19, 100, 71, 44, 48, 79, 57, 89, 27, 8, 10, 41, 33, 68, 56, 93, 90]

CLEFIA implementation in Python (improved version)

As commented in my last post, CLEFIA implementation in Python, the performance figures for key setup and block encryption of my first coding attempt were poor when compared to those of AES.

Nevertheless, in order to establish a fair comparison between the two implementations, I included in CLEFIA’s code the same optimization technique used in AES, memoization, for speeding up the multiplication in the binary finite field GF(2^8). For this purpose, the generic multiplication function used in the first version was replaced by five memoized functions, x2(), x4(), x6(), x8() and x10(), thanks to the fact that there are only five different constants (besides the multiplicative identity of the field, 1) in the matrices m0 and m1.

The payback of this small effort was remarkable: the new code is 7 times(!) faster than the old, as the time required for the 128-bit key setup plus block encryption (in ECB mode) went down from 3.8 ms to 0.5 ms. A similar reduction was observed for 192- and 256-bit keys: from 5.2ms to 0.7 ms and 5.7 ms to 0.77 ms. All numbers were observed in a Core i7-2600K @ 3.5 GHz desktop running 64-bit Windows 7. The Python interpreter used was ActivePython-3.2.2.3-win64-x64.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: CLEFIA implementation in Python 3
#
# Date: 2012-02-03
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================
import sys
from time import time

# Key sizes supported
ksTable = {"SIZE_128": 16,
           "SIZE_192": 24,
           "SIZE_256": 32}

# Number of rounds related to key size
nrTable = {"SIZE_128": 18,
           "SIZE_192": 22,
           "SIZE_256": 26}

# Number of round keys related to key size
nrkTable = {"SIZE_128": 36,
            "SIZE_192": 44,
            "SIZE_256": 52}

# Number of rounds
nr = None

# Number of round keys effectively used
nrk = None

# Number of whitening keys
nwk = 4

# Round keys vector
rk = [None] * 2 * nrTable[max(nrTable)]

# Whitening keys
wk = [None] * 4

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

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

m0 = [0x01, 0x02, 0x04, 0x06, 0x02, 0x01, 0x06, 0x04,
      0x04, 0x06, 0x01, 0x02, 0x06, 0x04, 0x02, 0x01]

m1 = [0x01, 0x08, 0x02, 0x0a, 0x08, 0x01, 0x0a, 0x02,
      0x02, 0x0a, 0x01, 0x08, 0x0a, 0x02, 0x08, 0x01]

con128 = [0xf56b7aeb, 0x994a8a42, 0x96a4bd75, 0xfa854521,
          0x735b768a, 0x1f7abac4, 0xd5bc3b45, 0xb99d5d62,
          0x52d73592, 0x3ef636e5, 0xc57a1ac9, 0xa95b9b72,
          0x5ab42554, 0x369555ed, 0x1553ba9a, 0x7972b2a2,
          0xe6b85d4d, 0x8a995951, 0x4b550696, 0x2774b4fc,
          0xc9bb034b, 0xa59a5a7e, 0x88cc81a5, 0xe4ed2d3f,
          0x7c6f68e2, 0x104e8ecb, 0xd2263471, 0xbe07c765,
          0x511a3208, 0x3d3bfbe6, 0x1084b134, 0x7ca565a7,
          0x304bf0aa, 0x5c6aaa87, 0xf4347855, 0x9815d543,
          0x4213141a, 0x2e32f2f5, 0xcd180a0d, 0xa139f97a,
          0x5e852d36, 0x32a464e9, 0xc353169b, 0xaf72b274,
          0x8db88b4d, 0xe199593a, 0x7ed56d96, 0x12f434c9,
          0xd37b36cb, 0xbf5a9a64, 0x85ac9b65, 0xe98d4d32,
          0x7adf6582, 0x16fe3ecd, 0xd17e32c1, 0xbd5f9f66,
          0x50b63150, 0x3c9757e7, 0x1052b098, 0x7c73b3a7]

con192 = [0xc6d61d91, 0xaaf73771, 0x5b6226f8, 0x374383ec,
          0x15b8bb4c, 0x799959a2, 0x32d5f596, 0x5ef43485,
          0xf57b7acb, 0x995a9a42, 0x96acbd65, 0xfa8d4d21,
          0x735f7682, 0x1f7ebec4, 0xd5be3b41, 0xb99f5f62,
          0x52d63590, 0x3ef737e5, 0x1162b2f8, 0x7d4383a6,
          0x30b8f14c, 0x5c995987, 0x2055d096, 0x4c74b497,
          0xfc3b684b, 0x901ada4b, 0x920cb425, 0xfe2ded25,
          0x710f7222, 0x1d2eeec6, 0xd4963911, 0xb8b77763,
          0x524234b8, 0x3e63a3e5, 0x1128b26c, 0x7d09c9a6,
          0x309df106, 0x5cbc7c87, 0xf45f7883, 0x987ebe43,
          0x963ebc41, 0xfa1fdf21, 0x73167610, 0x1f37f7c4,
          0x01829338, 0x6da363b6, 0x38c8e1ac, 0x54e9298f,
          0x246dd8e6, 0x484c8c93, 0xfe276c73, 0x9206c649,
          0x9302b639, 0xff23e324, 0x7188732c, 0x1da969c6,
          0x00cd91a6, 0x6cec2cb7, 0xec7748d3, 0x8056965b,
          0x9a2aa469, 0xf60bcb2d, 0x751c7a04, 0x193dfdc2,
          0x02879532, 0x6ea666b5, 0xed524a99, 0x8173b35a,
          0x4ea00d7c, 0x228141f9, 0x1f59ae8e, 0x7378b8a8,
          0xe3bd5747, 0x8f9c5c54, 0x9dcfaba3, 0xf1ee2e2a,
          0xa2f6d5d1, 0xced71715, 0x697242d8, 0x055393de,
          0x0cb0895c, 0x609151bb, 0x3e51ec9e, 0x5270b089]

con256 = [0x0221947e, 0x6e00c0b5, 0xed014a3f, 0x8120e05a,
          0x9a91a51f, 0xf6b0702d, 0xa159d28f, 0xcd78b816,
          0xbcbde947, 0xd09c5c0b, 0xb24ff4a3, 0xde6eae05,
          0xb536fa51, 0xd917d702, 0x62925518, 0x0eb373d5,
          0x094082bc, 0x6561a1be, 0x3ca9e96e, 0x5088488b,
          0xf24574b7, 0x9e64a445, 0x9533ba5b, 0xf912d222,
          0xa688dd2d, 0xcaa96911, 0x6b4d46a6, 0x076cacdc,
          0xd9b72353, 0xb596566e, 0x80ca91a9, 0xeceb2b37,
          0x786c60e4, 0x144d8dcf, 0x043f9842, 0x681edeb3,
          0xee0e4c21, 0x822fef59, 0x4f0e0e20, 0x232feff8,
          0x1f8eaf20, 0x73af6fa8, 0x37ceffa0, 0x5bef2f80,
          0x23eed7e0, 0x4fcf0f94, 0x29fec3c0, 0x45df1f9e,
          0x2cf6c9d0, 0x40d7179b, 0x2e72ccd8, 0x42539399,
          0x2f30ce5c, 0x4311d198, 0x2f91cf1e, 0x43b07098,
          0xfbd9678f, 0x97f8384c, 0x91fdb3c7, 0xfddc1c26,
          0xa4efd9e3, 0xc8ce0e13, 0xbe66ecf1, 0xd2478709,
          0x673a5e48, 0x0b1bdbd0, 0x0b948714, 0x67b575bc,
          0x3dc3ebba, 0x51e2228a, 0xf2f075dd, 0x9ed11145,
          0x417112de, 0x2d5090f6, 0xcca9096f, 0xa088487b,
          0x8a4584b7, 0xe664a43d, 0xa933c25b, 0xc512d21e,
          0xb888e12d, 0xd4a9690f, 0x644d58a6, 0x086cacd3,
          0xde372c53, 0xb216d669, 0x830a9629, 0xef2beb34,
          0x798c6324, 0x15ad6dce, 0x04cf99a2, 0x68ee2eb3]

def _8To32(x32):
    """Convert a 4-byte list to a 32-bit integer"""
    return (((((x32[0] << 8) + x32[1]) << 8) + x32[2]) << 8) + x32[3]    

def _32To8(x32):
    """Convert a 32-bit integer to a 4-byte list"""
    return [(x32 >> 8 * i) & 0xff for i in reversed(range(4))]

def _32To128(x32):
    """Convert a 32-bit 4-element list to a 128-bit integer"""
    return (((((x32[0] << 32) + x32[1]) << 32) + x32[2]) << 32) + x32[3]

def _128To32(x128):
    """Convert a 128-bit integer to a 32-bit 4-element list"""
    return [(x128 >> 32 * i) & 0xffffffff for i in reversed(range(4))]

def _192To32(x192):
    """Convert a 192-bit integer to a 32-bit 6-element list"""
    return [(x192 >> 32 * i) & 0xffffffff for i in reversed(range(6))]

def _256To32(x256):
    """Convert a 256-bit integer to a 32-bit 8-element list"""
    return [(x256 >> 32 * i) & 0xffffffff for i in reversed(range(8))]

def sigma(x128):
    """The double-swap function sigma (used in key scheduling)"""
    return [(x128[0] << 7) & 0xffffff80  | (x128[1] >> 25),
            (x128[1] << 7) & 0xffffff80  | (x128[3] & 0x7f),
            (x128[0] & 0xfe000000)       | (x128[2] >> 7),
            (x128[2] << 25) & 0xfe000000 | (x128[3] >> 7)]

def memoize(f):
    """Memoization function"""
    memo = {}
    def helper(x):
        if x not in memo:
            memo[x] = f(x)
        return memo[x]
    return helper

def mult(p1, p2):
    """Multiply two polynomials in GF(2^8)
       (the irreducible polynomial used in this
       field is x^8 + x^4 + x^3 + x^2 + 1)"""
    p = 0
    while p2:
        if p2 & 0b1:
            p ^= p1
        p1 <<= 1
        if p1 & 0x100:
            p1 ^= 0b11101
        p2 >>= 1
    return p & 0xff

# Auxiliary one-parameter functions defined for memoization
# (to speed up multiplication in GF(2^8))

@memoize    
def x2(y):
    """Multiply by 2 in GF(2^8)"""
    return mult(2, y)

@memoize
def x4(y):
    """Multiply by 4 in GF(2^8)"""    
    return mult(4, y)

@memoize
def x6(y):
    """Multiply by 6 in GF(2^8)"""        
    return mult(6, y)

@memoize
def x8(y):
    """Multiply by 8 in GF(2^8)"""    
    return mult(8, y)

@memoize
def x10(y):
    """Multiply by 10 in GF(2^8)"""    
    return mult(10, y)

def multm0(t32):
    """Multiply the matrix m0 by a 4-element transposed vector in GF(2^8)"""
    return [   t32[0]  ^ x2(t32[1]) ^ x4(t32[2]) ^ x6(t32[3]),
            x2(t32[0]) ^    t32[1]  ^ x6(t32[2]) ^ x4(t32[3]),
            x4(t32[0]) ^ x6(t32[1]) ^    t32[2]  ^ x2(t32[3]),
            x6(t32[0]) ^ x4(t32[1]) ^ x2(t32[2]) ^    t32[3]]    
    
def multm1(t32):
    """Multiply the matrix m1 by a 4-element transposed vector in GF(2^8)"""
    return [    t32[0]  ^  x8(t32[1]) ^  x2(t32[2]) ^ x10(t32[3]),
             x8(t32[0]) ^     t32[1]  ^ x10(t32[2]) ^  x2(t32[3]),
             x2(t32[0]) ^ x10(t32[1]) ^     t32[2]  ^  x8(t32[3]),
            x10(t32[0]) ^ x2(t32[1])  ^  x8(t32[2]) ^     t32[3]]      

def f0(rk, x32):
    """F0 function"""
    t8 = _32To8(rk ^ x32)
    t8 = [s0[t8[0]], s1[t8[1]], s0[t8[2]], s1[t8[3]]]
    return _8To32(multm0(t8))

def f1(rk, x32):
    """F1 function"""
    t8 = _32To8(rk ^ x32)
    t8 = s1[t8[0]], s0[t8[1]], s1[t8[2]], s0[t8[3]]
    return _8To32(multm1(t8))

def gfn4(x32, n):
    """4-branch Generalized Feistel Network function"""
    t32 = x32[:]
    for i in range(0, n << 1, 2):
        t32[1] ^= f0(rk[i], t32[0])
        t32[3] ^= f1(rk[i + 1], t32[2])
        t32 = t32[1:] + t32[:1]
    return t32[3:] + t32[:3]

def gfn4i(x32, n):
    """4-branch Generalized Feistel Network inverse function"""
    t32 = x32[:]
    for i in reversed(range(0, n << 1, 2)):
        t32[1] ^= f0(rk[i], t32[0])
        t32[3] ^= f1(rk[i + 1], t32[2])
        t32 = t32[3:] + t32[:3]
    return t32[1:] + t32[:1]

def gfn8(x32, n):
    """8-branch Generalized Feistel Network function"""
    t32 = x32[:]
    for i in range(0, n << 2, 4):
        t32[1] ^= f0(rk[i], t32[0])
        t32[3] ^= f1(rk[i + 1], t32[2])
        t32[5] ^= f0(rk[i + 2], t32[4])
        t32[7] ^= f1(rk[i + 3], t32[6])     
        t32 = t32[1:] + t32[:1]
    return t32[7:] + t32[:7]

def setKey128(k128):
    """Generate round/whitening keys from a 128-bit key"""
    k32 = _128To32(k128)
    for i in range(len(con128) - nrk):
        rk[i] = con128[i]
    l = gfn4(k32, 12)
    for i in range(nwk):
        wk[i] = k32[i]
    for i in range(0, nrk, 4):
        t32 = [r ^ s for r, s in zip(l, con128[i + 24:i + 28])]
        l = sigma(l)
        if i & 0b100:
            rk[i:i + 4] = [r ^ s for r, s in zip(t32, k32)]
        else:
            rk[i:i + 4] = t32

def setKey192(k192):
    """Generate round/whitening keys from a 192-bit key"""
    k32 = _192To32(k192)
    kl = k32[:4]
    kr = k32[4:6] + [k32[0] ^ 0xffffffff] + [k32[1] ^ 0xffffffff]
    for i in range(len(con192) - nrk):
        rk[i] = con192[i]        
    l = gfn8(kl + kr, 10)
    ll, lr = l[:4], l[4:]
    kk = [r ^ s for r, s in zip(kl, kr)]
    for i in range(nwk):
        wk[i] = kk[i]
    for i in range(0, nrk, 4):
        if i & 0b1100 < 8:
            t32 = [r ^ s for r, s in zip(ll, con192[i + 40:i + 44])]
            ll = sigma(ll)
            if i & 0b100:
                t32 = [r ^ s for r, s in zip(t32, kr)]
        else:
            t32 = [r ^ s for r, s in zip(lr, con192[i + 40:i + 44])]
            lr = sigma(lr)
            if i & 0b100:
                t32 = [r ^ s for r, s in zip(t32, kl)]            
        rk[i:i + 4] = t32     

def setKey256(k256):
    """Generate round/whitening keys from a 256-bit key"""
    k32 = _256To32(k256)
    kl, kr = k32[:4], k32[4:]
    for i in range(len(con256) - nrk):
        rk[i] = con256[i]    
    l = gfn8(kl + kr, 10)
    ll, lr = l[:4], l[4:]
    kk = [r ^ s for r, s in zip(kl, kr)]    
    for i in range(nwk):
        wk[i] = kk[i]
    for i in range(0, nrk, 4):
        if i & 0b1100 < 8:
            t32 = [r ^ s for r, s in zip(ll, con256[i + 40:i + 44])]
            ll = sigma(ll)
            if i & 0b100:
                t32 = [r ^ s for r, s in zip(t32, kr)]
        else:
            t32 = [r ^ s for r, s in zip(lr, con256[i + 40:i + 44])]
            lr = sigma(lr)
            if i & 0b100:
                t32 = [r ^ s for r, s in zip(t32, kl)]                
        rk[i:i + 4] = t32

def setKey(key, keySize):
    """Generate round/whitening keys from the given key"""
    global nr, nrk
    try:
        assert keySize in ksTable
    except AssertionError:
        print("Key size identifier not valid")
        sys.exit("ValueError")
    try:
        assert isinstance(key, int)
    except AssertionError:
        print("Invalid key")
        sys.exit("ValueError")
    try:
        assert key.bit_length() // 8 <= ksTable[keySize]
    except AssertionError:
        print("Key size mismatch")
        sys.exit("ValueError")
    nr = nrTable[keySize]
    nrk = nrkTable[keySize]
    if keySize == "SIZE_128":
        setKey128(key)
    elif keySize == "SIZE_192":
        setKey192(key)
    elif keySize == "SIZE_256":
        setKey256(key)
    else:
        sys.exit("Invalid key size identifier")

def encrypt(ptext):
    """Encrypt a block"""
    t32 = _128To32(ptext)
    t32[1] ^= wk[0]
    t32[3] ^= wk[1]
    t32 = gfn4(t32, nr)
    t32[1] ^= wk[2]
    t32[3] ^= wk[3]
    return _32To128(t32)

def decrypt(ctext):
    """Decrypt a block"""
    t32 = _128To32(ctext)
    t32[1] ^= wk[2]
    t32[3] ^= wk[3]
    t32 = gfn4i(t32, nr)
    t32[1] ^= wk[0]
    t32[3] ^= wk[1]
    return _32To128(t32)

if __name__ == "__main__":

    def checkTestVector(key, keySize, plaintext, ciphertext, nIter = 1000):
        testSuccess = True
        setKey(key, keySize)
        ks = ksTable[keySize] * 8
        ctext = encrypt(plaintext)
        ptext = decrypt(ctext)
        try:
            assert ctext == ciphertext
        except AssertionError:
            print("Error in encryption")
            print("Resulting ciphertext: {:02x}".format(ctext))
            print("Expected ciphertext: {:02x}".format(ciphertext))
            testSuccess = False
        try:
            assert ptext == plaintext
        except AssertionError:
            print("Error in decryption:")
            print("Recovered plaintext: {:02x}".format(ptext))
            print("Expected plaintext: {:02x}".format(plaintext))
            testSuccess = False         
        if not testSuccess:
            return False
        t1 = time()        
        for i in range(nIter):
            setKey(key, keySize)
            ctext = encrypt(plaintext)
        t2 = time()
        avg_elapsed_time = (t2 - t1) * 1000 / nIter        
        print("{:3d}-bit key test ok!".format(ksTable[keySize] * 8))
        print("Average elapsed time for 16-byte block ", end="")
        print("ECB-{0:3d} encryption: {1:0.3f}ms".format(ks, avg_elapsed_time))
        t3 = time()   
        for i in range(nIter):
            setKey(key, keySize)
            ptext = decrypt(ctext)
        t4 = time()
        avg_elapsed_time = (t4 - t3) * 1000 / nIter        
        print("{:3d}-bit key test ok!".format(ksTable[keySize] * 8))
        print("Average elapsed time for 16-byte block ", end="")
        print("ECB-{0:3d} decryption: {1:0.3f}ms".format(ks, avg_elapsed_time))        
        return True

    # The test vectors below are described in document "The 128-bit Blockcipher
    # CLEFIA Algorithm Specification" rev.1, June 1, 2007, Sony Corporation.

    ptext = 0x000102030405060708090a0b0c0d0e0f

    # Test vector for 128-bit key
    key1 = 0xffeeddccbbaa99887766554433221100
    ctext1 = 0xde2bf2fd9b74aacdf1298555459494fd

    # Test vector for 192-bit key
    key2 = 0xffeeddccbbaa99887766554433221100f0e0d0c0b0a09080
    ctext2 = 0xe2482f649f028dc480dda184fde181ad

    # Test vector for 256-bit key
    key3 = 0xffeeddccbbaa99887766554433221100f0e0d0c0b0a090807060504030201000
    ctext3 = 0xa1397814289de80c10da46d1fa48b38a

    try:
        assert checkTestVector(key1, "SIZE_128", ptext, ctext1) and \
               checkTestVector(key2, "SIZE_192", ptext, ctext2) and \
               checkTestVector(key3, "SIZE_256", ptext, ctext3)
    except AssertionError:
        print("At least one test failed")
        sys.exit(1)
    print("All tests passed!")
    sys.exit()

CLEFIA implementation in Python

As Wikipedia puts it, “CLEFIA is a new block cipher algorithm developed by Sony. Its name is derived from the French word clef, meaning “key”. The block size is 128 bits and the key size can be 128 bit, 192 bit or 256 bit [exactly like NIST’s Advanced Encryption Standard (AES)]. It is intended to be used in DRM systems.”

CLEFIA has just been standardized in ISO/IEC 29192 (Information technology – Security techniques – Lightweight cryptography – Part 2: Block ciphers). The rationale for “lightweight cryptography” is described in the short paper “Lightweight Cryptography for the Internet of Things” (Katagi & Moriai). CLEFIA’s architecture is presented in RFC 6114.

According to Sony, “when implemented in hardware [CLEFIA] achieves the world’s highest hardware gate efficiency. When implemented in software it can realize high speed performance on a wide variety of processors.” My AES Python implementation is, however, an order of magnitude faster than CLEFIA’s (see below) with respect to key setup & single block encryption. I think, though, that there is room for some improvement, as this version of CLEFIA is the very first.

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: CLEFIA implementation in Python 3
#
# Date: 2012-02-02
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================
import sys
from time import time

# Key size supported
ksTable = {"SIZE_128": 16,
           "SIZE_192": 24,
           "SIZE_256": 32}

# Number of rounds related to key size
nrTable = {"SIZE_128": 18,
           "SIZE_192": 22,
           "SIZE_256": 26}

# Number of round keys related to key size
nrkTable = {"SIZE_128": 36,
            "SIZE_192": 44,
            "SIZE_256": 52}

# Number of rounds
nr = None

# Number of round keys effectively used
nrk = None

# Number of whitening keys
nwk = 4

# Round keys vector
rk = [None] * 2 * nrTable[max(nrTable)]

# Whitening keys
wk = [None] * 4

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

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

m0 = [0x01, 0x02, 0x04, 0x06, 0x02, 0x01, 0x06, 0x04,
      0x04, 0x06, 0x01, 0x02, 0x06, 0x04, 0x02, 0x01]

m1 = [0x01, 0x08, 0x02, 0x0a, 0x08, 0x01, 0x0a, 0x02,
      0x02, 0x0a, 0x01, 0x08, 0x0a, 0x02, 0x08, 0x01]

con128 = [0xf56b7aeb, 0x994a8a42, 0x96a4bd75, 0xfa854521,
          0x735b768a, 0x1f7abac4, 0xd5bc3b45, 0xb99d5d62,
          0x52d73592, 0x3ef636e5, 0xc57a1ac9, 0xa95b9b72,
          0x5ab42554, 0x369555ed, 0x1553ba9a, 0x7972b2a2,
          0xe6b85d4d, 0x8a995951, 0x4b550696, 0x2774b4fc,
          0xc9bb034b, 0xa59a5a7e, 0x88cc81a5, 0xe4ed2d3f,
          0x7c6f68e2, 0x104e8ecb, 0xd2263471, 0xbe07c765,
          0x511a3208, 0x3d3bfbe6, 0x1084b134, 0x7ca565a7,
          0x304bf0aa, 0x5c6aaa87, 0xf4347855, 0x9815d543,
          0x4213141a, 0x2e32f2f5, 0xcd180a0d, 0xa139f97a,
          0x5e852d36, 0x32a464e9, 0xc353169b, 0xaf72b274,
          0x8db88b4d, 0xe199593a, 0x7ed56d96, 0x12f434c9,
          0xd37b36cb, 0xbf5a9a64, 0x85ac9b65, 0xe98d4d32,
          0x7adf6582, 0x16fe3ecd, 0xd17e32c1, 0xbd5f9f66,
          0x50b63150, 0x3c9757e7, 0x1052b098, 0x7c73b3a7]

con192 = [0xc6d61d91, 0xaaf73771, 0x5b6226f8, 0x374383ec,
          0x15b8bb4c, 0x799959a2, 0x32d5f596, 0x5ef43485,
          0xf57b7acb, 0x995a9a42, 0x96acbd65, 0xfa8d4d21,
          0x735f7682, 0x1f7ebec4, 0xd5be3b41, 0xb99f5f62,
          0x52d63590, 0x3ef737e5, 0x1162b2f8, 0x7d4383a6,
          0x30b8f14c, 0x5c995987, 0x2055d096, 0x4c74b497,
          0xfc3b684b, 0x901ada4b, 0x920cb425, 0xfe2ded25,
          0x710f7222, 0x1d2eeec6, 0xd4963911, 0xb8b77763,
          0x524234b8, 0x3e63a3e5, 0x1128b26c, 0x7d09c9a6,
          0x309df106, 0x5cbc7c87, 0xf45f7883, 0x987ebe43,
          0x963ebc41, 0xfa1fdf21, 0x73167610, 0x1f37f7c4,
          0x01829338, 0x6da363b6, 0x38c8e1ac, 0x54e9298f,
          0x246dd8e6, 0x484c8c93, 0xfe276c73, 0x9206c649,
          0x9302b639, 0xff23e324, 0x7188732c, 0x1da969c6,
          0x00cd91a6, 0x6cec2cb7, 0xec7748d3, 0x8056965b,
          0x9a2aa469, 0xf60bcb2d, 0x751c7a04, 0x193dfdc2,
          0x02879532, 0x6ea666b5, 0xed524a99, 0x8173b35a,
          0x4ea00d7c, 0x228141f9, 0x1f59ae8e, 0x7378b8a8,
          0xe3bd5747, 0x8f9c5c54, 0x9dcfaba3, 0xf1ee2e2a,
          0xa2f6d5d1, 0xced71715, 0x697242d8, 0x055393de,
          0x0cb0895c, 0x609151bb, 0x3e51ec9e, 0x5270b089]

con256 = [0x0221947e, 0x6e00c0b5, 0xed014a3f, 0x8120e05a,
          0x9a91a51f, 0xf6b0702d, 0xa159d28f, 0xcd78b816,
          0xbcbde947, 0xd09c5c0b, 0xb24ff4a3, 0xde6eae05,
          0xb536fa51, 0xd917d702, 0x62925518, 0x0eb373d5,
          0x094082bc, 0x6561a1be, 0x3ca9e96e, 0x5088488b,
          0xf24574b7, 0x9e64a445, 0x9533ba5b, 0xf912d222,
          0xa688dd2d, 0xcaa96911, 0x6b4d46a6, 0x076cacdc,
          0xd9b72353, 0xb596566e, 0x80ca91a9, 0xeceb2b37,
          0x786c60e4, 0x144d8dcf, 0x043f9842, 0x681edeb3,
          0xee0e4c21, 0x822fef59, 0x4f0e0e20, 0x232feff8,
          0x1f8eaf20, 0x73af6fa8, 0x37ceffa0, 0x5bef2f80,
          0x23eed7e0, 0x4fcf0f94, 0x29fec3c0, 0x45df1f9e,
          0x2cf6c9d0, 0x40d7179b, 0x2e72ccd8, 0x42539399,
          0x2f30ce5c, 0x4311d198, 0x2f91cf1e, 0x43b07098,
          0xfbd9678f, 0x97f8384c, 0x91fdb3c7, 0xfddc1c26,
          0xa4efd9e3, 0xc8ce0e13, 0xbe66ecf1, 0xd2478709,
          0x673a5e48, 0x0b1bdbd0, 0x0b948714, 0x67b575bc,
          0x3dc3ebba, 0x51e2228a, 0xf2f075dd, 0x9ed11145,
          0x417112de, 0x2d5090f6, 0xcca9096f, 0xa088487b,
          0x8a4584b7, 0xe664a43d, 0xa933c25b, 0xc512d21e,
          0xb888e12d, 0xd4a9690f, 0x644d58a6, 0x086cacd3,
          0xde372c53, 0xb216d669, 0x830a9629, 0xef2beb34,
          0x798c6324, 0x15ad6dce, 0x04cf99a2, 0x68ee2eb3]

def _8To32(x32):
    return ((x32[0] * 256 + x32[1]) * 256 + x32[2]) * 256 + x32[3]    

def _32To8(x32):
    return [(x32 >> 8 * i) & 0xff for i in reversed(range(4))]

def _32To128(x32):
    return ((x32[0] * 256 ** 4 + x32[1]) * 256 ** 4 + x32[2]) * 256 ** 4 + x32[3]

def _128To32(x128):
    return [(x128 >> 32 * i) & 0xffffffff for i in reversed(range(4))]

def _192To32(x192):
    return [(x192 >> 32 * i) & 0xffffffff for i in reversed(range(6))]

def _256To32(x256):
    return [(x256 >> 32 * i) & 0xffffffff for i in reversed(range(8))]

def sigma(x128):
    return [(x128[0] << 7) & 0xffffff80  | (x128[1] >> 25),
            (x128[1] << 7) & 0xffffff80  | (x128[3] & 0x7f),
            (x128[0] & 0xfe000000)       | (x128[2] >> 7),
            (x128[2] << 25) & 0xfe000000 | (x128[3] >> 7)]

def mMult(m, t):
    """Multiply a 4x4 matrix by a transposed 1x4 vector in GF(2^8)"""

    def mult(p1, p2):
        """Multiply two polynomials in GF(2^8)"""
        p = 0
        while p2:
            if p2 & 1: p ^= p1
            p1 <<= 1
            if p1 & 256: p1 ^= 0x1d
            p2 >>= 1
        return p & 255    

    s = [0] * 4
    for i in range(4):
        for j in range(4):
            s[i] ^= mult(m[4 * i + j], t[j])
    return s

def f0(rk, x32):
    """F0 function"""
    t8 = _32To8(rk ^ x32)
    return _8To32(mMult(m0, [s0[t8[0]], s1[t8[1]], s0[t8[2]], s1[t8[3]]]))

def f1(rk, x32):
    """F1 function"""
    t8 = _32To8(rk ^ x32)
    return _8To32(mMult(m1, [s1[t8[0]], s0[t8[1]], s1[t8[2]], s0[t8[3]]]))

def gfn4(x32, n):
    """4-branch Generalized Feistel Network function"""
    t32 = x32[:]
    for i in range(n):
        t32[1] ^= f0(rk[2 * i], t32[0])
        t32[3] ^= f1(rk[2 * i + 1], t32[2])
        t32 = t32[1:] + t32[:1]
    return t32[3:] + t32[:3]

def gfn4i(x32, n):
    """4-branch Generalized Feistel Network inverse function"""
    t32 = x32[:]
    for i in reversed(range(n)):
        t32[1] ^= f0(rk[2 * i], t32[0])
        t32[3] ^= f1(rk[2 * i + 1], t32[2])
        t32 = t32[3:] + t32[:3]
    return t32[1:] + t32[:1]

def gfn8(x32, n):
    """8-branch Generalized Feistel Network function"""
    t32 = x32[:]
    for i in range(n):
        t32[1] ^= f0(rk[4 * i], t32[0])
        t32[3] ^= f1(rk[4 * i + 1], t32[2])
        t32[5] ^= f0(rk[4 * i + 2], t32[4])
        t32[7] ^= f1(rk[4 * i + 3], t32[6])     
        t32 = t32[1:] + t32[:1]
    return t32[7:] + t32[:7]

def setKey128(k128):
    """Generate round/whitening keys from 128-bit key"""
    k32 = _128To32(k128)
    for i in range(len(con128) - nrk):
        rk[i] = con128[i]
    l = gfn4(k32, 12)
    for i in range(nwk):
        wk[i] = k32[i]
    for i in range(nrk // 4):
        t32 = [r ^ s for r, s in zip(l, con128[4 * i + 24: 4 * i + 28])]
        l = sigma(l)
        if i % 2:
            t32 = [r ^ s for r, s in zip(t32, k32)]
        rk[4 * i:4 * i + 4] = t32

def setKey192(k192):
    """Generate round/whitening keys from 192-bit key"""
    k32 = _192To32(k192)
    kl = k32[:4]
    kr = k32[4:6] + [k32[0] ^ 0xffffffff] + [k32[1] ^ 0xffffffff]
    for i in range(len(con192) - nrk):
        rk[i] = con192[i]        
    l = gfn8(kl + kr, 10)
    ll, lr = l[:4], l[4:]
    kk = [r ^ s for r, s in zip(kl, kr)]
    for i in range(nwk):
        wk[i] = kk[i]
    for i in range(nrk // 4):
        if i % 4 == 0 or i % 4 == 1:
            t32 = [r ^ s for r, s in zip(ll, con192[4 * i + 40: 4 * i + 44])]
            ll = sigma(ll)
            if i % 2:
                t32 = [r ^ s for r, s in zip(t32, kr)]
        else:
            t32 = [r ^ s for r, s in zip(lr, con192[4 * i + 40: 4 * i + 44])]
            lr = sigma(lr)
            if i % 2:
                t32 = [r ^ s for r, s in zip(t32, kl)]                
        rk[4 * i:4 * i + 4] = t32     

def setKey256(k256):
    """Generate round/whitening keys from 256-bit key"""
    k32 = _256To32(k256)
    kl = k32[:4]
    kr = k32[4:]
    for i in range(len(con256) - nrk):
        rk[i] = con256[i]    
    l = gfn8(kl + kr, 10)
    ll, lr = l[:4], l[4:]
    kk = [r ^ s for r, s in zip(kl, kr)]    
    for i in range(nwk):
        wk[i] = kk[i]
    for i in range(nrk // 4):
        if i % 4 == 0 or i % 4 == 1:
            t32 = [r ^ s for r, s in zip(ll, con256[4 * i + 40: 4 * i + 44])]
            ll = sigma(ll)
            if i % 2:
                t32 = [r ^ s for r, s in zip(t32, kr)]
        else:
            t32 = [r ^ s for r, s in zip(lr, con256[4 * i + 40: 4 * i + 44])]
            lr = sigma(lr)
            if i % 2:
                t32 = [r ^ s for r, s in zip(t32, kl)]                
        rk[4 * i:4 * i + 4] = t32

def setKey(key, keySize):
    """Generate round/whitening keys"""
    global nr, nrk
    try:
        assert keySize in ksTable
    except AssertionError:
        print("Key size identifier not valid")
        sys.exit("ValueError")
    try:
        assert isinstance(key, int)
    except AssertionError:
        print("Invalid key")
        sys.exit("ValueError")
    try:
        assert key.bit_length() // 8 <= ksTable[keySize]
    except AssertionError:
        print("Key size mismatch")
        sys.exit("ValueError")
    nr = nrTable[keySize]
    nrk = nrkTable[keySize]
    if keySize == "SIZE_128":
        setKey128(key)
    elif keySize == "SIZE_192":
        setKey192(key)
    elif keySize == "SIZE_256":
        setKey256(key)
    else:
        sys.exit("Invalid key size identifier")

def enc(ptext):
    t32 = _128To32(ptext)
    t32[1] ^= wk[0]
    t32[3] ^= wk[1]
    t32 = gfn4(t32, nr)
    t32[1] ^= wk[2]
    t32[3] ^= wk[3]
    return _32To128(t32)

def dec(ctext):
    t32 = _128To32(ctext)
    t32[1] ^= wk[2]
    t32[3] ^= wk[3]
    t32 = gfn4i(t32, nr)
    t32[1] ^= wk[0]
    t32[3] ^= wk[1]
    return _32To128(t32)

if __name__ == "__main__":
    
    def checkTestVector(key, keySize, plaintext, ciphertext, nIter = 1000):
        testSuccess = True
        setKey(key, keySize)
        ks = ksTable[keySize] * 8
        try:
            assert  enc(plaintext) == ciphertext
        except AssertionError:
            print("Error in encryption")
            print("Resulted ciphertext: {:s}".format(ctext))
            print("Expected ciphertext: {:s}".format(ciphertext))
            testSuccess = False
        try:
            assert dec(enc(plaintext)) == plaintext
        except AssertionError:
            print("Error in decryption:")
            print("Recovered plaintext: {:s}".format(ptext))
            print("Expected plaintext: {:s}".format(plaintext))
            testSuccess = False         
        if not testSuccess:
            return False
        t1 = time()        
        for i in range(nIter):
            setKey(key, keySize)
            ctext = enc(plaintext)
        t2 = time()
        avg_elapsed_time = (t2 - t1) * 1000 / nIter        
        print("{:3d}-bit key test ok!".format(ksTable[keySize] * 8))
        print("Average elapsed time for 16-byte block ", end="")
        print("ECB-{0:3d} encryption: {1:0.3f}ms".format(ks, avg_elapsed_time))
        t3 = time()   
        for i in range(nIter):
            setKey(key, keySize)
            ptext = dec(ctext)
        t4 = time()
        avg_elapsed_time = (t4 - t3) * 1000 / nIter        
        print("{:3d}-bit key test ok!".format(ksTable[keySize] * 8))
        print("Average elapsed time for 16-byte block ", end="")
        print("ECB-{0:3d} decryption: {1:0.3f}ms".format(ks, avg_elapsed_time))        
        return True
    
    # The test vectors below are described in document "The 128-bit Blockcipher
    # CLEFIA Algorithm Specification" rev.1, June 1, 2007, Sony Corporation.

    ptext = 0x000102030405060708090a0b0c0d0e0f
    
    # Test vector for 128-bit key
    key1 = 0xffeeddccbbaa99887766554433221100
    ctext1 = 0xde2bf2fd9b74aacdf1298555459494fd
    
    # Test vector for 192-bit key
    key2 = 0xffeeddccbbaa99887766554433221100f0e0d0c0b0a09080
    ctext2 = 0xe2482f649f028dc480dda184fde181ad
    
    # Test vector for 256-bit key
    key3 = 0xffeeddccbbaa99887766554433221100f0e0d0c0b0a090807060504030201000
    ctext3 = 0xa1397814289de80c10da46d1fa48b38a
    
    try:
        assert checkTestVector(key1, "SIZE_128", ptext, ctext1) and \
               checkTestVector(key2, "SIZE_192", ptext, ctext2) and \
               checkTestVector(key3, "SIZE_256", ptext, ctext3)
    except AssertionError:
        print("At least one test failed")
        sys.exit(1)
    print("All tests passed!")
    sys.exit()

RSA implementation in Python

This Python script below implements the basic RSA encryption and decryption operations without any concern about padding or character encoding. Nevertheless, it has all the primitive machinery needed to encrypt and decrypt messages using the RSA public-key algorithm. The getprime() function is in charge of generating primes of the required size and, for this purpose, uses a probabilistic algorithm, the Miller-Rabin primality test. The inv() function, which calculates modular multiplicative inverses, gets a helping hand from the Extended Euclidean algorithm xgcd() function. The genRSA() function is pretty straightforward: it generates the standard public-key small exponent, 3 or 65537 (hex 10001), depending on the size of modulus (a small public exponent speeds up both the encryption and signature verification operations).

#!/usr/bin/python3
#
# Author: Joao H de A Franco (jhafranco@acm.org)
#
# Description: RSA implementation in Python 3
#
# Date: 2012-01-30
#
# License: Attribution-NonCommercial-ShareAlike 3.0 Unported
#          (CC BY-NC-SA 3.0)
#===========================================================
from random import randrange, getrandbits
from itertools import repeat
from functools import reduce

def getPrime(n):
    """Get a n-bit pseudo-random prime"""
    def isProbablePrime(n, t = 7):
        """Miller-Rabin primality test"""
        def isComposite(a):
            """Check if n is composite"""
            if pow(a, d, n) == 1:
                return False
            for i in range(s):
                if pow(a, 2 ** i * d, n) == n - 1:
                    return False
            return True
    
        assert n > 0
        if n < 3:
            return [False, False, True][n]
        elif not n & 1:
            return False
        else:
            s, d = 0, n - 1
            while not d & 1:
                s += 1
                d >>= 1
        for _ in repeat(None, t):
            if isComposite(randrange(2, n)):
                return False
        return True    
    
    p = getrandbits(n)
    while not isProbablePrime(p):
        p = getrandbits(n)
    return p

def inv(p, q):
    """Multiplicative inverse"""
    def xgcd(x, y):
        """Extended Euclidean Algorithm"""
        s1, s0 = 0, 1
        t1, t0 = 1, 0
        while y:
            q = x // y
            x, y = y, x % y
            s1, s0 = s0 - q * s1, s1
            t1, t0 = t0 - q * t1, t1
        return x, s0, t0      

    s, t = xgcd(p, q)[0:2]
    assert s == 1
    if t < 0:
        t += q
    return t

def genRSA(p, q):
    """Generate public and private keys"""
    phi, mod = (p - 1) * (q - 1), p * q
    if mod < 65537:
        return (3, inv(3, phi), mod)
    else:
        return (65537, inv(65537, phi), mod)    

def text2Int(text):
    """Convert a text string into an integer"""
    return reduce(lambda x, y : (x << 8) + y, map(ord, text))

def int2Text(number, size):
    """Convert an integer into a text string"""
    text = "".join([chr((number >> j) & 0xff)
                    for j in reversed(range(0, size << 3, 8))])
    return text.lstrip("\x00")

def int2List(number, size):
    """Convert an integer into a list of small integers"""
    return [(number >> j) & 0xff
            for j in reversed(range(0, size << 3, 8))]

def list2Int(listInt):
    """Convert a list of small integers into an integer"""
    return reduce(lambda x, y : (x << 8) + y, listInt)

def modSize(mod):
    """Return length (in bytes) of modulus"""
    modSize = len("{:02x}".format(mod)) // 2
    return modSize

def encrypt(ptext, pk, mod):
    """Encrypt message with public key"""
    size = modSize(mod)
    output = []
    while ptext:
        nbytes = min(len(ptext), size - 1)
        aux1 = text2Int(ptext[:nbytes])
        assert aux1 < mod
        aux2 = pow(aux1, pk, mod)
        output += int2List(aux2, size + 2)
        ptext = ptext[size:]
    return output

def decrypt(ctext, sk, p, q):
    """Decrypt message with private key
    using the Chinese Remainder Theorem"""
    mod = p * q
    size = modSize(mod)
    output = ""
    while ctext:
        aux3 = list2Int(ctext[:size + 2])
        assert aux3 < mod
        m1 = pow(aux3, sk % (p - 1), p)
        m2 = pow(aux3, sk % (q - 1), q)
        h = (inv(q, p) * (m1 - m2)) % p
        aux4 = m2 + h * q
        output += int2Text(aux4, size)
        ctext = ctext[size + 2:]
    return output

if __name__ == "__main__":

    from math import log10
    from time import time

    def printHexList(intList):
        """Print ciphertext in hex"""
        for index, elem in enumerate(intList):
            if index % 32 == 0:
                print()            
            print("{:02x}".format(elem), end = "")
        print()

    def printLargeInteger(number):
        """Print long primes in a formatted way"""
        string = "{:02x}".format(number)
        for j in range(len(string)):
            if j % 64 == 0:
                print()
            print(string[j], end = "")
        print()

    def testCase(p, q, msg, nTimes = 1):
        """Execute test case: generate keys, encrypt message and
           decrypt resulting ciphertext"""
        print("Key size: {:0d} bits".format(round(log10(p * q) / log10(2))))
        print("Prime #1:", end = "")
        printLargeInteger(p)
        print("Prime #2:", end = "")
        printLargeInteger(q)
        print("Plaintext:", msg)
        pk, sk, mod = genRSA(p, q)
        ctext = encrypt(msg, pk, mod)
        print("Ciphertext:", end = "")
        printHexList(ctext)
        ptext = decrypt(ctext, sk, p, q)
        print("Recovered plaintext:", ptext, "\n")

    # First test: RSA-129 (see http://en.wikipedia.org/wiki/RSA_numbers#RSA-129)
    p1 = 3490529510847650949147849619903898133417764638493387843990820577
    p2 = 32769132993266709549961988190834461413177642967992942539798288533
    testCase(p1, p2, "The Magic Words are Squeamish Ossifrage", 1000)
  
    # Second test: random primes (key size: 512 to 4096 bits)
    for n in [256, 512, 1024, 2048]:    
        t1 = time()
        p5 = getPrime(n)
        t2 = time()
        print("Elapsed time for {:0d}-bit prime ".format(n), end = "")
        print("generation: {:0.3f} s".format(round(t2 - t1, 3)))
        t3 = time()
        p6 = getPrime(n)
        t4 = time()
        print("Elapsed time for {:0d}-bit prime ".format(n), end = "")
        print("generation: {:0.3f} s".format(round(t4 - t3, 3)))
        testCase(p5, p6, "It's all greek to me")

The first test uses the famous 425-bit key pair belonging to RSA-129, a challenge published in Martin Gardner’s column in Scientific American back in 1977. RSA-129 owes its name to the 129-digit number that expresses the product of two (unknown at the time) primes.

RSA-129 = 1143816257578888676692357799761466120102182967212423625625618429 35706935245733897830597123563958705058989075147599290026879543541

This challenge, posed by the inventors of the RSA algorithm (Rivest, Shamir & Adleman), was solved only in 1994. The decrypted plaintext was the phrase “The Magic Words are Squeamish Ossifrage”. The second test was done with pairs of randomly generated primes. Their key sizes (as a matter of fact, the size of the modulus) are 512, 1024 and 2048 bits (the key size most used in SSL certificates is 1024 bits). The output of this Python script is shown below.

Key size: 425 bits
Prime #1:
87c296ed480f9ab17885decd31197d617779c0dac70c3234996e1
Prime #2:
4fa84812157119acc8ecca98c404b2e5ee24ce18f60ea818091895
Plaintext: The Magic Words are Squeamish Ossifrage
Ciphertext:
000100d00ec048028f8de9998bfcdb271ad5d34ab70a3990ebdfdd30a32ae15e
1ae01ccda1945493cc02be7d739ded0c56d8cd9c996ee8
Recovered plaintext: The Magic Words are Squeamish Ossifrage 

Elapsed time for 256-bit prime generation: 0.108 s
Elapsed time for 256-bit prime generation: 0.018 s
Key size: 511 bits
Prime #1:
7a29e838f6afe351c24b0dfee0b12d32e092c6ec129382b75b224f749e9785c9
Prime #2:
caffe728fa460b4f009a7647bdbb0763868206942a4ced79aefd36a6b05f539b
Plaintext: It's all greek to me
Ciphertext:
00005741fdf3f6f20fc16524c9a1c80f62e7a293aaff5238f319692e4c469443
10d2272ed76d6a9a884372c2eae03e2d3b8c10df04353d459659a3f8c3e63be3
2c6d
Recovered plaintext: It's all greek to me 

Elapsed time for 512-bit prime generation: 0.162 s
Elapsed time for 512-bit prime generation: 0.576 s
Key size: 1015 bits
Prime #1:
37f8ab0bbdeca6af2a8a8874844472add93745cb98a6224175c6528c1e3c2f26
5ef0fcba9c6044ecd3c98aa22981e1f75f8779ba5f72fa33ab4da89d4abfc357
Prime #2:
2aabce36a668c2f8e161521238d2be8541d286e612b6572e4c1531f8218dfe0e
08de2e177b02c2ab1a76056eca13920e2847bfdffffdb6b2b9fe874849c3637
Plaintext: It's all greek to me
Ciphertext:
00006d1c01974161e2c4710ec423076bbab8107cf31d2d2ecfd6127ca6dc142e
aa7f37d93281222684227892c82575ef93491ffa802003efac9085e3e8076000
63780dd7ea03b1418e7cdb17dec38b533913d2e11f37de1f364dd96e559add5b
a205c367ed6d7270e1ab02da80522264238e3cb93aff81e53eeec7910d5013c8
a4
Recovered plaintext: It's all greek to me 

Elapsed time for 1024-bit prime generation: 7.189 s
Elapsed time for 1024-bit prime generation: 8.921 s
Key size: 2047 bits
Prime #1:
7843d991e4fa29466f4e0fa385f0b8612db5d51c0f09f0fca801339155cd1ca2
fde3b810569de371e99e8d1e3d127a9b4d7f07944d64c7d9da973252b73dd1ed
847183c1a855b65817411cfea22c05b58764d5bb016770feef93d6cb5e6274e3
fce2c5fb251bb41a8ba879cd5a2d6755c591e921c6aa08327ef6b6e8f1cf1a1b
Prime #2:
e751af1ed509106fdd814d52778da0cf6998b768bafd16a10ec91c4becf05856
9e96007661d4525cb0b40b20247aa879657c4550d3d58a74ab23666b504febaf
5de5b40f3cc7b8c766e5d426fa37d9832b91fef718e4f117984168a2eeb087bc
8ad09364d9f57c6534d5332eb78eacf9f5cdaa12f95aee44de0d1df475ec52ff
Plaintext: It's all greek to me
Ciphertext:
000056ce1bab323cc03876dc786cdc9f53e406a4e321a2e8233ed648f93d0026
3b9504faa72bfbf2d97b166f320344596bb2634addd2766ef9e4816841a76b95
ba3c5b899bd559222162d5580924ef073ba7118023710cebc901a20f51603acb
5a797e84ea7cf853d3b2031f3263879d3afb3d05bd74f15219b877753b3182ef
a93cdfc015db601b8fd93540e1e50d1b9676e152769268039484cdd821f9291b
93f408d1e8c851fb79b143a80d7ef03629f3c11c81f8188a35076df9a88df1ba
5da77de10857a9c7238122f3818c4011a6cbe4e702c8b481789edca59ab74606
3608ca5c98f6126ddec1a3eb972021b9297c7d4d2c38f50f8105745eb3110977
f76a
Recovered plaintext: It's all greek to me

I will comment on code details in another post.