  We consider the electrostatic Born-Infeld energy \begin{equation*}
\int_{\mathbb{R}^N}\left(1-{\sqrt{1-|\nabla u|^2}}\right)\, dx
-\int_{\mathbb{R}^N}\rho u\, dx, \end{equation*} where $\rho \in
L^{m}(\mathbb{R}^N)$ is an assigned charge density, $m \in [1,2_*]$,
$2_*:=\frac{2N}{N+2}$, $N\geq 3$. We prove that if $\rho \in L^q(\mathbb{R}^N)
$ for $q>2N$, the unique minimizer $u_\rho$ is of class
$W_{loc}^{2,2}(\mathbb{R}^N)$. Moreover, if the norm of $\rho$ is sufficiently
small, the minimizer is a weak solution of the associated PDE
  \begin{equation}\label{eq:BI-abs} \tag{$\mathcal{BI}$}
-\operatorname{div}\left(\displaystyle\frac{\nabla u}{\sqrt{1-|\nabla
u|^2}}\right)= \rho \quad\hbox{in }\mathbb{R}^N, \end{equation} with the
boundary condition $\lim_{|x|\to\infty}u(x)=0$ and it is of class
$C^{1,\alpha}_{loc}(\mathbb{R}^N)$, for some $\alpha \in (0,1)$.
