Constructing 3x3x3 magic cube requires solving 49 linear equations for 27 unknowns. On top of that, we are looking for an all integer solution. The number of equations is given by 5 equations per 3x3 layer, times 3 layers per dimension, times 3 dimensions, for a total of 5x3x3 = 45. Add 4 equations from the body diagonals for a total of 49. The 49 equations may be considered as 49 vectors in a 27-dimensional space. Maximum of 27 can be linarly independent. The corresponding homogenous system (right sides of all equations equal to zero) obviously has the trivial solution (x1 = x2 = ... = x27 = 0) and the non-homogenous system (all right sides equal to the magic number M) also has the trivial solution (x1 = x2 = ... = x27 = M/3). Therefore, no linearly dependent equations have incommpatible right sides (such as 2x1 + x2 = 3, 2x1 + x2 = 4). If we find 27 linearly independent vectors (equations), the trivial solution is the only solution and it obviously does not satisfy the magic cube conditions (all numbers must be mutually different). Number the small 1x1x1 cubes on any 3x3x3 cube face (for example, the top horizontal face) as follows: 1 2 3 4 5 6 7 8 9 Write down the 6 vectors a1, a2, ..., a6 corresponding to the 6 possible equations for the left-right and forward-backward rows, including the right sides (the magic number M): 1 0 0 1 0 0 1 0 0 M -> a1 0 1 0 0 1 0 0 1 0 M -> a2 0 0 1 0 0 1 0 0 1 M -> a3 1 1 1 0 0 0 0 0 0 M -> a4 0 0 0 1 1 1 0 0 0 M -> a5 0 0 0 0 0 0 1 1 1 M -> a6 The remaining 18 coordinates of these vectors are all zeros. It is obvious that these 6 vectors are linearly dependent, because: a1 + a2 + a3 - (a4 + a5 + a6) = 0 So lets kick out the vector a4. 1 0 0 1 0 0 1 0 0 M -> a1 0 1 0 0 1 0 0 1 0 M -> a2 0 0 1 0 0 1 0 0 1 M -> a3 0 0 0 1 1 1 0 0 0 M -> a5 0 0 0 0 0 0 1 1 1 M -> a6 The remaining 5 vectors are linearly independent because of the triangular shape of the matrix (we may exchange the 6th and 7th columns to see the triangular shape). Add the 2 possible diagonals: 1 0 0 0 1 0 0 0 1 M -> a7 1 0 0 1 0 0 1 0 0 M -> a1 0 1 0 0 1 0 0 1 0 M -> a2 0 0 1 0 0 1 0 0 1 M -> a3 0 0 1 0 1 0 1 0 0 M -> a8 0 0 0 1 1 1 0 0 0 M -> a5 0 0 0 0 0 0 1 1 1 M -> a6 Transform the matrix to the triangular shape: 1 0 0 1 0 0 1 0 0 M -> a1 0 1 0 0 1 0 0 1 0 M -> a2 0 0 1 0 0 1 0 0 1 M -> a3 0 0 0 1 -1 0 1 0 -1 0 -> a1-a7 0 0 0 0 1 -1 1 0 -1 0 -> a8-a3 0 0 0 0 0 3 -3 0 3 M -> a5-(a1-a7)-2(a8-a3) 0 0 0 0 0 0 1 1 1 M -> a6 Although it is not essential, this matrix can be further reduced to: 1 0 0 0 0 0 0 0 1 2/3M -> a7-(a8-a3)-1/3[a5-(a1-a7)-2(a8-a3)] 0 1 0 0 0 0 0 1 0 2/3M -> a2-(a8-a3)-1/3[a5-(a1-a7)-2(a8-a3)] 0 0 1 0 0 0 1 0 0 2/3M -> a3-1/3[a5-(a1-a7)-2(a8-a3)] 0 0 0 1 0 0 1 0 -1 1/3M -> a1-a7+(a8-a3)+1/3[a5-(a1-a7)-2(a8-a3)] 0 0 0 0 1 0 0 0 0 1/3M -> a8-a3+1/3[a5-(a1-a7)-2(a8-a3)] 0 0 0 0 0 1 -1 0 1 1/3M -> 1/3[a5-(a1-a7)-2(a8-a3)] 0 0 0 0 0 0 1 1 1 M -> a6 We have 7 linearly independent vectors from the first (top) 3x3 layer that span a 7-dimensional subspace of the 27 dimensional space. We can get 7 linearly independent vectors from each of the three parallel horizontal 3x3 layers, for a total of 21 linearly independent vectors, because the three 7-dimensional subspaces are mutually ortogonal. We have to find 6 more vectors without causing any linear dependence. We have 9 columns, 12 vertical 3x3 layer diagonals, and 4 body diagonals of the 3x3x3 cube to work with. All these 25 remaining 27-dimensional vectors take one (and only one) non-zero coordinate from each of the three 7-dimensional subspaces. The remaining linearly independent vectors either have the 8th and/or 9th coordinates in each 7-dimensional subspace non-zero or we can reduce them to such vectors by separate linear combinations within the subspaces. Let's see if we can pick up the 6 remaining linearly independent vectors: 8 9 17 18 26 27 <- coordinate number ----------------------- 1 0 1 0 1 0 <- column 8 0 1 0 1 0 1 <- column 9 0 1 1 0 0 0 <- face diagonal from the top of column 9 across column 8 0 0 1 0 0 1 <- face diagonal from the bottom of column 9 across column 8 0 1 0 0 0 0 <- face diagonal from the top of column 9 across column 6 0 0 0 0 0 1 <- face diagonal from the bottom of column 9 across column 6 0 1 0 0 0 0 <- body diagonal from the top of column 9 0 0 0 0 0 1 <- body diagonal from the bottom of column 9 1 0 0 0 0 0 <- face diagonal from the top of column 8 0 0 0 0 1 0 <- face diagonal from the bottom of column 8 We can eliminate either the face diagonal from column 9 across column 6 of the body diagonals, because they have the same 8th and 9th coordinates in each subspace. Let's eliminate the former. 1 0 0 0 0 0 <- v1: FD top C8 1 0 1 0 1 0 <- v2: C8 0 1 0 0 0 0 <- v3: BD top C9 0 1 0 1 0 1 <- v4: C9 0 1 1 0 0 0 <- v5: FD top C9 0 0 1 0 0 1 <- v6: FD bottom C9 0 0 0 0 0 1 <- v7: BD bottom C9 Transform the matrix to the triangular shape: 1 0 0 0 0 0 <- v1 0 1 0 0 0 0 <- v3 0 0 1 0 0 0 <- v5-v3 0 0 0 1 0 0 <- v4-v3-v7 0 0 0 0 1 0 <- v2-v1-(v5-v3) 0 0 0 0 0 1 <- v6-(v5-v3) 0 0 0 0 0 1 <- v7 Kicking out either the vector v6 (face diagonal from the bottom of column 9) or the vector v7 (body diagonal from the bottom of column 9) leaves us with the required 6 linearly independent vectors. Together with the 21 linearly independent equation vectors obtained earier, we have a total of 27 linarly independent vectors or equations, and therefore the trivial solution (x1 = x2 = ... = x27 = M/3) is the only solution. We do not even have to restrict ourselves to the sequential numbers 1, 2, ... , 27.