Deciding EA-equivalence via invariants

We define a family of efficiently computable invariants for (n,m)-functions under EA-equivalence, and observe that, unlike the known invariants such as the differential spectrum, algebraic degree, and extended Walsh spectrum, in the case of quadratic APN functions over F2n\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathbb {F}_{2^n}$\end{document} with n even, these invariants take on many different values for functions belonging to distinct equivalence classes. We show how the values of these invariants can be used constructively to implement a test for EA-equivalence of functions from F2n\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathbb {F}_{2}^{n}$\end{document} to F2m\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\mathbb {F}_{2}^{m}$\end{document}; to the best of our knowledge, this is the first algorithm for deciding EA-equivalence without resorting to testing the equivalence of associated linear codes.


Introduction
Let F 2 n denote the finite field with 2 n elements for some positive integer n, and let F * 2 n denote its multiplicative group. The vector space of dimension n over F 2 will be denoted by F n 2 . An (n, m)-function, or vectorial Boolean function, is any mapping from F n 2 to F m 2 or, equivalently, from F 2 n to F 2 m . We typically assume that m = n, i.e. we concentrate on functions from a finite field of characteristic two to itself; however, the approach given in the present paper can be applied to an arbitrary pair of dimensions (n, m).
In the particular case of n = m, we can conveniently represent (n, n)-functions by univariate polynomials over F 2 n ; more precisely, any (n, n)-function F can be written as F (x) = 2 n −1 i=0 c i x i for some coefficients c i ∈ F 2 n . This form, called the univariate representation of F , always exists, and is unique. In the general case, any (n, m)-function F can be represented uniquely as a multivariate polynomial of the form F (x) = u∈F n 2 a u n i=1 x u i i for a u ∈ F m 2 , where v i denotes the i-th component of the vector v = (v 1 , v 2 , . . . , v n ) ∈ F n 2 . We will not concern ourselves too deeply with the choice of representation here, as it is typically only important when defining and classifying (n, m)functions. However, when illustrating some of the concepts of the EA-equivalence test with examples, we will mostly use (n, n)-functions and the univariate representation.
As suggested above, the finite field F 2 n can be identified with the vector space F n 2 , and so the elements of F 2 n can be identified with binary vectors of n bits; thus, (n, m)-functions can be understood as transformations that take an n-bit sequence as input, and produce an m-bit sequence as output. Thanks to this interpretation, vectorial Boolean functions naturally find applications in the theory and practice of various areas of mathematics and computer science. In particular, vectorial Boolean functions are used in modern block ciphers in the role of so-called "S-boxes", or "substitution boxes", and typically constitute the only non-linear part of the cipher. Consequently, the security of the cipher directly depends on the properties of the underlying S-boxes, which motivates the study of vectorial Boolean functions with respect to their cryptographic properties. It is also intuitively clear that (n, n)-functions constitute one of the most practically significant cases, since in cryptography one typically wants to replace a bit sequence with a different bit sequence of the same length.
A basic property of any (n, m)-function, which also has cryptographic implications, is its algebraic degree. Given an (n, m)-function with ANF F (x) = u∈F n 2 u n i=1 x u i i , the algebraic degree of F is defined as the largest degree of any term x u i i that has a non-zero coefficient a u . In other words, the algebraic degree of F is the multivariate degree of its ANF. Thus, for instance, F (x) = x 1 x 2 x 4 + x 2 x 3 would have algebraic degree 3. In the case of an (n, n)-function given by its univariate representation, the algebraic degree also has a natural interpretation. Given a positive integer i, the binary weight of i is the number of non-zero entries in its binary representation; for example, 19 is written as 10011 in binary, and hence has binary weight 3. The algebraic degree of F (x) = 2 n −1 i=0 c i x i is the largest binary weight of any exponent i with c i = 0. Functions of algebraic degree 1, 2, and 3, are called affine, quadratic, and cubic, respectively. An affine function A with A(0) = 0 is called linear. An affine (n, n)-function A satisfies A(x) + A(y) + A(z) = A(x + y + z) for any x, y, z ∈ F 2 n ; similarly, a linear (n, n)-function L satisfies L(x) + L(y) = L(x + y) for any x, y ∈ F 2 n . It is desirable for vectorial Boolean functions used as S-boxes to have a high algebraic degree, since the latter indicates a good resistance to higher-order differential attacks [10,14].
Two of the most important cryptographic properties of (n, m)-functions are the differential uniformity and the nonlinearity. Suppose that F is an (n, m)-function for some positive integers n and m. Let δ F (a, b) denote the number of solutions x ∈ F n 2 to the equation } is called the differential spectrum of F . The largest value in the differential spectrum is denoted by δ F and is called the differential uniformity of F .
The existence of a large number of solutions x to some equation of the form F (x)+F (a+ x) = b for some a = 0 and b makes the function F vulnerable to differential cryptanalysis [2]. The value of δ F should thus be as low as possible. Since x + a is a solution to the aforementioned equation whenever x is, the minimum possible value of δ F is 2; the class of (n, n)functions attaining this optimal value is called the class of almost perfect nonlinear (APN), and has been an object of intense study since its introduction by Nyberg in the 90's [16].
The nonlinearity N L(F ) of F is simply the minimum Hamming distance between any component function of F and any affine (n, 1)-function, the component functions of F being the (n, 1)-functions F c of the form F c (x) = Tr m (cF (x)), where Tr m : F 2 m → F 2 is the absolute trace defined by Tr m (x) = m−1 i=0 x 2 i for any x ∈ F 2 n ; when the dimension m is clear from context, we will sometimes write just Tr instead of Tr m . The nonlinearity should be high in order to resist linear cryptanalysis [15].
When studying "linear properties" of functions, such as their nonlinearity, it is useful to adapt some linear-algebraic notions from the vector space F n 2 even in the case when we are working with the finite field F 2 n . The linear span of a set S ⊆ F 2 n is simply the set of all possible linear combinations of the elements in S, i.e. if S = {s 1 , s 2 , . . . , s k } for some positive integer k and for s i ∈ F 2 n , then Span(S) = {c 1 s 1 + c 2 s 2 + · · · + c k s k : c 1 , c 2 , . . . , c k ∈ F 2 }; obviously, the span can be defined in the same way in the case of the vector space F n 2 . Having formalized the linear span, it is straightforward to carry over other notions from F n 2 , such as that of linear independence, and that of a basis of F 2 n (being a linearly independent set B ⊆ F 2 n with Span(B) = F 2 n ).
A useful tool for analyzing vectorial Boolean functions is the Walsh transform, which is an integer valued function W F : F n 2 × F m 2 → Z associated with F : F n 2 → F m 2 , and given by the prescription In the case of (n, n)-functions, we can more succinctly write Due to the huge number of (n, m)-functions even for small values of n and m, their classification is typically only performed up to some equivalence relation that preserves the properties being studied. In the case of cryptographically optimal vectorial Boolean functions, the most general equivalence relation preserving both the differential uniformity and the nonlinearity is the so-called Carlet-Charpin-Zinoviev-equivalence, or CCZ-equivalence [9]. Two (n, m)-functions F and G are said to be CCZ-equivalent if there is an affine per- Testing whether two given functions F and G are CCZ-equivalent is usually done by means of the equivalence of linear codes [5,13]. More precisely, a particular linear code C F is associated with F , and a particular linear code C G is associated with G; F and G are then CCZ-equivalent if and only if C F and C G are equivalent. Testing whether two given linear codes are equivalent has the advantage that it is usually already implemented in most mathematical software, such as the Magma programming language that we use for most of our computations [4].
Unfortunately, computationally testing CCZ-equivalence in this way can reliably be performed only when the dimensions m and n are relatively small; in the case of (n, n)functions, this means n ≤ 9, since for higher values of n, the memory consumption becomes overwhelming, and the test cannot be performed in a lot of cases. Furthermore, the current implementation of Magma can give false negatives due to insufficient memory; in other words, if the equivalence test outputs "false", we have no reliable way of determining whether this is due to insufficient memory, or due to a successfully completed exhaustive search proving the inequivalence of the linear codes (and hence, vectorial Boolean functions) in question.
Ruling out e.g. CCZ-equivalence can be facilitated by means of invariants, i.e. properties or statistics that are preserved under CCZ-equivalence. The differential spectrum and the extended Walsh spectrum are invariant under CCZ-equivalence, i.e. if two (n, m)-functions F and G are CCZ-equivalent, then their differential spectra and extended Walsh spectra are the same. Unfortunately, the Walsh spectra and differential spectra of all known APN functions are the same (with some rare exceptions in the case of the Walsh spectrum), rendering these invariants nearly useless in practice. Other invariants, such as the -rank and -rank have been introduced [12], that can take different values for distinct CCZ-classes of functions, and can therefore be used to rule out CCZ-equivalence in some cases. The major drawback of these invariants is that they require significant computational resources, meaning that, on the one hand, their calculation takes a long time, e.g. around ten days for a single -rank over F 2 10 , and, on the other hand, computing these invariants for F 2 n with n > 10 is impossible at the moment due to overwhelming memory requirements.
A special case of CCZ-equivalence is extended affine equivalence, or EA-equivalence. Two (n, m)-functions F and G are said to be EA-equivalent if there exist affine functions with A 1 and A 2 being bijective. We will refer to A 1 as the outer permutation and to A 2 as the inner permutation throughout the paper. CCZ-equivalence is strictly more general than EAequivalence combined with taking inverses [7], but in certain cases, such as for quadratic and monomial functions, checking whether two functions (or, potentially, their inverses) are EA-equivalent is enough to decided CCZ-equivalence: two quadratic APN functions are CCZ-equivalent if and only if they are EA-equivalent [18]; and two power functions are CCZ-equivalent if and only if they are cyclotomic equivalent [11]. Recall that two power functions F (x) = x d and G(x) = x e over F 2 n are said to be cyclotomic equivalent if e ≡ 2 k d mod (2 n − 1) or e −1 ≡ 2 k d mod (2 n − 1); furthermore, cyclotomic equivalence is a special case of EA-equivalence and taking inverses. This is particularly interesting when one takes into account that all known APN functions (which fall into more than 20000 distinct CCZ-equivalence classes) are CCZ-equivalent to a monomial or quadratic function, with only a single exception in dimension n = 6. Thus, from a practical point of view, being able to test functions for EA-equivalence is virtually as useful as being able to test them for CCZ-equivalence, at least as far as the classification of APN functions is concerned. Surprisingly, despite its simple definition, the only known algorithm to date for computationally testing the EA-equivalence of two given functions is one by means of associated linear codes, much like in the case of CCZ-equivalence [13]; the associated codes used in this approach are of a somewhat more complicated form than the ones used in the CCZ-equivalence test, and so this approach is even more restrictive with respect to computational resources and memory requirements. Indeed, testing EA-equivalence for quadratic functions (which coincides with CCZ-equivalence) is typically done by testing them for CCZ-equivalence. Algorithms for testing the EA-equivalence of two functions in some other special cases (grouped under the umbrella term "restricted EA-equivalence") have previously been studied in [3], [8], and [17].
Since EA-equivalence is a special case of CCZ-equivalence, any CCZ-invariant is also an EA-invariant. As mentioned above, the extended Walsh spectrum and differential spectrum are practically useless in the case of APN functions, as they almost always take the same value in the APN case, while the -and -rank involve somewhat laborious computations. EA-equivalence being less general than CCZ-equivalence, it is natural to expect to have properties that are EA-invariant but not CCZ-invariant. One such property is the algebraic degree, which is preserved by EA-equivalence, but not by CCZ-equivalence. Unfortunately, this is not terribly useful for classifying APN function either since, as mentioned above, nearly all known instances of APN functions are quadratic.
In this paper, we present an approach for computationally testing the EA-equivalence of two (n, m)-functions by first guessing the outer permutation A 1 , applying its inverse to (1) to obtain a relation of the form F • A 2 + A = G , and then solving the latter for A 2 and A . In the case of (n, n)-functions with n even, our approach allows the set of possible affine permutations A 1 to be drastically reduced (as opposed to exhaustive search), which makes the entire procedure computationally feasible. Our approach has the advantage that it can be broken down into a multitude of small independent steps, which makes the resulting algorithm easily parallelizable. Unlike the CCZ-equivalence test and EA-equivalence test described in [13], which rely on testing the equivalence of a pair of linear codes (and therefore require specialized and rather complex algorithms), our approach uses only basic arithmetics and linear algebra, and can be easily implemented in any general-purpose programming language, and ran on any computer. Furthermore, each of the individual steps comprising the algorithm has a concrete and meaningful input and output that can be monitored and verified. This precludes the possibility of false positives or negatives as in the case of the current CCZ-equivalence test.

A family of EA-invariants
Let m, n, k be positive integers, and t be an element of F n 2 . We denote by T k (t) the set of all k-tuples of elements from F n 2 that add up to t, i.e.
If A is an affine (n, n)-permutation, then the image of any k-tuple ( , the sum of whose elements is For any (n, m)-function F , let F k (t) denote the multiset of all sums of the form F ( The multiplicities of F k (0) are then an EA-invariant for any even value of k.
In particular, the multiplicities of F k (0) and G k (0) (that is, the number of times that each element occurs in each multiset) are the same when k is even.
for odd values of k, and A 2 (x 1 + x 2 + · · · + x k ) + A 2 (0) for even values of k. Since x 1 + x 2 + · · · + x k = 0 by assumption, the sum A 2 (x 1 ) + · · · + A 2 (x k ) is then A 2 (0) for odd k, and 0 for even k. Thus, computing all sums of the form   Recall that the majority of known APN functions are quadratic, and that testing the equivalence of quadratic (n, n)-functions represents the case of highest practical interest. One very useful observation that we can make in the quadratic case is that we can assume A 1 (0) = A 2 (0) = 0, which (as we see later), greatly simplifies the complexity of the entire EA-equivalence test; and, in particular, means that the multiplicities of F k (0) are an EA-invariant for quadratic functions in the case of odd values of k as well.
Proposition 2 Let F, G be quadratic (n, n)-functions for some positive integer n, and sup- Proof We can assume that F is purely quadratic, i.e. of the form is an affine function. Then the composition As suggested above, Proposition 2 implies that the multiplicities of F k (0) are an invariant for both odd and even values of k in the quadratic case. Note that the condition of the function being quadratic is necessary, as witnessed by e.g. F (x) = x 15 and 15 over F 2 six, where α is a primitive element of F 2 six: the elements of the finite field in question fall into three distinct classes based on their multiplicities in F 3 (0), but into five distinct classes based on their multiplicities in G 3 (0). Corollary 3 Following the notation and hypothesis of Proposition 1, if F and G are in addition quadratic, then the multiplicities of F k (0) and G k (0) are the same for any value of k.
The complexity of computing the multiplicities of F k (t) for an (n, m)-function F increases exponentially with each increment of k. Fortunately, computing the multiplicities via the Walsh transform of F results in a complexity that does not depend on the value of k.

Proposition 4
Let F be an (n, m)-function, k be a positive integer, t ∈ F n 2 and s ∈ F m 2 . Let M F k (t, s) denote the number of k-tuples (x 1 , x 2 , . . . , x k ) such that x 1 + x 2 + · · · + x k = t and F (x 1 ) + F (x 2 ) + · · · + F (x k ) = s. Then Proof From the definition of the Walsh transform, the expression By changing the order of summation, this becomes The statement then follows by recalling that a∈F n 2 (−1) Tr n (ax) evaluates to 0 for any 0 = x ∈ F n 2 , and evaluates to 2 n for x = 0.
Finding the multiplicity of a given element s ∈ F m 2 in F k (t) now amounts to computing the Walsh coefficients W F (a, b) of F , raising them to the power k, and combining them according to (2). We note that for the purposes of testing EA-equivalence, we always assume t = 0, and hence (2) F (a, b). Furthermore, the Walsh coefficients W F (a, b) can be precomputed for all F from a known set of EA-representatives, allowing the computations to be sped up at the cost of storing the precomputed result.

Remark 5
We note that, in the case of APN functions, the multiset of the multiplicities of F 3 (0) is essentially the same as the multiset 0 F studied in [6]. The latter, given an (n, n)-function F , is defined as the multiset The equation F (x) + F (a + x) + F (a) = b has either 0 or 2 solutions for any 0 = a ∈ F 2 n and any b ∈ F 2 n if F is APN. Thus, an equivalent invariant would be the multiset and it is easy to see that this can equivalently be rewritten as which is essentially the same as F 3 (0). As pointed out in [6], the multiset 0 F is a CCZinvariant for quadratic APN functions.

Guessing the outer permutation
Suppose that we are given two EA-equivalent functions F and G from F n 2 to F m 2 for some positive integers n, m, so that In particular, we can always assume that c 1 = 0 without loss of generality, so that A 1 is linear. In the following, we will write simply G = L 1 • F • A 2 + A.
By Proposition 1, for even values of k, we have Besides justifying that the multiset of multiplicities of F k (0) is an EA-invariant for any positive even integer k, the above relation gives us some information about L 1 ; namely, it implies that if L 1 (x) = y for some x, y ∈ F m 2 , then the multiplicity of x in F k (0) should be the same as that of y in G k (0). If the elements of F m 2 are partitioned according to the multiplicities of F as for some positive integer s, so that all elements in K i for every 1 ≤ i ≤ s have the same multiplicity in F k (0), and elements in K i and K j have distinct multiplicities for i = j ; and, similarly, as according to the multiplicities of G, then we must have for all 1 ≤ i ≤ s. We will say that any permutation L 1 satisfying L 1 (K i ) = C i for all i respects the two partitions of F m 2 . Consequently, we obtain conditions that can be used to restrict the possible choices for L 1 . Intuitively, the larger the number of classes s in the partition of F m 2 , the fewer linear permutations L 1 can satisfy the conditions thus obtained. In particular, if all elements of F m 2 occur with the same multiplicity, we do not obtain any information on L 1 . This is clearly the case when F is a permutation. Furthermore, the same appears to be true for all APN (n, n)-functions with odd n (regardless of whether they are permutations or not), which is why we concentrate on fields of even extension degree in our work.
All linear permutations L respecting the partitions F m 2 = K 1 ⊕ · · · ⊕ K s and F m 2 = C 1 ⊕ · · · ⊕ C s can now be found by trying to guess the values of L on a basis of F m 2 , and backtracking whenever some assignment violates these partitions. An algorithmic description of this procedure is provided below under Algorithm 1 and Algorithm 2. The former presents the general framework for partitioning F m 2 according to the multiplicities of sums of values of F and G, while the latter describes the process of reconstructing all linear permutations that respect the constructed partitions. We remark that the algorithm is described for the particular case of k = 4 (which is what we have mostly used in practice for our computational experiments), but the principle trivially generalizes to any value of k. We also note that computing the number M F k (0, s) of k-tuples whose values under F add up to a given s ∈ F m 2 can be done via the values of the Walsh transform as described in Proposition 4; this is particularly useful if the selected value of k is large, or if a precomputed table of the Walsh coefficients for one (or both) of the tested functions is available.
Let us take a closer look at Algorithm 1. We first fix an even value of k, for instance k = 4. Given two (n, n)-functions, F and G, that we would like to test for equivalence, we begin by computing the multiplicities of the elements in the multisets F k (0) and G k (0). The number of times that the element s ∈ F m 2 appears in F k (0) is denoted by M F k (0, s) (this means that M F k (0, s) k-tuples (x 1 , x 2 , . . . x k ) with x 1 + x 2 + · · · + x k = 0 satisfy F (x 1 ) + F (x 2 ) + · · · + F (x k ) = s).
Using these multiplicities, we partition F m 2 in two ways: using the multiplicities {M F k (0, s) : s ∈ F m 2 }, and using the multiplicities {M G k (0, s) : s ∈ F m 2 }. More precisely, we write F m 2 as F m 2 = K 1 ⊕ K 2 ⊕ · · · ⊕ K s , with K 1 , K 2 , . . . , K s being disjoint sets of elements; two elements s 1 and s 2 are in the same block K i of the partition if and only if M F k (0, s 1 ) = M F k (0, s 2 ), i.e. if s 1 and s 2 occur with the same multiplicity in F k (0). Equiv-alently, we could say that the multiplicities M F k (0, s) induce an equivalence relation, in which two elements s 1 , s 2 ∈ F m 2 are equivalent precisely when M F k (0, s 1 ) = M F k (0, s 2 ); the blocks K 1 , K 2 , . . . , K s are then the equivalence classes of this equivalence relation. In the same way that K 1 , K 2 , . . . , K s is the partition induced by F k (0), C 1 , C 2 , . . . , C s is the partition induced by G k (0). If F and G are EA-equivalent, then the number of blocks in both partitions must be the same, and the individual blocks must have the same sizes. Thus, if s = s , or if the multiset {#K i : i = 1, 2, . . . , s} is not equal to {#C i : i = 1, 2, . . . , s}, we can immediately conclude that F and G are not EA-equivalent. Otherwise, we can rearrange the blocks K 1 , K 2 , . . . , K s and C 1 , C 2 , . . . , C s in such a way that #K i = #C i for i = 1, 2, . . . , s. At this point, we know that if F and G are equivalent via L 1 • F • A 2 + A = G, then L 1 must map K i to C i for i = 1, 2, . . . s. This additional information allows us to significantly reduce the number of linear permutations L 1 that needs to be considered.
The set of all linear permutations preserving the partitions can be found using Algorithm 2. The latter is essentially an exhaustive search that tries to guess the values of L 1 on a basis For any such element x, we can find the indices j, j such that x ∈ K j and L 1 (x) ∈ C j . If j = j , then L 1 does not respect the partitions, and so we backtrack, attempting a different guess for b i . If we do not find any contradiction of this type, we proceed to guessing the value of b i+1 . We continue in this manner until we have exhausted all possibilities.
The partitions F m 2 = K 1 ⊕K 2 ⊕· · ·⊕K s can be precomputed for representatives from e.g. all known EA-classes of APN functions; in particular, we refer to our computational results described in Section 5 where we describe how we provide such pre-computed results for all currently known APN functions over F 2 n up to dimension n = 10. When using Algorithms 1 and 2 to find all possibilities for the outer permutation L 1 in L 1 •F •A 2 +A = G, however, we need to know the partitions according to both F and G, which makes the precomputation of the permutations L 1 impossible.
Nonetheless, we can observe that the set of linear permutations L 1 mapping K i to C i for every 1 ≤ i ≤ s is simply a coset in the symmetric group of F m 2 of the subgroup of linear permutations mapping K i to K i for 1 ≤ i ≤ s. The latter can be precomputed for known EA-representatives, and hence finding a single linear permutation mapping every K i to C i with 1 ≤ i ≤ s allows us to reconstruct all such permutations by composing it with the precomputed ones. This can be formalized as follows.

Proposition 6
Let n be a positive integer, and F n 2 = K 1 ⊕ K 2 ⊕ · · · ⊕ K s and F n 2 = C 1 ⊕ C 2 ⊕ · · · ⊕ C s be two partitions of the elements of F n 2 such that #K i = #C i for every 1 ≤ i ≤ s. Let K be the set of all linear permutations L of F n 2 such that L(K i ) = K i for all 1 ≤ i ≤ s, and let P be the set of all linear permutations L of F n 2 such that L(K i ) = C i for 1 ≤ i ≤ s. Then K is a subgroup of the symmetric group of F n 2 , and P is a coset of K.
Proof The composition of two linear permutations is clearly a linear permutation itself, and so is the inverse of a linear permutation. Furthermore, if L 1 and L 2 are linear permutations that permute some set K i ⊆ F n 2 , then their composition and their inverses do so as well. Thus, K is closed under composition and taking inverses, and is a subgroup of the symmetric group of F n 2 . Now, suppose that L is a linear permutation of F n 2 mapping some subset K i ⊆ F n 2 onto some C i ⊆ F n 2 . Then K • L is also a linear permutation mapping K i onto C i for any K ∈ K. Thus, K → K • L maps K to P, and is clearly invertible since L is a permutation. Consequently, P is a coset of K represented by L.
Besides delegating a large portion of the work in constructing P to the precomputation of K, Proposition 6 allows us to estimate the complexity of testing EA-equivalence between a function F (which we can assume is a known EA-representative) and another function G inducing a partition of F m 2 compatible with the one induced by F . Furthermore, it is clear that the size of the group K of linear permutations that preserve the partition F m 2 = K 1 ⊕ · · · ⊕ K s induced by the multiplicities in k F (0) is an EA-invariant. What makes this interesting, is that it is more discriminating than the sizes of the partition classes: for instance, the APN functions F (x) = x 3 and G(x) = x 3 + α 11 x 6 + αx 9 over F 2 six (where α is primitive in F 2 six) both partition F 2 six into three classes of size 1, 21, and 42, respectively; but the group of linear permutations preserving the partition of F (x) contains 1008 elements, while the group of linear permutations preserving the partition of G has 336 elements. Thus, precomputing the groups K of linear permutations preserving the partition for representatives from the known classes of APN functions has the additional advantage that it allows us to rule out equivalence in more cases (using a stronger invariant). We note that the actual elements of the group K are not, in general, invariant under EAequivalence.
To give some basic idea of how efficient these processes are, we have computed the groups K for representatives from all switching classes of APN functions over F 2 n with n ∈ {6, 8} [12]. The results are presented in Table 1 below. The first column gives the dimension n of F 2 n . The functions are indexed in the second column in the same way as in [12]. The next two columns give the time in seconds for computing the partition of F 2 n according to the quadruple sums of F directly and using the Walsh transform, respectively (including the time in seconds for precomputing the Walsh coefficients). The following column gives the time for computing all linear permutations preserving the corresponding partition. The last column gives the number of linear permutations found in each case, which is a direct measure of the complexity of an EA-equivalence test by our method, as the approach for guessing the inner permutation (described in the following Section 4) has to be applied to every possible choice of the outer permutation.
We note that the running times are highly dependent on the programming language, implementation, and computational equipment used, and the ones presented in the paper are given only for illustrative purposes.
For comparison, there are 27998208 linear permutations of F 2 6 , and 132640470466560 linear permutations of F 2 8 .

Guessing the inner permutation
If, in addition to the (n, m)-functions F and G, we know the linear permutation L 1 in the relation L 1 • F • A 2 + A = G, we can apply its inverse, L −1 1 to both sides, obtaining where A = L −1 1 • A and G = L −1 1 • G. A pair of affine (n, m)-functions A 2 , A satisfying the above relation then exists if and only if F is EA-equivalent to G.
Once again, let us write c = A (0) and c 2 = A 2 (0), and L 2 = A 2 + c 2 and A = L + c for the linear parts of A 2 and A, respectively. Substituting 0 for x in (3) yields F (c 2 ) + c = G (0). Since we know G , and hence also G (0), this means that any choice of c 2 uniquely determines c. It is thus enough to loop over all possible choices of c 2 ∈ F n 2 and take c = F (c 2 ) + G (0) in order to exhaust all possibilities for (c 2 , c ). As observed in Proposition 2, if F and G are quadratic, then we can assume that c 2 = 0 and c = G (0) without loss of generality; for functions of higher algebraic degree, we have to consider all possible values of c 2 . In the following, we assume that we have guessed the constants c 2 and c, and rewrite (3) as where G (x) = G (x +c 2 )+c. It now remains to look for a pair of functions L 2 : F m 2 → F m 2 and L : F n 2 → F m 2 satisfying (4). To guess the permutation L 2 , we observe that, given some k-tuple (x 1 , x 2 , · · · , x k ) ∈ T k (0), by Proposition 1, we have G (x 1 ) + · · · + G (x k ) = F (L 2 (x 1 )) + · · · + F (L 2 (x k )), and thus, if some element x i ∈ F n 2 is part of a k-tuple whose sum under G is t, then its image L 2 (x i ) under L 2 must be part of some k-tuple whose sum under F is t. We state this formally as follows.

Proposition 7
Let F and G be (n, m)-functions for some positive integers n, m such that F • L 2 + L = G for L 2 , L linear and L 2 bijective. Let k be a positive integer, and, for any Then, if (x 1 , x 2 , · · · x k ) ∈ T k (0) with G(x 1 ) + · · · + G(x k ) = t, we must have Using (5) for k = 3, we can significantly reduce the domains of L 2 (x) for x ∈ F m 2 , i.e. the ranges of possible values that L 2 (x) can take. A large number of the domains end up consisting of three elements (although we do obtain larger domains in some cases). Since guessing L 2 amounts to guessing its values on a basis of F m 2 , the elements of the basis can be chosen in such a way that the Cartesian product of the respective domains is small. In most cases, we can indeed choose the basis elements in such a way that all domains consist of three elements, and thus end up with only 3 n possibilities for L 2 .
In addition, assuming that e.g. F is a known representative, some precomputations are possible; namely, the sets O F k (0, t) can be precomputed for k = 3 and all values of t. Alternatively, the roles of F and G in (4) can be swapped by composing both sides with the inverse of L 2 from the right, in which case the sets F k (0, x) can be precomputed for k = 3 and for all x ∈ F m 2 . We note that for values of k greater than 3, it seems to always be possible to express all elements in F 2 n as the sum of four values of F for an APN function F , in which case Algorithm 3 describes the approach for reconstructing L 2 from (3) suggested by the above considerations. The first part of the algorithm computes the domains D(x) for all elements x ∈ F n 2 ; we then know that for any x ∈ F n 2 we must have L 2 (x) ∈ D(x). All domains are initially set to F n 2 , i.e. no restrictions on the value of L 2 (x) is made. We then compute the sets O F For any element y 1 belonging to a triple (y 1 , y 2 , y 1 +y 2 ) with G(y 1 )+G(y 2 )+G(y 1 +y 2 ) = t, we know that L 2 (y 1 ) must belong to O F 3 (t); we use this to reduce the domain D(y 1 ) of y 1 .
Having computed the domains, the second part of the algorithm consists of finding a basis B = {b 1 , b 2 , . . . , b n } of F n 2 , and constructing all linear permutations L 2 for which L 2 (b i ) ∈ D(b i ) for i = 1, 2, . . . , n. Since we assume that F • A 2 + A = G (with A 2 = L 2 +c 2 , where we also guess the value of c 2 by going through all possibilities), if the choice of L 2 is correct, then A = F • A 2 + G must be affine. For every possible of choice of L 2 and c 2 , we thus compute A and check whether its algebraic degree is at most 1; if so, then we have found the equivalence between F and G.
Recall that by Proposition 2 we can assume that c 2 = 0 if the functions being tested for equivalence are quadratic, which significantly reduces the computation time.
We note that Algorithm 3 will return all affine permutations A 2 for which A = F •A 2 +G is affine. If our goal is to check whether such a permutation exists (which is all that we need for the purposes of the EA-equivalence test), we can immediately terminate as soon as a single such permutation is found. Furthermore, we remark that if (3) is obtained by applying the inverse L −1 1 of the outer permutation, and a solution (A 2 , A) of (3) is found, then this already witnesses that F and G are EA-equivalent.
In order to get an idea of the efficiency of this method, we once again run a number of experiments on representatives from the known APN functions for n = 6 and n = 8. For every pair (F, G) of representatives from the switching classes in [12], we generate a random affine permutation A 2 and a random affine function A, and use Algorithm 3 to attempt to reconstruct A 2 and A from F and G. In the cases when F and G are not EAequivalent this, of course, will fail; in the remaining cases (when F and G do belong to the same EA-equivalence class), we stop as soon as we find the first pair of affine functions (A 2 , A) solving F •A 2 +A = G. For each combination of F and G, we generate 10 pairs of (A 2 , A). Table 2 gives the average running time for solving F • A 2 + A = G for dimensions n = 8. There are 23 switching APN representatives in F 2 8 , and we index them from 1 to 23 in Table 2 in the same order that they are listed in [12]. In the case of n = 6, the running time does not exceed 0.2 seconds in the worst case; we omit a detailed table of the running times for the sake of brevity. The running times are given in seconds, multiplied by a factor of 100; e.g. deciding that F • A 2 + A = G is unsolvable when F is 1.1 and G is 1.2 from [12] takes 7.51 seconds.

Computational results
A recent paper [1] introduces 12 923 new APN functions over F 2 8 , in addition to the more than 8000 instances previously found and documented in [19]. For the purposes of measuring how efficient the multiplicities of the elements in F k (0) are as an invariant, and for speeding-up potential EA-equivalence tests, we have computed the exact partitions for k = 4 for all of these functions. We also perform similar computations for all known APN functions up to dimension n = 10. A complete list of these partitions is available online at https://boolean.h.uib.no/mediawiki. Here, we give a summary of the computed data.
In total, we have computed the partition induced by F k (0) for 21105 CCZ-inequivalent functions F . From these, we have obtained 19300 distinct partitions. Of these, the "Goldlike" partition (which splits the field into three partition classes, of size 1, 70, and 185, respectively) is the most frequently occurring, and is induced by 21 functions including, of course, the Gold function x 3 . The number of partitions that occur only once is 18103; and the remaining partitions occur between two and eleven times.
Most of the partitions contain a large number of classes: indeed, only the "Gold-like" partition described above has three classes, while all other observed partitions have at least 6 classes; the vast majority of functions induce a partition having between 12 and 16 classes, while the largest number of classes, 22, is achieved by only two functions. We recall that a large number of classes intuitively corresponds to a small number of linear permutations respecting the corresponding partition, and consequently to a faster test for EA-equivalence.
In the case of n = 10, we only observe the "Gold-like" partition for all the ten known representatives from the infinite families. However, among the five new functions given in the dataset accompanying [1], we find three that have different (and pairwise distinct) spectra.
For odd dimensions (n = 7 and n = 9), we also compute the partitions induced by the known APN representatives, but these always yield a trivial partition of F 2 n into a zero and non-zero elements (even when we take into account the newly discovered APN classes from [1]).

Conclusion
We have introduced a family of invariants under EA-equivalence, and have shown how their values can be efficiently computed using the Walsh transform. We have experimentally observed that over F 2 n with even n, these invariants can be used to partition quadratic APN functions into small subclasses, thereby significantly facilitating their classification up to EA-and CCZ-equivalence. We have demonstrated how the values of these invariants can be used to restrict the values of the outer permutation A 1 in the relation A 1 • F • A 2 + A = G for two given (n, m)-functions F and G, and have ran experiments in order to measure how much this approach reduces the search space. We have described how a variation of the same invariants can be used to restrict the values of the images of F 2 n under the inner permutation, A 2 , and have combined the above into a computational test for deciding the EA-equivalence of any two (n, m)-functions F and G. Although slower than the standard test for CCZ-equivalence via the permutation equivalence of linear codes, our approach has the advantage that it is easily implementable on any programming language, and can be separated into a multitude of small, independent steps with concrete output, the majority of which can be naturally parallelized and run in different processes or on different computers. Furthermore, this is, to the best of our knowledge, the first efficient algorithm for directly testing the EA-equivalence of two given functions.
One direction for future work would be to investigate the invariants described in Section 2 more closely, and see whether they can be modified in order to provide more efficient restrictions. In the same vein, it would be interesting to investigate the functions for which our experimental results show a large number of choices for the outer permutation A 1 following the restriction described in Section 3, and to see whether some of these choices can be ruled out using some other criterion; this would directly impact the efficiency of the entire EA-equivalence test for these functions.
So far, we have implemented the algorithms described in Sections 3 and 4 in the Magma programming language [4] due to the ease of implementation. As pointed out above, our approach is quite simple, and does not depend on anything more complicated than computing linear combinations of binary vectors, and so it should be readily implementable in any general-purpose programming language. We expect that a careful implementation in an efficient language would further reduce the computational time needed for testing EA-equivalence, and make the method even more useful in practice.

Funding Open access funding provided by University of Bergen (incl Haukeland University Hospital).
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.