I am looking for a straightforward way to upper bound the covering number of a $d$-dimensional euclidean ball by $\ell\_\infty$-balls of radius $\varepsilon$, which I will call cubes of sidelength $2\varepsilon$ for clarity. Let us denote this number by $\mathcal N(\varepsilon)$.


An elementary upper bound is to say $\mathcal N(\varepsilon) \leqslant 1 / \varepsilon^d$, as we can cover the whole cube of sidelength 2 with this many small cubes.
However, I know that the precise asymptotic behavior, up to multiplicative and additive constants, is the following
$$
\log \mathcal N(\varepsilon) \approx \left\{\begin{split}
\frac{1}{\varepsilon^2} \log d \varepsilon^2\quad \text{if} \quad \varepsilon \geqslant 1/\sqrt{d}\\
d \log \frac{1}{d\varepsilon^2} \quad \text{if} \quad \varepsilon \leqslant 1/\sqrt{d} \, .
\end{split}\right. \
$$
Therefore the easy upper bound of $d \log 1/ \varepsilon$ is loose when $\varepsilon$ is large. The change of regime happens when the diagonal of the smaller cubes become comparable to the radius of the larger ball.


The only proof I found, in online teaching notes [1], is convoluted. It consists in first bounding the covering number of $\ell\_1$ balls by $\ell\_2$ balls, thanks to Maurey's empirical method, and then to appeal to a duality result of [2].


While the proof is elegant and sophisticated, the inconvenient is that it yields non-explicit constants. Also, it does not seem to take advantage of the ``easiness'' of the problem: covering with hypercubes should not be too hard as the cubes fit well together. I tried counting the number of cubes needed to cover the ball in the natural covering (cutting the cube of sidelength $2$ into small cubes of length $2\varepsilon$ ) but I have no idea how one would do so. For all I know, it might even be that the optimal covering has some overlap between the cubes.


EDIT: A nice improvement on the trivial bound by Rémi Peyre.
Considering the natural covering, it suffices to count the cubes that lie inside the ball of radius $1 + 2 \varepsilon \sqrt{d}$. Therefore
\begin{equation}
\mathcal N(\varepsilon)
\leq \frac{1}{2 \epsilon^d} \mathrm{Vol}\big(B(1 + 2\varepsilon\sqrt{d}) \big)
\leq \frac{1}{\sqrt{d\pi}} (2 \pi e)^{d/2} \Big( 1 + \frac{1}{2\varepsilon \sqrt{d}} \Big)^d \, .
\end{equation}
using a non-asymptotic version of Stirling. Taking logs we obtain something of order $d ( c + \log(1 + 1/(\varepsilon \sqrt{d}) )$. This is an improvement but is still far from the refined bounds above, when $\varepsilon \sqrt{d}$ is large.


[1] [http://www.stat.yale.edu/~yw562/teaching/598/lec15.pdf](http://www.stat.yale.edu/%7Eyw562/teaching/598/lec15.pdf)


[2] Duality of metric entropy, S. Artstein, V. Milman, and S. J. Szarek
<https://annals.math.princeton.edu/wp-content/uploads/annals-v159-n3-p07.pdf>



Elementary precise estimate of the covering number of euclidean balls by hypercubes
Once you know the answer, the proof is a trivial induction on $d$. We will show that for some positive constants $A,B,K$ (to be chosen in the end) one can cover the ball of radius $r$ by
$$
F(d,R^2)=\left(1+\frac{d}{R^2}\right)^{BR^2}+\left(\frac{AR^2}{d}\right)^{d/2}=F\_1(d,R^2)+F\_2(d,R^2).
$$
cubes with sidelength $2K$ where $R^2=r^2+K^2$. If you are not too picky about constants, this is equivalent to what you have requested.


We will notice first of all that $F(d,R^2)=0$ for $R<K$ ($r$ is killed for $R=K$ already, so we do not need to take any special care of too small values of $R$ in any dimension. We shall also notice that everything is fine in dimensions $1,2$, so we'll take $d\ge 2$ for the induction step.


Now take the $d+1$-dimensional ball of radius $r$, cover the central slice of width $2K$ by $F(d,R^2)$ cubes and then go up and down by $K$ looking at how the radii of the cross-sections change and using the optimal cover each time. We'll get
$$
F(d,R^2)+2\sum\_{m\ge 1} F(d,R^2-m^2K^2)
$$  

All we need to do is to compare it to $F(d+1,R^2)$. We want the domination to be overwhelming (in a sense), so let us look at each term separately first.


*Term 1:* We have
$$
F\_1(d+1,R^2)-F\_1(d,R^2)\ge B\left(1+\frac{d}{R^2}\right)^{BR^2-1}
$$
(mean value theorem).
On the other hand, the derivative of the function $f(x)=Bx\log(1+\frac dx)$ is $B[\log(1+\frac dx)-\frac{d}{d+x}]\ge B[\log(1+\frac d{R^2})-\frac d{d+R^2}]=B\lambda$ on $[0,R^2]$, so the sum in $m$ for $F\_1$ is dominated by
$$
\left(1+\frac{d}{R^2}\right)^{BR^2}\int\_{K^2/2}^\infty e^{-B\lambda t}dt=\left(1+\frac{d}{R^2}\right)^{BR^2}e^{-BK^2\lambda/2}\frac 1{B\lambda}
\\
\le
\left(1+\frac{d}{R^2}\right)^{BR^2-1}e^B e^{-(BK^2/2-1)\lambda}\frac 1{B\lambda}\le \left(1+\frac{d}{R^2}\right)^{BR^2-1}e^B e^{-BK^2\lambda/4}\frac 1{B\lambda}
$$
If $\lambda\ge\frac 1K$, we have the factor
$$
e^{-(K/4-1)B}\frac KB\ll B
$$
provided that $K$ is large enough (even in terms of $B$, if you want; we'll choose $K$ last).


So, we are in danger here only if $\lambda<\frac 1K$, in which case, since $\log(1+u)-\frac u{u+1}\approx u^2/2$ near $u=0$, we have $\frac{d}{R^2}\le\sqrt{\frac 6K}$, say, and we shall use the trivial estimate
$R^2e^{Bd}$ for the sum (each term is at most $e^d$ and there are not more than $R^2$ terms).


*Term 2:*
$$
F\_2(d+1,R^2)=\left(\frac{AR^2}{d+1}\right)^{(d+1)/2}\ge
\sqrt{\frac{AR^2}{d+1}}\left(\frac d{d+1}\right)^{d/2}\left(\frac{AR^2}{d}\right)^{d/2}
\\
\ge \sqrt{\frac{AR^2}{2de}}F(d,R^2)
$$
If $\frac{d}{R^2}<\sqrt{\frac 6K}$, this is at least
$$
\sqrt{\frac{AR^2}{2de}}\left(\frac{AR^2}{d}\right)^{d/2}\ge
\frac{AR^2}{2de}\left(\frac{AR^2}{d}\right)^{(d-1)/2}
\ge \frac{AR^2}{2de}(A\sqrt K/3)^{d/4}\gg R^2e^{Bd}
$$
if $A$ is not too small and $A\sqrt K/(3e^4)\gg e^{4B}$, which can be achieved by increasing $K$ for any fixed $A,B$.


Now it is time to estimate the sum in $m$. Again, comparing to the exponent, we see that we get
$$
\le F\_2(d,R^2)\sum\_{m\ge 1} e^{-\frac{m^2K^2d}{2R^2}}\le \sqrt{\frac{R^2}d} 
$$
so we need to compare
$
\sqrt{\frac{AR^2}{2de}}
$
and $1+2\sqrt{\frac{R^2}d}$.
If $d\le R^2$, the first expression wins big time. If not, we have at most $3\left(\frac{AR^2}{d}\right)^{d/2}$. We can play it against $B\left(1+\frac{d}{R^2}\right)^{BR^2-1}$ now. $B$ certainly beats $3$. If $d>AR^2$, the comparison is trivial. Finally, if $R^2\le d\le AR^2$, we have at most $A^{AR^2/2}$ on the left and at least $2^{BR^2/2}$ on the right, so if $A^A\le 2^B$, we are fine. This tells the relation between $A$ and $B$ (both large) and $K$ is chosen last.


The end.


P.S. If you think a bit, you see that such induction is doomed to work: no fancy cover is possible because for the cubes of twice smaller size the covers of our slices are disjoint (the cube just cannot stretch enough in the vertical direction), so the stupid algorithm (cover slice by slice reducing the dimension by 1, then do the same for each slice of slice, etc.) is actually the best.


