Exact analytical expressions for the matrix elements of the Uehling potential in a basis of explicitly correlated exponential wave functions are presented. The obtained formulas are then used to compute with an improved accuracy the vacuum polarizati
on correction to the binding energy of muonic and pionic molecules, both in a first-order perturbative treatment and in a nonperturbative approach. The first resonant states lying below the n=2 threshold are also studied, by means of the stabilization method with a real dilatation parameter.
