The following question arises when I attempt to understand the modular parameterization of the elliptic curve $$E:y^2-y=x^3-x$$


In [Mazur-Swinnerton-Dyer](https://eudml.org/doc/142281) and [Zagier](http://people.mpim-bonn.mpg.de/zagier/files/doi/10.1007/BFb0084592/fulltext.pdf)'s construction, a theta function associated with a positive definite quadratic form is induced:


$$\theta(q)=\sum\_{x\in\mathbb{Z}^4}q^{\frac{1}{2}x^{T}Ax}$$


where $$A=\left(\begin{matrix}2 & 0 & 1 & 1\\
0 & 4 & 1 &2\\
1 & 1 & 10 & 1\\
1 & 2 & 1 & 20
\end{matrix}\right)$$


$A$ is a positive definite matrix of determinant $37^2$, and we have $37A^{-1}=K^TAK$ where $K$ is an integral matrix of determinant $\pm 1$.


**Question:** Suppose $A$ is a positive definite $4\times 4$ matrix with integral entries. All diagonal entries are even numbers. The determinant of $A$ is a square number $N^2$. Is it true that for every $N=p$ ($p>2$ is a prime number), there is *at least one* $A$ that $NA^{-1}=K^TAK$, where $K$ is an integral matrix of determinant $\pm 1$?



Integral quaternary forms and theta functions
The answer is true, using the following construction.
Let $B$ be the quaternion algebra of discriminant $p$ and let $O$ be a maximal order with an element $x$ satisfying $x^2 = -p$. The reduced norm is a quadratic form on $O$, with positive definite Gram matrix $A$ of determinant $p^2$. The matrix $A^{-1}$ then represents the reduced norm on the dual lattice $O^{\sharp}$. Recall that $\text{nrd}(\text{diff}(O)) = \text{discrd}(O) = p$ and since $O$ is maximal, $\text{diff}(O)$ is invertible and $O = \text{diff}(O) O^{\sharp}$ (see e.g.
*Voight, John*, [**Quaternion algebras**](http://dx.doi.org/10.1007/978-3-030-56694-4), [ZBL07261776](https://zbmath.org/?q=an:07261776). Section 16.8), hence
$$
p O^{\sharp} \subseteq \text{diff}(O) O^{\sharp} = O
$$
which implies that $p \cdot (O^{\sharp} / O) = 0$ (therefore $O^{\sharp} / O$ is an abelian group of exponent $p$ and size $p^2$, so isomorphic to $(\mathbb{Z} / p \mathbb{Z})^2$). Next, we note that the matrix $A^{-1}$ is also the change of basis matrix between the chosen basis of $O$ and its dual, therefore $pA^{-1}$ is integral.


This is not enough, but so far we have not used the element $x$. The matrix $pA^{-1}$ is the matrix representing the norm form on the ideal $xO^{\sharp}$, since $x^2 = -p$. But $nrd(xO^{\sharp})=p$ is a maximal order with an element $x$ such that $x^2 = -p$. By a theorem of Ibukiyama (see reference below), it is isomorphic (hence isometric as lattices) to $O$.


This could also be seen directly (referring to some of the answers above - it is possible to do it using a finite number of cases) by explicitly writing down the maximal orders. I will list below explicit constructions, which are based on
*Ibukiyama, Tomoyoshi*, [**On maximal orders of division quaternion algebras over the rational number field with certain optimal embeddings**](http://dx.doi.org/10.1017/S002776300002016X), Nagoya Math. J. 88, 181-195 (1982). [ZBL0473.12012](https://zbmath.org/?q=an:0473.12012). -


If $p \equiv 3 \bmod 4$, then $B = (-p,-1)$, and $O = \mathbb{Z}<(1+i)/2,j>$. In this case, we have
$$
A = \left( \begin{array}{cccc} 
2 & 1 & 0 & 0 \\ 
1 & \frac{p+1}{2} & 0 & 0 \\
0 & 0 & 2 & 1 \\
0 & 0 & 1 & \frac{p+1}{2}
\end{array} \right)
,
pA^{-1} = \left( \begin{array}{cccc} 
\frac{p+1}{2} & -1 & 0 & 0 \\ 
-1 & 2 & 0 & 0 \\
0 & 0 & \frac{p+1}{2} & -1 \\
0 & 0 & -1 & 2
\end{array} \right),
K = \left( \begin{array}{cccc} 
-1 & 1 & 0 & 0 \\ 
1 & 0 & 0 & 0 \\
0 & 0 & -1 & 1 \\
0 & 0 & 1 & 0
\end{array} \right)
$$


If $p \equiv 1 \bmod 4$, then choosing $q \equiv 3 \bmod 4$ such that
$\left( \frac{q}{p} \right) = -1$, we can find $c$ such that
$c^2 \equiv -p \bmod q$, and then $B = (-p,-q)$ and
$$
O = \mathbb{Z} \oplus \mathbb{Z} \frac{1+j}{2} \oplus
\mathbb{Z} \frac{i(1+j)}{2} \oplus \mathbb{Z} \frac{(c+i)j}{q}
$$
In this case, we compute that
$$
A = \left( \begin{array}{cccc} 
2 & 1 & 0 & 0 \\ 
1 & \frac{q+1}{2} & 0 & c \\
0 & 0 & \frac{p(q+1)}{2} & p \\
0 & c & p & \frac{2(p+c^2)}{q}
\end{array} \right)
,
pA^{-1} = \left( \begin{array}{cccc} 
\frac{(q+1)(c^2+p)}{2q} & -c-\frac{c^2+p}{q} & -c & \frac{c(q+1)}{2} \\ 
-c-\frac{c^2+p}{q} & 2(c^2 + \frac{c^2+p}{q}) & 2c & -c(q+1) \\
-c & 2c & 2 & -q \\
\frac{c(q+1)}{2} & -c(q+1) & -q & \frac{q(q+1)}{2}
\end{array} \right)
$$


and if $\{e\_1, e\_2, e\_3, e\_4 \}$ is the above basis for $O$, then we see that $ \{e\_2 e\_4, ce\_1 - e\_4, e\_1, j e\_2 \} $ is a basis for the resulting module, which written in terms of the original basis yields the matrix
$$
K = \left( \begin{array}{cccc}
0 & c & 1 & -\frac{q+1}{2} \\
-c & 0 & 0 & 1 \\
-1 & 0 & 0 & 0 \\
\frac{q+1}{2} & -1 & 0 & 0
\end{array}
 \right).
$$


Referring to Will Jagy's wondering in the first answer, the reason that this does not work for most lattices is that they do not correspond to a maximal order (as in the answer by few\_reps, the quotient of the lattices is actually cyclic) and there are only one or two (depending on $p \bmod 4$ isomorphism classes of maximal orders which contain a root of $-p$. For example, for $p = 37$, there are only two isomorphism classes of maximal orders in the quaternion algebra, and only one of them contains a root of $-p$.


