  The orbits of two individual planets in two known binary star systems, \gamma
Cephei and HD 196885 are numerically integrated using various numerical
techniques to assess the chaotic or quasi-periodic nature of the dynamical
system considered. The Hill stability (HS) function which measures the orbital
perturbation of a planet around the primary star due to the secondary star is
calculated for each system. The maximum Lyapunov exponent (MLE) time series are
generated to measure the divergence/convergence rate of stable manifolds, which
are used to differentiate between chaotic and non-chaotic orbits. Then, we
calculate dynamical Mean Exponential Growth factor of Nearby Orbits (MEGNO)
maps from solving the variational equations along with the equations of motion.
These maps allow us to accurately differentiate between stable and unstable
dynamical systems. The results obtained from the analysis of HS, MLE, and MEGNO
maps are analysed for their dynamical variations and resemblance. The HS test
for the planets shows stability and quasi-periodicity for at least ten million
years. The MLE and the MEGNO maps have also indicated the local
quasi-periodicity and global stability in a relatively short integration
period. The orbital stability of the systems is tested using each indicator for
various values of planet inclinations (i_{pl} \le 25^\circ) and binary
eccentricities. The reliability of HS criterion is also discussed based on its
stability results compared with the MLE and MEGNO maps.
