def _R2deriv(self,R,z,phi=0.,t=0.):
        """
        NAME:
           _Rderiv
        PURPOSE:
           evaluate the second radial derivative for this potential
        INPUT:
           R - Galactocentric cylindrical radius
           z - vertical height
           phi - azimuth
           t - time
        OUTPUT:
           the second radial derivative
        HISTORY:
           2018-08-04 - Written - Bovy (UofT)
        """
        Raz2= (R+self.a)**2+z**2
        Raz= nu.sqrt(Raz2)
        m= 4.*R*self.a/Raz2
        R2ma2mz2o4aR1m= (R**2-self.a2-z**2)/4./self.a/R/(1.-m)
        return (2*R**2+self.a2+3*R*self.a+z**2)/R/Raz2*self._Rforce(R,z)\
            +2.*self.a/R/Raz*(m*(R**2+self.a2+z**2)/4./(1.-m)/self.a/R**2\
                                  *special.ellipe(m)\
              +(R2ma2mz2o4aR1m/(1.-m)*special.ellipe(m)
                +0.5*R2ma2mz2o4aR1m*(special.ellipe(m)-special.ellipk(m))
                +0.5*(special.ellipe(m)/(1.-m)-special.ellipk(m))/m)\
                                  *4*self.a*(self.a2+z**2-R**2)/Raz2**2)
NAME:
           _Rderiv
        PURPOSE:
           evaluate the second radial derivative for this potential
        INPUT:
           R - Galactocentric cylindrical radius
           z - vertical height
           phi - azimuth
           t - time
        OUTPUT:
           the second radial derivative
        HISTORY:
           2018-08-04 - Written - Bovy (UofT)