Modeling the Evolution and Distribution of the Frequency's Second Derivative and Braking Index of Pulsar Spin with Simulations

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_{\rm 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_{\rm c}$ and $\log|n|-\log\tau_{\rm 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.