1.1. General setting

When we discuss integrable equations, continuous or discrete, there is always the question of which definition of integrability should be used. Unfortunately at this time it is still not possible to give a universal definition of integrability. This is because although integrability can be discussed for various different classes of equations, it manifests itself in different form depending on the class. In place of an elusive universal definition, we do have many different specific properties that are associated with equations that are considered integrable. These properties depend on the class of equations.^a However, there are some properties that can be applied to a variety of classes of equations, such as behavior around singularities, existence of a Lax pair with spectral parameter, or the existence of sufficient number of conservation laws, and these have sometimes been elevated as definitions of integrability, while some other proposed definitions can only be taken as indicators, or necessary conditions.

In this paper we consider difference equations defined on an elementary square of the Cartesian lattice, with the dynamical variables located on the corners of the square, so called quadrilateral equations. The following basic assumptions are made:

1.2. Detailed formulation

For consistency analysis we consider a 3D cube and assign equations on each side of the cube. Furthermore we will extend the consideration from one cube to the full 𝕑³ lattice and embed the equation as follows:

1. 1.

   On a given 𝕑² plane all elementary squares have equations of the same form (that is, only the corner variables change corresponding to the location). The coefficients of the equation may depend on the two lattice parameters associated with the plane but not on the location.

2. 2.

   When extended to the 𝕑³ lattice, the quadrilaterals on parallel planes all carry the same equation but intersecting planes may have different equations.

In 𝕑² the corner variables are usually notated as x_n,m where (n, m) gives the locations on the Cartesian plane, while in 𝕑³ the variables are indexed as x_n,m,k. When dealing with a specific quadrilateral or cube we use shorthand notation

An acceptable quadratic quadrilateral equation can then be written as

where the c_j may depend on the lattice variables p,q, which are associated with the n, m or tilde and hat directions, respectively. The lattice variable for k or bar direction is r. For the consistency cube we need Q12(x,x˜,x^,x˜^,p,q)=0, Q23(x,x^,x¯,x^¯,q,r)=0, and Q31(x,x¯,x˜,x˜¯;r,p)=0, that is respectively. We also need their shifts. When writing this set of equations we have used cyclic convention p → q → r → p and ∼ → ^ → − → ∼, see also Figure 1. Note that the names of the lattice variables in c_i are also used for function names. That means, for example, that there is no c₁(x, y) from which c₁(p, q) and c₁(q, r) could be obtained, rather c₁(p, q) and c₁(q, r) are independent functions.

Fig. 1:

The consistency cube. The picture on the right gives the location of the variables, lattice directions and parameters. The picture on the left explains the coefficients: the variables of a product are connected by a line, adjacent to which we give the name of the corresponding coefficient function. The six functions within a quadrilateral depend on the lattice parameters given at the center of the quadrilateral.

Under the assumptions above the consistency conditions arise as follows: The equations on the faces of the consistency cube are given by:

We need precisely four initial values to compute the other values, three of the initial values should be on one quadrilateral and the fourth somewhere else. The standard way (there are others) is to take x, x˜, x^, x¯ as the initial values and solve the LHS equations of (1.3) for x˜^, x^¯ and x˜¯ and substitute them to the RHS. Then from each of the three RHS equations we can solve for x˜^¯ and the results should be the same, therefore we get two conditions for consistency.

1.3. Organization of the search

The condition that the different ways to compute x˜^¯ provide the same result leads to a set of equations which are polynomials in the initial values x, x˜, x^, x¯. The coefficient of every monomial xαx˜βx^γx¯δ has to vanish separately and this yields polynomial equations on c_i(p, q), c_i(q, r),c_i(r, p), ∀i = 1,...,6.

In analyzing the ensuing equations the solution process branches depending on whether some coefficient function is zero or not. It is not easy to keep track of such branching and therefore we choose here a different approach.

We will do a pre-analysis by setting some coefficients in (1.1) to zero. For that purpose let us write (1.1) as

where ι_j ∈ {0, 1} while c_j are all nonzero. By various choices of ι_j we can have 2⁶ = 64 equations, but by the conditions in Definition 1.1 this actually yields 37 equations. In order to keep track of the different cases it is convenient to name the equations by the list {ι₁, ι₂, ι₃, ι₄, ι₅, ι₆} and if we consider this as a binary number then we can construct also the corresponding decimal number: N = ι₁ + 2ι₂ + 4ι₃ + 8ι₄ + 16ι₅ + 32ι₆. Then we have for example

With this coding the list of the 37 equations passing the conditions 1–3 of Definition 1.1 is given by

• {5,{1,0,1,0,0,0}},{ 7,{1,1,1,0,0,0}},{10,{0,1,0,1,0,0}},{11,{1,1,0,1,0,0}},

• {13,{1,0,1,1,0,0}},{14,{0,1,1,1,0,0}},{15,{1,1,1,1,0,0}},{21,{1,0,1,0,1,0}},

• {23,{1,1,1,0,1,0}},{26,{0,1,0,1,1,0}},{27,{1,1,0,1,1,0}},{29,{1,0,1,1,1,0}},

• {30,{0,1,1,1,1,0}},{31,{1,1,1,1,1,0}},{37,{1,0,1,0,0,1}},{39,{1,1,1,0,0,1}},

• {42,{0,1,0,1,0,1}},{43,{1,1,0,1,0,1}},{45,{1,0,1,1,0,1}},{46,{0,1,1,1,0,1}},

• {47,{1,1,1,1,0,1}},{48,{0,0,0,0,1,1}},{49,{1,0,0,0,1,1}},{50,{0,1,0,0,1,1}},

• {51,{1,1,0,0,1,1}},{52,{0,0,1,0,1,1}},{53,{1,0,1,0,1,1}},{54,{0,1,1,0,1,1}},

• {55,{1,1,1,0,1,1}},{56,{0,0,0,1,1,1}},{57,{1,0,0,1,1,1}},{58,{0,1,0,1,1,1}},

• {59,{1,1,0,1,1,1}},{60,{0,0,1,1,1,1}},{61,{1,0,1,1,1,1}},{62,{0,1,1,1,1,1}},

• {63,{1,1,1,1,1,1}}

Note that 5, 10 and 48 yield two-term equations.

In order to populate the consistency cube we need triplets of equations, leading in principle to 37³ = 50653 cases, fortunately this number can be reduced by applying symmetries.

2.1. Rotation

Let us first consider rotating the cube around the axis (0, 0, 0) − (1, 1, 1), counterclockwise when looking from (1, 1, 1). This moves the quadrilateral planes by (n, m) → (m, k) → (n, k) → (n, m) and certainly does not change the consistency of the equations. Thus we rotate cyclically the shifts ∼ → ^ → ¯ → ∼ accompanied with parameter change p → q → r → p. Due to the cyclic convention in writing the triplet (1.3) it means that the LHS equations become

This rotation does not change the roles of c_i and therefore it corresponds to the rotation of triplet codes by

Clearly ℛ³ = 𝕀. Using rotation symmetry we can reduce the number of equation triples to be analyzed, for example if we always rotate so that X ≤ Y, Z the remaining number of cases is 17575 (in practice we use a slightly different method).

2.2. Tilde-hat reflection

In addition to rotations we need to consider some reflections. Let us first take the reflection across the plane defined by points (0, 0, 0) − (1, 1, 0) − (1, 1, 1), −(0, 0, 1), i.e. the tilde-hat reflection ∼ ↔ ^, accompanied with p ↔ q. From the point of view of the equations, this rule is all we need to know to get another triplet of CAC equations. For the search process we need to know how this and other symmetries change the triplet codes of the equations, because we need to analyze only one triplet from the orbit of a triplet.

Applying tilde-hat reflection on (2.1) we get

Comparing this with (2.1) we find that the back and left equations get exchanged, and in addition the equations change by exchanging the ι or c subscripts 1 ↔ 5 and 3 ↔ 6 (compare, e.g., (2.3c) with (2.1b)). Thus the tilde-hat reflection is given by the operator Θ

where the θ operator changes the binary codes by

The above was for tilde-hat reflection, but we can of course also have hat-bar and bar-tilde reflections, these can be obtained using rotations: ℛΘℛ² and ℛ²Θℛ, respectively.

We have θ² = 𝕀 and hence Θ² = 𝕀, furthermore ℛ²Θ = Θℛ, (Θℛ)² = (ℛΘ)² = 𝕀.

2.3. Tilde reversal

Another important reflection is the shift reversal. Let us consider the tilde reversal, which means reflection across the plane defined by the points (12,0,0)−(12,1,0)−(12,1,1)−(12,0,1), or alternatively, changing all tildes to down-tildes and then applying an overall tilde shift. This will exchange the back and front equations, which were related by shift, and this we can ignore. To see how it changes the bottom and left equation we apply the operation on (2.1) and obtain

Thus for the bottom equation we have 2 ↔ 4 and 5 ↔ 6, the back equation is unchanged, and for the left equation we have 1 ↔ 3 and 2 ↔ 4. We can write the action of the reversal operation 𝒯 as

where and

We have T_aT_b = T_bT_a, Ta2=𝕀, Tb2=𝕀 and therefore 𝒯² = 𝕀. The main relation between T_x and θ is the conjugation θT_aθ = T_b. We also have ℛ𝒯ℛ² = Θ𝒯Θ, or alternatively (Θℛ𝒯)² = 𝕀 or ℛ𝒯ℛ = Θ𝒯ℛΘ, also (𝒯Θ)⁴ = 𝕀.

The above was for tilde-reversal, but we can of course also have hat-reversal and bar-reversal, these can be obtained considering: ℛ𝒯ℛ² and ℛ²𝒯ℛ.

2.4. Inversion

One important transformation that preserves affine linearity of the quad equation (up to overall factor) is the Moebius or affine linear transformation. Now that we are dealing with homogeneous equations only scaling x ↦ const. × x and inversion x ↦ 1/x remain.

From (2.1) we can read that inversion x ↦ 1/x implies 1 ↔ 3, 2 ↔ 4 and 5 ↔ 6. Similar exchanges were also seen in tilde reversal above. Indeed one finds that (ℛ𝒯)³ implies 1 ↔ 3 and 5 ↔ 6, but not 2 ↔ 4. Thus inversion would be an additional symmetry operation in cases where one of c₂, c₄ vanishes but not both. It will turn out that such an asymmetric possibility arises only with solutions in which c₂ and c₄ remain arbitrary functions, because in that case there is also solution in which one of c₂, c₄ vanishes. Such solutions are obtained by simple reductions and therefore we do not actually need inversion.

2.5. Orbit of the transformation group

By composing the basic operators of rotation ℛ and two reflections 𝒯 and Θ in various orders one can generate the orbit of a given initial triplet and for integrability we only need to check one representative of each orbit.

We have already mentioned some relations between the operators (and there may be others) and these can be used to eliminate some combinations as redundant. However, the full group generated by these three operations is probably rather large, but in our case we only need their representation when operating on three strings each containing six zeros or ones. It turns out that (after rotating the triplets into some canonical form) we only have orbits of length 1, 2, 4 or 8. Below we will give the orbits separately, if they have length > 1.

2.6. Gauge transformation

Finally, we can change the coefficients in the equations by a gauge transformation of the form

where we have also included terms with quadratic exponent. However, there are strong conditions on the allowed values for a, b, c, A, B, C, 𝒜, ℬ, 𝒞 because the transformation must not break the uniformity of the 3D lattice, i.e., the translation invariance of the quad equation. (Note that for non-homogeneous equations the gauge freedom would be considerably smaller, perhaps only leaving sign changes.)

In order to easily check the uniformity we have introduced a fixed point n₀, m₀, k₀, and when the gauge transformation is applied to a given equation the result cannot depend on the fixed point, except possibly through an overall multiplier. The restrictions this imposes on the parameters of the gauge transformation will depend on the particular triplet of equations.

We will use gauge transformations to simplify the final result, for example by eliminating some free coefficient(s) from the equation(s), when it simplifies the system. However, even if we can eliminate some free function c_i this way, the mere possibility that such a function can exist is by itself important for classification.

4.1. {63,10,10}

For the first group the most general triplet is

Thus the bottom equation is a completely free homogeneous quadratic equation while the side equations are simple. This is Equation (A.3) in Appendix A.1. The triply shifted x is given by

We see that this equation has the tetrahedron property if c₂ = c₄ = 0. The organization of the 36 sub-cases, obtained by setting some c_j = 0, into orbits is given in Appendix A.2. The simplest equation in this category is (A.2) in Appendix A.1.

For {63,10,10} one cannot change the back or left equations by gauge, but as the number of terms in the bottom equation decreases there is more freedom in the side equations. This is illustrated by a CAC result in the {58,10,10} category:

It has two free functions in the back equation. However it is transformed into a sub-case of (4.1) by the gauge transformation

The possibility of extra free functions that can be gauged away is irrelevant for MDC, but they will have information value in the analysis of BTs in Section 5.4.

4.2. {10,58,15}

This triplet is Equation (A.4) in Appendix A.2. Two of the equations have four free functions each:

Note that the last two equations are not connected by a cyclic variable change but by just tilde-hat exchange (and cyclic parameter change). The triply shifted quantity is given by

The sub-cases of (4.3) are listed in Appendix A.3 and are obtainable by setting some c_j = 0. For each of the four term equations there are 7 acceptable sub-cases but one of them leads to the {63,10,10} category, therefore there are 36 cases in this category. Missing c_j terms sometimes allow a more general bottom equation, but that freedom can be gauged away. TET is possible only if c₂(r, p) = c₄(r, p) = c₅(q, r) = c₆(q, r) = 0 or if c₁(r, p) = c₃(r, p) = c₂(q, r) = c₄(q, r) = 0, but both reduce (4.3) to the {X,10,10} category.

4.3. {53,15,58}

Another chain of triplets that can have free functions is obtained starting with the highest equation {53,15,58} or Equation (A.11) of Appendix A.5. After a gauge transformation the result can be written as

The triply shifted variable is

and clearly has the TET property. The sub-cases obtained by setting some c_i = 0 change the first number of the triplet code and are listed in Appendix A.5. Other types of reductions are given in Figure 2, they are due to limits and are explained below.

Fig. 2:

On the left some sub-case relations from {53,15,58}, on the right the reductions from {63,48,05}. Each box gives the code of the triplet and under it on the left its equation number. The ratio on the right gives above the number of terms and below the number of c₂, c₄ terms.

1. 1)

   {53,11,26} is obtained from {53,15,58} by scaling x¯→rx¯ and then taking the leading term as r → 0.

2. 2)

   {53,05,48} is obtained from {53,15,58} by taking the leading term as r → ∞.

3. 3)

   If r = 0 the triplet (4.4) reduces to a sub-case {53,10,10} of (4.1).

4. 4)

   If in {21,05,48} we take c₅ = 0, we get {05,05,48} which is rotated reflected (A.19).

4.4. {63,48,05}

Here the highest equation is Equation (A.21) of Appendix A.8:

The sign σ₁ cannot be eliminated by gauge.

It has TET if c₄ = 0, which has code {48,48,05}. The diagram on the RHS of Figure 2 describes the reductions of (4.5), they are all by setting some c_i = 0.

4.5. {63,63,63}-2, a triplet without free functions

The highest level triplet of this set has code {63,63,63}-2, it is Equation (A.29) in Appendix. It is in fact Q3(δ = 0) of the ABS list [2]:

The triply shifted variable is

Furthermore, all other triplets that cannot have free functions are obtained by some limit from this equation. Figure 3 describes the chain of limits. From each orbit one element is mentioned and the data in the box is as in Figure 2.

Fig. 3:

Sub-case relations, same notation as in Figure 2. The starting Equation {63,63,63}-2 (4.6) is Q3(δ = 0), {63,63,63}-1 (A.28) is Q1(δ = 0) and {53,53,53} (A.16) is H3(δ = 0). Every reduction simplifies the triplet of equations and eventually we reach equations that can contain free functions c_j.

The limits given in Figure 3 are obtained as follows:

1. 1)

   In (4.6) scale x˜→x˜/p, x^→−qx^ and change r → −1/r and then take the leading terms as p → ∞, q → 0.

2. 2)

   In (A.23) scale x^→qx^ and then take the leading term as q → ∞.

3. 3)

   In (A.13) scale x˜→px˜ and then take the leading term as p → ∞. In addition replace r ↦ 1/r and change signs by gauge.

4. 4)

   In (4.6) scale x¯→x¯/r and then take the leading term as r → ∞. Also change signs of p,q.

5. 5)

   In (A.26) scale x¯→rx¯ and then take the leading term as r → ∞.

6. 6)

   In (4.6) let p ↦ 1 + εp, q ↦ 1 + εq, r ↦ 1 + εr and take the leading term as ε → 0.

7. 7)

   In (A.28) scale x¯→x¯/r, x˜→−x˜, x^→−x^ and then take leading term in r → ∞.

8. 8)

   In (4.6) take the leading term (linear) as p, q, r → 0.

9. 9)

   In (A.16) scale x¯→rx¯ and then take the leading terms as r → 0.

In this way all equations to which free functions of two variables cannot be introduced by a gauge are reductions of Q3(δ = 0) (4.6).

It should be noted, however, that in Figure 3 the bottom box {10,15,58} (A.9) is an equation that does allow free functions, and this happens also on further reductions from the other bottom boxes: From {63,58,15} (A.10) the limit p, q → 0 leads to {53,10,10} while p, q → ∞ leads to {53,48,05} (A.14); from {10,15,58} (A.9) the limit r → leads to {10,05,10} while r → ∞ leads to {10,10,48}, which are in the same orbit; from {53,49,21} (A.17) the limit p, q → ∞ leads to {53,48,05} (A.14). All of these further limits yield equations which allow free functions. It should be noted also that some potential limits actually lead to factorizable equations, for example if in {63,59,31}-1 (A.25) one takes p, q → 0 or if in {63,58,15} (A.10) one takes p ↦ 1 + εp, q ↦ 1 + εq, r ↦ 1 + εr and then take the leading term as ε → 0.

We also note that (4.6) still allows gauge transformations of the form

where σj2=1. This changes the bottom equation into and corresponding results hold for the other equations.

4.6. Comparison with Boll’s results

As was mentioned before, Boll has done a classification of triplets of CAC equations in [3–5]. His setup was both more general (allowing flips between opposing sides and not restricting to homogeneous quadratic equations) and less general (by requiring the tetrahedron property). In [3] the results were collected into Theorems 3.4–3.11. Restricting those equations to quadratic irreducible equations without flips implies some restrictions on the parameters, and yield the following:

1. 1)

   For (3.6),(3.4) we must take δ = 0,ε = 0, then it gives (A.16), i.e. H3(δ = 0). If γ = 0 we get (A.14)

2. 2)

   For (3.10),(3.7) we must take δ = 0, then it gives {47,63,62}, which is reflected rotated (A.26). If furthermore ε = 0 we get rotated (A.10).

3. 3)

   For (3.14),(3.11) we must take δ = 0, ε = 0, then it gives (4.4). If furthermore we take γ = 0 we get (A.14). (Also for γ² = 1 there are no flips but then some equations factorize.)

4. 4)

   In (3.19–20) one additional parameter dependence is missing in B, this is corrected in [4]. In the corrected form it is necessary to take δ₃ = 0, which yields reflected (A.23). If after this we take δ₁ = 0 or δ₂ = 0 we get a reflected forms of (A.13) and if both coefficients vanish we get (A.9).

5. 5)

   For (3.21–22) (of [3], not [4]) we must take δ₁ = 0 or change δx₁x₂ into δx₁x, which is again a special case of (4.4). If α² = 1 some equations factorize.

6. 6)

   For (3.29–30) we must take δ₂ = 0, then it gives reflected rotated (A.17). If furthermore δ₁ = 0 we get rotated (A.14).

7. 7)

   For (3.31–32) we must take δ₁ = δ₂ = 0, then it gives rotated (A.15). In that case all coefficients in B of (3.31) can in fact be free. If in (3.31) we take α = 0 we get rotated (A.19).

8. 8)

   For (3.33–34) we must take δ₁ = 0, then it gives reflected (A.18). In B all coefficients can be free.

9. 9)

   For (3.41–42) we must take δ₁ = 0, then it corresponds to {42,21,14}, a rotated form of which appears in the orbit of (A.24). Next if δ₃ = 0 we get rotated reflected (A.22) and if δ₂ = 0 we get subcase {53,10,10} of (4.1).

10. 10)

    For (3.43–44) we must take δ₁ = δ₂ = 0, then it becomes rotated form of (A.20). (Note that in the formula for C the first variable must be x₁₃, not x₁₂.)

Thus all homogeneous quadratic shift invariant equations in Boll’s classification [3] are also included in our classification. However, there is one result we could not identify among Boll’s results, namely (A.25), which is related to Q1 the same way that (A.26) is related to Q3, and (A.16) is related to H3.

4.7. Triplets with two term equations

As another selection of our results we would like to collect here the simplest CAC triplets, those for which all equations have just two terms. Below on the LHS they are given in the general form after all equations from the CAC condition have been solved, and then on the RHS we present the simplest form obtainable by a gauge transform.

• {10,10,10}

  {c2(p,q)x˜x^+c4(p,q)xx˜^=0,c2(q,r)x^x¯+c4(q,r)xx^¯=0,c2(r,p)x¯x˜+c4(r,p)xx˜¯=0,→{x˜x^−xx˜^=0,x^x¯−xx^¯=0,x¯x˜−xx˜¯=0.(4.8)

• {5,10,10} or {48,10,10} by tilde-hat reflection

  {c1(p,q)xx˜+c3(p,q)x^x˜^=0,x^x¯+σxx^¯=0,c2(r,p)x¯x˜+c4(r,p)xx˜¯=0,→{c1(p,q)xx˜+c3(p,q)x^x˜^=0,x^x¯−xx^¯=0,x¯x˜−xx˜¯=0.(4.9)

• {10,48,5}

  {x˜x^+σxx˜^=0,c5(q,r)xx¯+c6(q,r)x^x^¯=0,c1(r,p)xx¯+c3(r,p)x˜x˜¯=0,→{x˜x^−xx˜^=0,c5(q,r)xx¯+c6(q,r)x^x^¯=0,c1(r,p)xx¯+c3(r,p)x˜x˜¯=0.(4.10)

• {5,5,48} and {48,5,48} by tilde-hat reflection

  {c1(p,q)xx˜+c3(p,q)x^x˜^=0,c(r)xx^+x¯x^¯=0,σc(r)xx˜+x¯x˜¯=0,→{c1(p,q)xx˜+c3(p,q)x^x˜^=0,xx^−x¯x^¯=0,xx˜−x¯x˜¯=0.(4.11)

These cases are easily distinguished by the number of c₂, c₄ terms, which have also been rotated into a particular position. These triplets can be obtained from several results listed in the Appendix by various limits.

5.1. General setup

Let us fix the equations and some concepts. We have six equations, which we can arrange as follows:

This is the same set as (1.3), except that we have named x¯=y, which has no effect on the algebra. We see that the top equation (5.1a) depends only on the y variables, the bottom equation (5.1d) only on the x variables, while the side equations (5.1b),(5.1c) depend on two x and on two y variables.

5.2. Computations for eqBT

In order to further analyze the eqBT, it is useful to develop the formulae to some extent without specifying the form of the Q-polynomials. For this purpose it is often useful to isolate one of the variables of Q, but since Q is multilinear this is easy. We can write, for example,

where

Using this notation we can solve the LHS equations of (5.1b), (5.1c),

and then from the RHS equations we can solve

Let us furthermore expand these in y since it cannot appear in final equation. We have

where

We use the symbol ≗ for equalities that should hold modulo an acceptable equation. For a set of equations the equality should be modulo the same equation, which we call Q₁₂.

Since Q₁₂ does not contain y, the Bäcklund condition y˜^23≗y˜^31 implies three equations

For a genuine eqBT these equations should not all be satisfied automatically but rather for at least one equation the LHS should factor with an acceptable equation as a factor.