Category: Matrix Cipher

Blog Entry (c) Wednesday, August 13, 2025, by James Pate Williams, Jr. Exercises from an Online Textbook

Exercises Szabo OstlundDownload
#include <complex>
#include <vector>

class CmpLinearAlgebra
{
public:
    static void CmpPrintMatrix(
        int m, int n,
        std::vector<std::vector<std::complex<double>>>& Ac);
    static void CmpAddition(
        size_t m, size_t n,
        std::vector<std::vector<std::complex<double>>>& A,
        std::vector<std::vector<std::complex<double>>>& B,
        std::vector<std::vector<std::complex<double>>>& C);
    static void CmpSubtraction(
        size_t m, size_t n,
        std::vector<std::vector<std::complex<double>>>& A,
        std::vector<std::vector<std::complex<double>>>& B,
        std::vector<std::vector<std::complex<double>>>& C);
    static void CmpMultiply(
        size_t m, size_t n, size_t p,
        std::vector<std::vector<std::complex<double>>>& A,
        std::vector<std::vector<std::complex<double>>>& B,
        std::vector<std::vector<std::complex<double>>>& C);
    static void CmpAnticommutator(
        size_t n,
        std::vector<std::vector<std::complex<double>>>& A,
        std::vector<std::vector<std::complex<double>>>& B,
        std::vector<std::vector<std::complex<double>>>& C);
    static void CmpCommutator(
        size_t n,
        std::vector<std::vector<std::complex<double>>>& A,
        std::vector<std::vector<std::complex<double>>>& B,
        std::vector<std::vector<std::complex<double>>>& C);
    static void CmpAdjoint(
        size_t m, size_t n,
        std::vector<std::vector<std::complex<double>>>& Ac,
        std::vector<std::vector<std::complex<double>>>& Ad);
    static std::complex<double> CmpDeterminant(
        bool& failure, int n,
        std::vector<std::vector<std::complex<double>>>& A);
    static bool CmpGaussianElimination(
        int m, int n,
        std::vector<std::vector<std::complex<double>>>& A,
        std::vector<std::complex<double>>& b,
        std::vector<std::complex<double>>& x,
        std::vector<size_t>& pivot);
    static bool CmpGaussianFactor(
        int n, std::vector<std::vector<std::complex<double>>>& M,
        std::vector<size_t>& pivot);
    static bool CmpGaussianSolution(
        int n, std::vector<std::vector<std::complex<double>>>& M,
        std::vector<std::complex<double>>& b,
        std::vector<std::complex<double>>& x,
        std::vector<size_t>& pivot);
    static bool CmpSubstitution(
        int m, int n, std::vector<std::vector<std::complex<double>>>& M,
        std::vector<std::complex<double>>& b,
        std::vector<std::complex<double>>& x,
        std::vector<size_t>& pivot);
    static bool CmpInverse(
        int n, std::vector<std::vector<std::complex<double>>>& M,
        std::vector<std::vector<std::complex<double>>>& Mi);
    static void CmpCharPolyAndAdjoint(
        int n,
        std::vector<std::vector<std::complex<double>>>& C,
        std::vector<std::vector<std::complex<double>>>& I,
        std::vector<std::vector<std::complex<double>>>& M,
        std::vector<std::vector<std::complex<double>>>& adjoint,
        std::vector<std::complex<double>>& a);
    static void CmpMatrixKernel(
        int m, int n,
        std::vector<std::vector<std::complex<double>>>& M,
        std::vector<std::vector<std::complex<double>>>& X,
        size_t& r);
    static void CmpMatrixImage(
        int m, int n,
        std::vector<std::vector<std::complex<double>>>& M,
        std::vector<std::vector<std::complex<double>>>& N,
        std::vector<std::vector<std::complex<double>>>& X,
        int rank);
};

#include <vector>

class DblLinearAlgebra
{
public:
    static void DblPrintMatrix(
        int m, int n, std::vector<std::vector<double>>& A);
    static void DblAddition(
        size_t m, size_t n,
        std::vector<std::vector<double>>& A,
        std::vector<std::vector<double>>& B,
        std::vector<std::vector<double>>& C);
    static void DblSubtraction(
        size_t m, size_t n,
        std::vector<std::vector<double>>& A,
        std::vector<std::vector<double>>& B,
        std::vector<std::vector<double>>& C);
    static void DblMultiply(
        size_t m, size_t n, size_t p,
        std::vector<std::vector<double>>& A,
        std::vector<std::vector<double>>& B,
        std::vector<std::vector<double>>& C);
    static void DblAnticommutator(
        size_t n,
        std::vector<std::vector<double>>& A,
        std::vector<std::vector<double>>& B,
        std::vector<std::vector<double>>& C);
    static void DblCommutator(
        size_t n,
        std::vector<std::vector<double>>& A,
        std::vector<std::vector<double>>& B,
        std::vector<std::vector<double>>& C);
    static double DblDeterminant(
        bool& failure, int n,
        std::vector<std::vector<double>>& A);
    static bool DblGaussianElimination(
        int m, int n, std::vector<std::vector<double>>& A,
        std::vector<double>& b, std::vector<double>& x,
        std::vector<size_t>& pivot);
    static bool DblGaussianFactor(
        int n, std::vector<std::vector<double>>& M,
        std::vector<size_t>& pivot);
    static bool DblGaussianSolution(
        int n, std::vector<std::vector<double>>& M,
        std::vector<double>& b, std::vector<double>& x,
        std::vector<size_t>& pivot);
    static bool DblSubstitution(
        int n, std::vector<std::vector<double>>& M,
        std::vector<double>& b, std::vector<double>& x,
        std::vector<size_t>& pivot);
    static bool DblInverse(
        int n, std::vector<std::vector<double>>& M,
        std::vector<std::vector<double>>& A);
    static void DblCharPolyAndAdjoint(
        int n,
        std::vector<std::vector<double>>& C,
        std::vector<std::vector<double>>& I,
        std::vector<std::vector<double>>& M,
        std::vector<std::vector<double>>& adjoint,
        std::vector<double>& a);
    static void DblMatrixKernel(
        int m, int n,
        std::vector<std::vector<double>>& M,
        std::vector<std::vector<double>>& X,
        size_t& r);
    static void DblMatrixImage(
        int m, int n,
        std::vector<std::vector<double>>& M,
        std::vector<std::vector<double>>& N,
        std::vector<std::vector<double>>& X,
        int rank);
};

// Exercises from "Modern Quantum Chemistry An Introduction to Advanced
// Electronic Structure Theory" by Attila Szabo and Neil S. Ostlund
// https://chemistlibrary.wordpress.com/wp-content/uploads/2015/02/modern-quantum-chemistry.pdf
// Program (c) Tuesday, August 12, 2025 by James Pate Williams, Jr.

#include <complex>
#include <iomanip>
#include <iostream>
#include <vector>
#include "DblLinearAlgebra.h"
#include "CmpLinearAlgebra.h"

int main()
{
	double AArcb[3][3] = { { 2, 3, -1 }, { 4, 4, -3 }, { -2, 3, -1 } };
	double AArso[3][3] = { { 1, 1, 0 }, { 1, 2, 2 }, { 0, 2, -1 } };
	double BBrso[3][3] = { { 1, -1, 1 }, { -1, 0, 0 }, { 1, 0, 1} };
	double BBr[3][3] = { { 1, -1, 1 }, { -1 , 0, 0 }, { 1, 0, 1 } };
	double AAcr[3][3] = { { 1, 2, 3 }, { 4, 5, 6 }, { 7, 8, 9 } };
	double AAci[3][3] = { { 1, 1, 2 }, { 3, 0, 1 }, { 0, 2, 4 } };
	double BBcr[3][3] = { { 1, 0, 1 }, { 1 , 1, 0 }, { 0, 1, 1 } };
	double BBci[3][3] = { { 1, 2, 3 }, { 4, 5, 6 }, { 7, 8, 9 } };
	int m = 3, n = 3, p = 3;

	std::vector<double> br(3);
	std::vector<size_t> pivot(3);
	std::vector<std::vector<double>> Arcb(3, std::vector<double>(3));
	std::vector<std::vector<double>> Arso(3, std::vector<double>(3));
	std::vector<std::vector<double>> Brso(3, std::vector<double>(3));
	std::vector<std::vector<double>> Br(3, std::vector<double>(3));
	std::vector<std::vector<double>> Cr(3, std::vector<double>(3));
	std::vector<std::vector<double>> Ai(3, std::vector<double>(3));
	std::vector<std::vector<double>> Ari(3, std::vector<double>(3));
	std::vector<std::vector<std::complex<double>>> Ac(3,
		std::vector<std::complex<double>>(3));
	std::vector<std::vector<std::complex<double>>> Bc(3,
		std::vector<std::complex<double>>(3));
	std::vector<std::vector<std::complex<double>>> Cc(3,
		std::vector<std::complex<double>>(3));
	std::vector<std::vector<std::complex<double>>> Dc(3,
		std::vector<std::complex<double>>(3));
	std::vector<std::vector<std::complex<double>>> Ec(3,
		std::vector<std::complex<double>>(3));
	std::vector<std::vector<std::complex<double>>> Fc(3,
		std::vector<std::complex<double>>(3));
	std::vector<std::vector<std::complex<double>>> Gc(3,
		std::vector<std::complex<double>>(3));
	
	for (int i = 0; i < m; i++)
	{
		for (int j = 0; j < p; j++)
		{
			Arcb[i][j] = AArcb[i][j];
			Arso[i][j] = AArso[i][j];
			Brso[i][j] = BBrso[i][j];
			Ac[i][j]._Val[0] = AAcr[i][j];
			Ac[i][j]._Val[1] = AAci[i][j];
		}
	}

	for (int i = 0; i < p; i++)
	{
		for (int j = 0; j < n; j++)
		{
			Br[i][j] = BBr[i][j];
			Bc[i][j]._Val[0] = BBcr[i][j];
			Bc[i][j]._Val[1] = BBci[i][j];
		}
	}
	
	DblLinearAlgebra::DblMultiply(3, 3, 3, Arcb, Br, Cr);
	std::cout << "Ar * Br = Cr Conte & de Boor" << std::endl;
	DblLinearAlgebra::DblPrintMatrix(3, 3, Cr);
	std::cout << std::endl;
	CmpLinearAlgebra::CmpMultiply(3, 3, 3, Ac, Bc, Cc);
	std::cout << "Ac * Bc = Cc" << std::endl;
	CmpLinearAlgebra::CmpPrintMatrix(3, 3, Cc);
	std::cout << std::endl;
	// Exercise 1.2
	std::cout << "Exercise 1.2 page 5 Commutator" << std::endl;
	DblLinearAlgebra::DblCommutator(3, Arso, Brso, Cr);
	DblLinearAlgebra::DblPrintMatrix(3, 3, Cr);
	std::cout << std::endl;
	std::cout << "Exercise 1.2 page 5 Anticommutator" << std::endl;
	DblLinearAlgebra::DblAnticommutator(3, Arso, Brso, Cr);
	DblLinearAlgebra::DblPrintMatrix(3, 3, Cr);
	std::cout << std::endl;
	CmpLinearAlgebra::CmpAdjoint(3, 3, Cc, Dc);
	std::cout << "Exercise 1.3 page 6 Cc adjoint" << std::endl;
	CmpLinearAlgebra::CmpPrintMatrix(3, 3, Dc);
	std::cout << std::endl;
	CmpLinearAlgebra::CmpAdjoint(3, 3, Ac, Ec);
	CmpLinearAlgebra::CmpAdjoint(3, 3, Bc, Fc);
	CmpLinearAlgebra::CmpMultiply(3, 3, 3, Fc, Ec, Gc);
	std::cout << "Exercise 1.3 page 6 Bc adjoint * Ac adjoint" << std::endl;
	CmpLinearAlgebra::CmpPrintMatrix(3, 3, Gc);
	std::cout << std::endl;
	std::cout << "Ar matrix" << std::endl;
	DblLinearAlgebra::DblPrintMatrix(3, 3, Arcb);
	bool inv = DblLinearAlgebra::DblInverse(n, Arcb, Ai);
	std::cout << std::endl;
	std::cout << "Ar Conte & de Boor inverse = " << inv << std::endl;
	DblLinearAlgebra::DblPrintMatrix(3, 3, Ai);
	std::cout << std::endl;
	std::cout << "Ar * Ar inverse" << std::endl;
	DblLinearAlgebra::DblMultiply(3, 3, 3, Arcb, Ai, Ari);
	DblLinearAlgebra::DblPrintMatrix(3, 3, Ari);
	std::cout << std::endl;
	std::cout << "Ac" << std::endl;
	CmpLinearAlgebra::CmpPrintMatrix(3, 3, Ac);
	std::cout << std::endl;
	inv = CmpLinearAlgebra::CmpInverse(3, Ac, Bc);
	std::cout << "Ac inverse = " << inv << std::endl;
	CmpLinearAlgebra::CmpPrintMatrix(3, 3, Bc);
	CmpLinearAlgebra::CmpMultiply(3, 3, 3, Ac, Bc, Cc);
	std::cout << std::endl;
	std::cout << "Ac * AC inverse" << std::endl;
	CmpLinearAlgebra::CmpPrintMatrix(3, 3, Cc);
}

Some Statistics for Another Matrix Cipher by James Pate Williams, Jr. BA, BS, MSwE, PhD

We first encipher the string “This is a test of the emergency broadcasting system!” which is a English language sample of length 52 ASCII characters.

New Matrix Cipher Encrypt 0

Below is a histogram of the plaintext characters.

New Matrix Cipher Encrypt 1

Here are the counts of the different plaintext characters and the statistic known as the index of coincidence. English has an index of coincidence of approximately 0.065, so this short sample is in that ballpark at 0.06067.

New Matrix Cipher Encrypt 2

Next we display part of the key material (upper triangular matrix elements), the ASCII encoded plaintext and the last column is the resulting ciphertext.

New Matrix Cipher Encrypt 3

We now display a histogram of the ciphertext.

New Matrix Cipher Decrypt 1

New Matrix Cipher Decrypt 2

The index of coincidence is 0. Each of the 52 plaintext characters map to a different ciphertext number in the range 2 to 977 inclusive. Now we show the decryption matrix, ciphertext, and plaintext.

New Matrix Cipher Decrypt 3

We now say a few more words about the cryptanalysis of this matrix cipher. Since the PRNG has a large key space, a direct brute force attack on the PRNG would probably be futile in a milieu without a quantum computer. We can guess values of N and we know that key matrix inverse has n * (n + 1) / 2 elements instead of n * n elements since it is an upper triangular matrix. In our sample the key matrix inverse has 52 * 53 / 2 = 26 * 53 = 1,378 matrix elements. Each matrix element is in the inclusive range 0 to 996 or 997 values. So we would need to brute force test N * n * (n + 1) / 2 possible key matrix inverses. In our case the number is 997 * 1,378 = 1,373,866 which is not a very large number by cryptanalytic standards but how many of those decryptions would make perfectly good sense? A reasonably adept adversary would not reuse the PRNG key and N for a new 52 ASCII character string so each message would require a new  cryptanalytic attack.

Using the upper triangular nature of the key matrix inverse and a further assumption that the plaintext is in the range 32 to 127 or 96 different values, we can reduce the brute force attack to 96 * 1,378 = 132,288 possibilities. Again the adversary could somewhat thwart these cryptanalytic efforts by choosing a much larger N and n.

However, all of this mental masturbation amounts to a mute point since no modern and sane adversary would utilize such an easy cipher to break. Just use a one time pad based on the ANSI X9.17 PRNG utilizing triple-AES instead of triple-DES with three 256-bit E-D-E keys or an astounding 768-bits of key material for the cipher core alone plus 256 bits in additional secret information. The final keyspace of such a scheme would be 1024-bits!

 

Another Matrix Cipher by James Pate Williams, Jr. BA, BS, MSwE, PhD

This is perhaps an improvement on the matrix cipher of a previous blog post of mine. In that post I introduced a matrix cipher whose keys were generated by selection of a seed such that 1 <= seed <= 2147483647, a number N such that 2 <= N <= 1000, and plaintext of length n such that 1 <= n <= N -1.

This matrix cipher relies on the ANSI X9.17 pseudorandom number generator (PRNG) of 5.11 Algorithm of the Handbook of Applied Cryptography by Alfred J. Menezes, et al. The PRNG uses triple-DES with a potential 168-bit (56 * 3 = 168) key space using E-D-E (Encryption key 1 – Decryption key 2 – Encryption key 3). Also, a 64-bit date related number and a 64-bit random seed are needed to initialize the PRNG.

The key space for the algorithm is (168 + 128) bits which is 296 bits. Here is the encryption and decryption of the ASCII ten characters string “ATTACK NOW”.

New Matrix Cipher 0

New Matrix Cipher 1

New Matrix Cipher 2

New Matrix Cipher 3

New Matrix Cipher 4

The first step in the cryptanalysis of this cipher would be to determine the modulus of the matrix and vector calculations N. I don’t know how many ciphertexts would be necessary to perform this task. From the preceding known ciphertext we find that N is at least 991. From traffic analysis we may have determined that the maximum value of N is 1000. That means would we only need to try 10 values of N.

 

A Matrix Cipher by James Pate Williams, Jr., BA, BS, MSwE, PhD

Suppose you have a n vector of ASCII encoded characters, arbitrarily choose 1 <= n <= 1000. Choose a modulus N such that 128 <= N <= 1000. Also choose a pseudo-random number seed 1 <= s <= 2147483647. Next find a random n x n matrix that is invertible by Gaussian elimination over the integer field consisting of N elements. Suppose this matrix is M and its tridiagonal form is M’. Now suppose the plaintext is the n vector P and the ciphertext is the n vector C then we have for encryption:

C = M’P

Further assume the inverse of M’ is N’. For decryption we use the equation:

P = N’C

Where

M’N’ = N’M’ = I

Such that I is the n x n identity matrix.

This cipher is related to the classic Hill Cipher. This cipher is polyalphabetic. We show the results of one encryption and decryption using 10 ASCII ‘A’ characters, N = 999, and s = 1. As you can see each occurrence of the letter ‘A’ which is encoded as the decimal number 65 leads to different integer in the range 0 to 998 which has a maximum of 10 bits. The key consists of the 100 integers in the original 10 x 10 matrix. The application was implemented in C# using a Gaussian elimination over a number field algorithm from Henri Cohen’s A Course in Computational Algebraic Number Theory.

Matrix Cipher 0

Matrix Cipher 1

Matrix Cipher 2

Machine Cryptanalysis of Basic Cryptosystems by James Pate Williams, Jr. BA, BS, MSwE, PhD

In Winter and Spring 2018 I wrote a simple C# computer program to perform machine cryptanalysis of the following basic (elementary and easily breakable) cryptosystems:

  1. Affine Cipher Operating on Monographs and Digraphs
  2. Matrix Cipher
  3. Mono-alphabetic Cipher
  4. n-Rotor with Shifting Polyalphabetic Cipher

The key ingredient in this program is a relatively extensive English language dictionary.

This slideshow requires JavaScript.