  Consider the multivariate Stein equation $\Delta f - x\cdot \nabla f = h(x) -
E h(Z)$, where $Z$ is a standard $d$-dimensional Gaussian random vector, and
let $f\_h$ be the solution given by Barbour's generator approach. We prove
that, when $h$ is $\alpha$-H\"older ($0<\alpha\leq1$), all derivatives of order
$2$ of $f\_h$ are $\alpha$-H\"older {\it up to a $\log$ factor}; in particular
they are $\beta$-H\"older for all $\beta \in (0, \alpha)$, hereby improving
existing regularity results on the solution of the multivariate Gaussian Stein
equation. For $\alpha=1$, the regularity we obtain is optimal, as shown by an
example given by Rai\v{c} \cite{raivc2004multivariate}. As an application, we
prove a near-optimal Berry-Esseen bound of the order $\log n/\sqrt n$ in the
classical multivariate CLT in $1$-Wasserstein distance, as long as the
underlying random variables have finite moment of order $3$. When only a finite
moment of order $2+\delta$ is assumed ($0<\delta<1$), we obtain the optimal
rate in $\mathcal O(n^{-\frac{\delta}{2}})$. All constants are explicit and
their dependence on the dimension $d$ is studied when $d$ is large.
