General relativity is a fully conservative theory, but there exist other possible metric theories of gravity. We consider non-conservative ones with a parameterized post-Newtonian (PPN) parameter, $\zeta_2$. A non-zero $\zeta_2$ induces a self-acceleration for the center of mass of an eccentric binary pulsar system, which contributes to the second time derivative of the pulsar spin frequency, $\ddot{\nu}$. In our work, using the method in Will (1992), we provide an improved analysis with four well-timed, carefully-chosen binary pulsars. In addition, we extend Will's method and derive $\zeta_2$'s effect on the third time derivative of the spin frequency, $\dddot{\nu}$. For PSR B1913+16, the constraint from $\dddot{\nu}$ is even tighter than that from $\ddot{\nu}$. We combine multiple pulsars with Bayesian inference, and obtain an upper limit, $\left|\zeta_{2}\right|<1.3\times10^{-5}$ at 95% confidence level, assuming a flat prior in $\log_{10} \left| \zeta_{2}\right|$. It improves the existing bound by a factor of three. Moreover, we propose an analytical timing formalism for $\zeta_2$. Our simulated times of arrival with simplified assumptions show binary pulsars' capability in limiting $\zeta_{2}$, and useful clues are extracted for real data analysis in future. In particular, we discover that for PSRs B1913+16 and J0737$-$3039A, $\dddot{\nu}$ can yield more constraining limits than $\ddot{\nu}$.

In Einstein's general relativity, gravity is mediated by a massless spin-2 metric field, and its extension to include a mass for the graviton has profound implication for gravitation and cosmology. In 2002, Finn and Sutton used the gravitational-wave (GW) back-reaction in binary pulsars, and provided the first bound on the mass of graviton. Here we provide an improved analysis using 9 well-timed binary pulsars with a phenomenological treatment. First, individual mass bounds from each pulsar are obtained in the frequentist approach with the help of an ordering principle. The best upper limit on the graviton mass, $m_{g}<3.5\times10^{-20} \, {\rm eV}/c^{2}$ (90% C.L.), comes from the Hulse-Taylor pulsar PSR B1913+16. Then, we combine individual pulsars using the Bayesian theorem, and get $m_{g}<5.2\times10^{-21} \, {\rm eV}/c^{2}$ (90% C.L.) with a uniform prior for $\ln m_g$. This limit improves the Finn-Sutton limit by a factor of more than 10. Though it is not as tight as those from GWs and the Solar System, it provides an independent and complementary bound from a dynamic regime.