Computing all "suboptimal" rational approximations to $\pi/2$

Question

I have an irrational number $\alpha$ ($\alpha=\frac\pi2$), and I would like to determine all integers $n\in[1,N]$ ($N=10^{16}$) that satisfy $$ n \epsilon(n)^2 \leq \tau $$ where $\tau$ is a known real number ($\tau=78$), and $\epsilon(n)$ is the distance between $n$ and the closest multiple of $\alpha$: $$ n = m \alpha + \epsilon(n), \qquad |\epsilon(n)| < \tfrac12\alpha. $$

I know that if I approximate $\alpha$ by a rational number $A/B$, then the problem can be solved with extended Euclid's algorithm applied to $$ Bn-Am = B\epsilon, $$ where $B\epsilon$ is an integer, with an upper bound derived from $n\epsilon^2<\tau$, $B\epsilon<B\sqrt{\tau}$. For each value of $B\epsilon$ I can compute corresponding solutions $n\in[1,N]$.

The problem is that if $A/B$ is a sufficiently accurate rational approximation to $\alpha$ (i.e., $B\sim N$ ensures that $n\epsilon^2$ is approximated accurately), then the bound $B\epsilon < B\sqrt\tau$ (e.g., $B\sqrt\tau=10^{17}$) leads me to consider prohibitively many different $n$.

The standard theory of rational approximations tells me how to compute the best rational approximations with continued fractions, but I need to compute all rational approximations $n/m$ that are good enough in the above sense. If I write $$ \alpha = \frac nm-\frac \epsilon m = \frac nm+\delta, $$ the equivalent problem is to compute all $(n,m)$ with $$ \big|\alpha - \frac nm\big| < \frac{\sqrt\tau}{m \sqrt{n}} \sim \frac{\sqrt{\tau/\alpha}}{m^{3/2}}, $$ but the right-hand side here is quite different from $1/m^2$ (or $1/m^\mu$ where $\mu\geq2$ is the irrationality measure of $\alpha$). Since $\tfrac32\leq\mu$, that tells me there are infinitely many solutions, but I don't know how to actually compute them.

Is there an efficient algorithm for computing $n$ with $n\epsilon^2<\tau$?

Answer 1

As noted in the comments (by Felipe Voloch and myself), one expects about $O_\tau(N^{1/2})$ solutions, because that's the area $|R|$ of the region $$ R = \{ (m,n) \in {\bf R}^2 \colon 1 \leq n \leq N, \ (\alpha m - n)^2 < \tau/n \} $$ whose intersection with ${\bf Z}^2$ we want to list. (In the comment I got the dependence on $\tau$ wrong: it's $(\tau N)^{1/2}$, not $\tau N^{1/2}$, because the upper bound on $|\alpha m - n|$ is $\sqrt{\tau/n}$, not $\tau\,/\sqrt n$ as I misread; so for $(N,\tau) = (10^{16},78)$ we expect about $10^9$ pairs $(m,n)$, not $10^{10}$.)

To list $R \cap {\bf Z}^2$ in $\,\tilde{\!O}(|R|)$ time:

i) Partition $R$ into about $\log_2 N$ parts $R_1,R_2,\ldots R_k$ where $R_j = \{(m,n) \in R \colon 2^{-j} N \leq n < 2^{1-j} N$. The $R_j$ are all equivalent under affine transformations (except for $R_k$ [unless $N$ is exactly a power of $2$], but $R_k$ is negligible anyway). Then do the remaining steps for each $j$.

ii) Inscribe $R_j$ in an ellipse $E_j$ whose area is $O(|R_j|)$. For example $E_j$ might be $$ \{ (x,y) \in {\bf R}^2 \colon \frac{4^{j+1}}{N^2}\Bigl(x - \frac32 2^{-j} N\Bigr)^2 + \frac {2^{j-1}\tau}{N} (\alpha x - y)^2 \leq 2 \} $$ if I did the algebra right.

iii) Choose an affine transformation $T_j$ that takes $E_j$ to the unit circle $D_1: X^2 + Y^2 \leq 1$.

iv) Apply $2$-dimensional lattice basis reduction (a.k.a.\ the Euclidean algorithm) to $T_j {\bf Z}^2 - T_j (0,0)$. Use the reduced basis to efficiently list $D_1 \cap T_j {\bf Z}^2$, and thus $E_j \cap {\bf Z}^2$.

v) Check whether each of the resulting integer vectors $(m,n)$ is actually in $R_j$, and output those that pass this test (which should be roughly the constant fraction $|R_j| \, / \, |E_j|$), QEF.