We model the evolution of spin frequency ’ s second derivative \ddot { \nu } and braking index n of radio pulsars with simulations within the phenomenological model of their surface magnetic field evolution , which contains a long-term decay modulated by short-term oscillations .
For the pulsar PSR B0329+54 , the model can reproduce the main characteristics of its \ddot { \nu } variation with oscillation periods , predicts another \sim 50 yr oscillation component and another recent swing of the sign of \ddot { \nu } .
We show that the “ averaged ” n is different from the instantaneous n , and its oscillation magnitude decreases abruptly as the time span increases , due to the “ averaging ” effect .
The simulation predicted timing residuals agree with the main features of the reported data .
We further perform Monte Carlo simulations for the distribution of the reported data in | \ddot { \nu } | versus characteristic age \tau _ { c } diagram .
The model with a power law index \alpha = 0.5 can reproduce the slope of the linear fit to pulsars ’ distributions in the diagrams of \log| \ddot { \nu } | - \log \tau _ { c } and \log|n| - \log \tau _ { c } , but the oscillations are responsible for the almost equal number of positive and negative values of \ddot { \nu } , in agreement with our previous analytical studies ; an oscillation period of about several decades is also preferred .
However the range of the oscillation amplitudes is -11.4 \lesssim \log f \lesssim - 10.2 , slightly lager than the analytical prediction , \log f \simeq - 11.85 , because the “ averaging ” effect was not included previously .