  Isolated, undamped geodesic-acoustic-mode (GAM) packets have been
demonstrated to obey a (focusing) nonlinear Schr\"odinger equation (NLSE) [E.
Poli, Phys. Plasmas 2021]. This equation predicts susceptibility of GAM packets
to the modulational instability (MI). The necessary conditions for this
instability are analyzed analytically and numerically using the NLSE model. The
predictions of the NLSE are compared to gyrokinetic simulations performed with
the global particle-in-cell code ORB5, where GAM packets are created from
initial perturbations of the axisymmetric radial electric field $E_r$. An
instability of the GAM packets with respect to modulations is observed both in
cases in which an initial perturbation is imposed and when the instability
develops spontaneously. However, significant differences in the dynamics of the
small scales are discerned between the NLSE and gyrokinetic simulations. These
discrepancies are mainly due to the radial dependence of the strength of the
nonlinear term, which we do not retain in the solution of the NLSE, and to the
damping of higher spectral components. The damping of the high-$k_r$
components, which develop as a consequence of the nonlinearity, can be
understood in terms of Landau damping. The influence of the ion Larmor radius
$\rho_i$ as well as the perturbation wavevector $k_{\text{pert}}$ on this
effect is studied. For the parameters considered here the aforementioned
damping mechanism hinders the MI process significantly from developing to its
full extent and is strong enough to stabilize some of the (according to the
undamped NLSE model) unstable wavevectors.
