  We study a random matching problem on closed compact $2$-dimensional
Riemannian manifolds (with respect to the squared Riemannian distance), with
samples of random points whose common law is absolutely continuous with respect
to the volume measure with strictly positive and bounded density. We show that
given two sequences of numbers $n$ and $m=m(n)$ of points, asymptotically
equivalent as $n$ goes to infinity, the optimal transport plan between the two
empirical measures $\mu^n$ and $\nu^{m}$ is quantitatively well-approximated by
$\big(\mathrm{Id},\exp(\nabla h^{n})\big)_\#\mu^n$ where $h^{n}$ solves a
linear elliptic PDE obtained by a regularized first-order linearization of the
Monge-Amp\`ere equation. This is obtained in the case of samples of correlated
random points for which a stretched exponential decay of the $\alpha$-mixing
coefficient holds and for a class of discrete-time Markov chains having a
unique absolutely continuous invariant measure with respect to the volume
measure.
