# Is this asymmetric (public key) cryptosystem based on a 16x16 s-table safe and useful?

May be this is absolutely off-topic, but here is. The cryptosystem description follows. Any hints of what to do with it, or flaws found are welcomed.

This description is here as well.

We will use a substitution 16x16 table, this is one of them:

 7  9 13 10 15  2  0  6  3 12  8  4  1  5 14 11
1 15  6  3  9  4 11 13 10  5 14  2  7 12  8  0
3  0 12  1 11  8  9  5  7 13  2 14 10  6  4 15
4  6 15  8 13  1  5  9 14 11 10  7  2  0  3 12
0  3  8 15 10 12  7 14  9  2 13  5 11  4  6  1
10 11  5  7  0 14 15 12  1  6  4  8  3 13  2  9
5 14 10 13  8 11  4  3  6  1 15  0 12  7  9  2
15  1  4  0  7  6 10  2 11 14  5 13  9  8 12  3
12  8  3  6 14  0  2 10 13  7  9 11  5  1 15  4
13  2  7  5  4  9  8  1 12  3  0 15  6 10 11 14
6  4  1 12  2 15 14  7  5 10 11  9 13  3  0  8
9  7  2 11  1 13  3  4  0  8 12  6 15 14  5 10
11 10 14  9  3  5  1  8 15  4  6 12  0  2 13  7
14  5 11  2 12 10  6  0  4 15  1  3  8  9  7 13
8 12  0  4  5  3 13 11  2  9  7 10 14 15  1  6
2 13  9 14  6  7 12 15  8  0  3  1  4 11 10  5


to define a function $$c=f(a,b)$$, where $$c$$ is the element in the $$a$$-th row and $$b$$-th column.

The following three properties hold:

$$f(f(a,b),c)\neq f(a,f(b,c))$$ -- non-associativity in general

$$f(a,b)\neq f(b,a)$$ -- non-commutativity in general

$$f(f(a,b),f(c,d))=f(f(a,c),f(b,d))$$ -- restricted commutativity

Next we define a list of $$N$$ integers in the range $$[0,15]$$ to meet the size required. For $$256$$ bits we need $$N=64$$. This list can be interpreted as a $$64$$ digit base-$$16$$ number.

Next we define a mixing procedure of elements of this kind, $$t$$ and $$k$$, N-element lists of numbers in the integer range $$[0,15]$$.

But before, we define a deterministic fixed sequence of numbers in the interval $$[0,N-1]$$, in a pseudorandom way. frist_seq gets the first element of the sequence and next_seq the next one.

The mixing procedure is:

function m(t,k) returns r

//one-to-one mixing of k and t
for i in 0..N-1
r[i] <- f(t[i],k[i])
end for
i <- first_seq
for M number of applications of f -- 4096 for example
// accumulative mixing of r with itself
j <- next_seq
r[j]<-f(r[j],r[i])
i <- j
end for
//one-to-one mixing of k and r
for i in 0..N-1
r[i] <- f(r[i],k[i])
end for

return r


The function m is neither associative nor commutative, and meets the restricted commutativity property:

$$m(m(a,b),m(c,d))=m(m(a,c),m(b,d))$$

The computationally hard problem proposed is:

in $$c=m(t,k)$$, knowing $$c$$ and $$t$$, find $$k$$.

About differential and linear cryptoanalysis. If we take each row and each column as a 4-bit lenght substitution table, i turns out that the probability of a linear equation holding is at most 14/16 and more or less is the same for differential characteristics. As we're doing 4096 iteration on the s-table, despite 14/16 is a high probability, $$log_2((14/16)^{4096})$$ is by far lesser than $$2^{-256}$$. So there's no apparent attack feasible here.

Now lets define the secret agreement and the digital signature using the mixing function $$m$$. To put it more clear we will use the following notation:

$$m(a,b)$$ is written as $$(a b)$$

$$m(m(a,b),m(c,d))$$ is written as $$(a b)(c d)$$

$$m(m(a,b),c\dots)$$ is written as $$(a b c\dots)$$

For the secret agreement the procedure is the following:

Both Alice and Bob agree on some constant $$C$$. Alice chooses a random key $$K$$, and Bob does the same choosing a random key $$Q$$. Alice sends to Bob $$(C K)$$, Bob sends to Alice $$(C Q)$$. Alice computes using bob sent value $$(C Q)(K C)$$, and Bob does the same and computes $$(C K)(Q C)$$.

By the property of restricted commutatitvity $$(C Q)(K C)=(C K)(Q C)$$

For the signature the procedure is the following:

Alice, the signer, chooses a public value $$C$$ and a two random keys $$K$$,$$Q$$. Its credentials are $$C$$, $$(C K)$$ and $$(K Q)$$. To sing a value, $$H$$, Alice computes $$S=(H Q)$$.

Bob needs to verify if Alice has signed $$H$$. Computes $$(C K)(H Q)$$ and $$(C H)(K Q)$$. Both values must be equal due to restricted commutativity if $$(H Q)$$ is a valid signature from Alice.

In this document there also a proof of correctness of the algorithm.

• A random substitution table would not have the property $f(f(a,b),f(c,d))=f(f(a,c),f(b,d))$. This one does. How was it constructed? This short C code verifies that the stated property is met.
– fgrieu
Mar 14, 2021 at 17:18
• I've built a search program. This is the reason is only 16x16. The number of random latin squares, the s-table is one is intractable, but with various strategies the process can be made faster. Discarding branches that already fails then property and filling cells that can have only a value due to the property. Mar 15, 2021 at 11:13
• @poncho I've claimed it. And it's a must. Indeed if somebody else finds a function like this that's safe the cryptosystem can be applied the same. I've failed several times, so I opted for the s-table. Mar 15, 2021 at 11:23
• Regardless of security, does it perform well? Does it perform well when it scales up to 128-bit security or 192-bit security? crypto.meta.stackexchange.com/a/1482/36960 Mar 15, 2021 at 11:36
• @DannyNiu : For the 256-bit size in the question, I think the computational cost is much lower than one point multiplication on curve25519. A convincing security argument is all I see needed to make this attractive in key exchange, signature, and asymmetric encryption.
– fgrieu
Mar 15, 2021 at 12:17

## 3 Answers

Your mixing function is an isomorphism of some group that's representable using matrices, and it's generally not difficult to find inverses of matrices (even ones with the determinant value 0) unless the dimension is too high. It took months for my hair to re-grow after this.

For example, the permutation representation of the 0'th row $$p(x) = f(0,x)$$ can be represented as the following matrix:

$$\begin{pmatrix} 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \\ 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \end{pmatrix}$$

This is simply setting the $$f(0,i)$$'th column of $$i$$'th row to 1 and every other cells to 0.

And as shown, this is a quite large matrix. It may be possible to simplify it, but the degree to which it can be simplified depends on the function $$f$$. And when $$f$$ is used in the mixing function $$m$$, then expanded matrix may be as unpractical to invert as XSL attack on AES.

The only matrix-related problem known to be hard with moderately low dimension is "lattice-reduction", which is used in post-quantum schemes such as Kyber, NTRU, Saber, Dilithium, Falcon, etc. If you do manage to build a cryptographically meaningful and efficient commutative function that require cryptanalysis at unpractically high matrix dimension, then you definitely should post it on IACR ePrint, and have it peer-reviewed.

• The matrix algebra I learned has coefficients in a field. If we could show that any $f$ with "latin square" and "restricted commutativity" must be of the form $f(a,b)=g^{-1}(g(a)u+g(b)v+w)$ in some instance of $\mathbb F_{2^4}$ for some fixed $u,v,w$ and bijection $g$ in that field, it would strengthen this answer's counter-argument of security.
– fgrieu
Mar 15, 2021 at 12:49
• Your answer is rather vague. But you do need to find the inverse of the isomorphism, or are you telling the isomorphism (mix function) can be represented by a matrix of some form? Anyway, perhaps with a more complex basic f the theorichal part is still useful. Mar 15, 2021 at 13:31
• @daniel I'm saying that the isomorphism can be represented by a matrix of some form. Right now, I'm not yet able to list some of the elements as matrix, and it's entirely possible that to represent an element of your mixing function require as much space as that needed for a sporadic "Monster" group. That is to say, my envised cryptanalysis may be as unpractical as XSL attack on AES. Mar 16, 2021 at 2:29
• @DannyNiu 64 4 bit values can be represented in a vector of 1024 bits. Each step is done by a matrix of size $1024^2$. The product of all these matrices represents the mixing. With gaussian elimination is manageable to recover the key knowing the plaintext and the result. Mar 16, 2021 at 10:49
• Actually I think pairs are needed $f(a,b)=(a',b)$ so the matrices are bigger, but not a lot. But with elements of, say, $2^64$ it can become unmanageable. Thank you. Mar 16, 2021 at 11:19

First, I've attempted to create a similar but quadratic sbox over $$\mathbb{Z}_{17}$$, and found that all non-trivial candidates are linear.

While I'm not yet able to find a proof that all functions meeting the "restricted-commutative" requirement are all isomorphically linear, I think it may have some serious implication on the security of your scheme.

Having no provable way to construct a non-linear 16x16 sbox is currently the greatest thing that undermines the confidence of the security.

Second, I've checked the sbox proposed by the question, and cannot yet find significantly meaningful linear relationship within.

C Program 1: quadratic sbox search program

Compile the following C program and then invoke it with 1 parameter specifying (in decimal) the number of threads to create; linking with libc and libpthread is required.

#include <stdio.h>
#include <stdlib.h>
#include <pthread.h>

#define mod 17
#define ANY 0
#define ALL 1

void long2vec(long s, int *v, int c)
{
for(int i=0; i<c; i++)
{
v[i] = s % mod;
s /= mod;
}
}

int f(int a, int b, int c[])
{
a %= mod;
b %= mod;
return (c[0] + c[1]*a + c[2]*b + c[3]*a*a + c[4]*a*b + c[5]*b*b) % mod;
}

int thrdcnt;
int cmax = mod*mod*mod*mod*mod*mod;

void *t(void *argp)
{
long ind = (long)argp;
int c[6], x[4], u, v;

for(; ind<cmax; ind+=thrdcnt)
{
long s;
long comm = 0, assoc = 0, fail = 0;
s = ind;
long2vec(s, c, 6);

if( (!c[1] && !c[3] && !c[4]) || (!c[2] && !c[4] && !c[5]) )
continue; // Exclude ones ignoring either operand.
if( c[1] == c[2] && !c[3] && !c[5] )
continue; // Exclude ones that are commutative non-trivials.
if( !c[3] && !c[4] && !c[5] )
continue; // Exclude all linear results.

for(s=0; s<mod*mod*mod*mod; s++)
{
long2vec(s, x, 4);

if( f(x[0], x[1], c) == f(x[1], x[0], c) && x[1] != x[0] )
comm++;

if( f(f(x[0], x[1], c), x[2], c) == f(x[0], f(x[1], x[2], c), c) )
assoc++;

u = f(f(x[0], x[1], c), f(x[2], x[3], c), c);
v = f(f(x[0], x[2], c), f(x[1], x[3], c), c);
if( u != v ) { fail++; break; }
}
if( fail == 0 )
dprintf(
1, "%d\t%d\t%d\t%d\t%d\t%d : %d, %d\n",
c[0], c[1], c[2], c[3], c[4], c[5], comm, assoc);
}

return NULL;
}

int main(int argc, char *argv[])
{
thrdcnt = argc >= 2 ? atoi(argv[1]) : 1;
pthread_t *threads = malloc(sizeof(pthread_t) * thrdcnt);
for(int i=0; i<thrdcnt; i++)
{
pthread_create(threads+i, NULL, t, (void *)(intptr_t)i);
}
for(int i=0; i<thrdcnt; i++)
{
pthread_join(threads[i], NULL);
}
return 0;
}


C Program 2: linear relationship search program

Compile and invoke the following C program; linking with libc is required.

#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>

int latinsquare[256] = {
+7,  9, 13, 10, 15,  2,  0,  6,  3, 12,  8,  4,  1,  5, 14, 11,
+1, 15,  6,  3,  9,  4, 11, 13, 10,  5, 14,  2,  7, 12,  8,  0,
+3,  0, 12,  1, 11,  8,  9,  5,  7, 13,  2, 14, 10,  6,  4, 15,
+4,  6, 15,  8, 13,  1,  5,  9, 14, 11, 10,  7,  2,  0,  3, 12,
+0,  3,  8, 15, 10, 12,  7, 14,  9,  2, 13,  5, 11,  4,  6,  1,
10, 11,  5,  7,  0, 14, 15, 12,  1,  6,  4,  8,  3, 13,  2,  9,
+5, 14, 10, 13,  8, 11,  4,  3,  6,  1, 15,  0, 12,  7,  9,  2,
15,  1,  4,  0,  7,  6, 10,  2, 11, 14,  5, 13,  9,  8, 12,  3,
12,  8,  3,  6, 14,  0,  2, 10, 13,  7,  9, 11,  5,  1, 15,  4,
13,  2,  7,  5,  4,  9,  8,  1, 12,  3,  0, 15,  6, 10, 11, 14,
+6,  4,  1, 12,  2, 15, 14,  7,  5, 10, 11,  9, 13,  3,  0,  8,
+9,  7,  2, 11,  1, 13,  3,  4,  0,  8, 12,  6, 15, 14,  5, 10,
11, 10, 14,  9,  3,  5,  1,  8, 15,  4,  6, 12,  0,  2, 13,  7,
14,  5, 11,  2, 12, 10,  6,  0,  4, 15,  1,  3,  8,  9,  7, 13,
+8, 12,  0,  4,  5,  3, 13, 11,  2,  9,  7, 10, 14, 15,  1,  6,
+2, 13,  9, 14,  6,  7, 12, 15,  8,  0,  3,  1,  4, 11, 10,  5,
};

#define P 0x11

int mulb(int a, int b) // multiplication in F_{2^4}
{
int x = 0;
for(int i=0; i<4; i++)
{
x ^= b&1 ? a : 0;
b >>= 1;
a <<= 1;
a ^= a&0x10 ? P : 0;
}
return x;
}

int mulf(int a, int b) // modular multiplication in Z_{16}
{
return (a * b) & 15;
}

int (*mul)(int a, int b) = mulf; // either mulb or mulf

void mksbox(int x[16], uint64_t s)
{
int i;

for(i=0; i<16; i++) x[i] = i;

for(i=16; i>1; i--)
{
int p = s % i;
int j = 16 - i;
int t = x[p+j];
x[p+j] = x[j];
x[j] = t;
s = (s - p) / i;
}
}

float misses(int ab, int map[256])
{
int a = ab >> 4, b = ab & 15;
int x, y;
int cand = 0;
float min = -1;
float num = 0, den = 0;

if( !a || !b || a == b ) return min;

for(x=0; x<256; x++) map[x] = 0;

for(x=0; x<16; x++)
{
for(y=0; y<16; y++)
{
int u = mul(a,x) ^ mul(b,y);
int v = latinsquare[x*16+y];
map[u*16+v]++;
}
}

for(x=0; x<16; x++)
{
for(y=0; y<16; y++)
{
num += map[x*16+y];
if( map[x*16+y] ) den++;
}
}

for(x=0; x<16; x++)
{
for(y=0; y<16; y++)
{
num += map[x*16+y];
if( map[x+16*y] ) den++;
}
}

dprintf(
1, "\na: % 2d, b: % 2d\n",
ab>>4, ab&15);
for(int i=0; i<256; i++)
dprintf(1, "% 3d%c", map[i], i%16==15?'\n':' ');

min = num / den;
return min;
}

int thrdcnt;

void *t(void *argp)
{
long ind = (long)argp;
int rec = -1;
int map[256];
float max = 0;
float subret;

for(; ind<256; ind+=thrdcnt)
{
subret = misses((int)ind, map);
if( subret > max )
{
max = subret;
rec = (int)ind;
}
}

return NULL;
}

int main()
{
thrdcnt = 1;
t(0);
return 0;
}

• I would like OP disclose the method through which he derives the 16x16 sbox so that we can verify its resistance against differential and linear analysis more rapidly. Mar 17, 2021 at 9:13
• The method is bruteforcing all latin squares and checking for the property, cutting off branches with yet undefined cells but failing the property on the defined ones. So bruteforcing at all. I've not been able to find a method based on increasing size of tables or something similar. Hope this answers your question. Mar 17, 2021 at 12:33
• @daniel I've updated the answer with a second part, and good news for you, that I find no significant linear relationship in your sbox using modular arithmetic and binary finite field. Mar 18, 2021 at 4:23

Alice, the signer, chooses a public value $$C$$ and a two random keys $$K$$, $$Q$$. Its credentials are $$C$$, $$(CK)$$ and $$(QK)$$. To sing a value, $$H$$, Alice computes $$S=(HQ)$$.

Bob needs to verify if Alice has signed $$H$$. Computes $$(HQ)(CK)$$ and $$(HC)(QK)$$. Both values must be equal due to restricted commutativity if $$(HQ)$$ is a valid signature from Alice.

The current $$m(t,k)$$ is invertible when knowing $$k$$, so with a verification transcript $$(HC)(QK)$$ it's easy to create a forgery of $$S=(HQ)$$ by inverting $$(HC)(QK)$$ using $$(CK)$$ as $$k$$

The fix is simple, change the public credential $$(QK)$$ to $$(KQ)$$, and verification checks $$(CK)(HQ) = (CH)(KQ)$$.

• But $S=(HQ)$ is already the signature you have. You are forging a signature for H using the signature for H. I don't see very clearly what do you gain. Your scheme is valid though. Mar 18, 2021 at 13:06
• @daniel I'm forging an S for a different H, whereby I create a verification transcript (HC)(QK) using the public C and (QK), then compute the forgery S=(HQ) again with the public (CK) as k to the inverse m function. It's always subtle to notice things like this, do please take your time. Mar 18, 2021 at 15:20
• More work. You're right I compute R=(HC)(QK) for the desired H and reverse R=(HQ)(CK). I'll take your suggestion. I think I'll put a file on contributors. Do you mind? Any idea on this? Mar 19, 2021 at 14:34
• @daniel I don't mind, go ahead. Mar 20, 2021 at 2:33