@staticmethod
def density(r, sigma0, r_core):
    """
        rho(r) =  2/pi * Sigma_crit R_c**3 * (R_c**2 + r**2)**(-2)

        :param r: radius (angular scale)
        :param sigma0: convergence in the core
        :param r_core: core radius
        :return: density at radius r
        """
    return 1.0 / 2 * sigma0 * r_core ** 2 * (r_core ** 2 + r ** 2) ** (-3.0 / 2
        )
