  The three-body decays $B^0_s \rightarrow \psi(2S,3S) \pi^+ \pi^-$ are studied
based on the perturbative QCD approach. With the help of the nonperturbative
two-pion distribution amplitudes, the analysis is simplified into the
quasi-two-body processes. Besides the traditional factorizable and
nonfactorizable diagrams at the leading order, the next-to-leading order vertex
corrections are also included to cancel the scale dependence. The $f_0(980)$,
$f_0(1500)$ resonance contributions as well as the nonresonant contributions
are taken into account using the presently known $\pi\pi$ time-like scalar form
factor for the $s\bar{s}$ component. It is found that the predicted $B^0_s
\rightarrow \psi(2S) \pi^+ \pi^-$ decay spectra in the pion pair invariant mass
shows a similar behavior as the experiment. The calculated S-wave contributions
to the branching ratio of $B^0_s \rightarrow \psi(2S) \pi^+ \pi^-$ is
$6.0\times 10^{-5}$, which is in agreement with the LHCb data $\mathcal
{B}(B^0_s \rightarrow \psi(2S) \pi^+ \pi^-)=(7.2\pm 1.2)\times 10^{-5} $ within
errors. The estimate of $\mathcal {B}(B^0_s \rightarrow \psi(3S) \pi^+ \pi^-)$
can reach the order of $10^{-5}$, pending for the corresponding measurements.
