diff --git a/src/sage/rings/number_field/order.py b/src/sage/rings/number_field/order.py index 7e2b95e3db5..025c8bd7a28 100644 --- a/src/sage/rings/number_field/order.py +++ b/src/sage/rings/number_field/order.py @@ -120,6 +120,82 @@ def quadratic_order_class_number(disc): return ZZ(h) +def quadratic_order_approximate_class_number(disc, *, bound=10**4): + r""" + Return *an approximation of* the class number of + the quadratic order of given discriminant. + + Currently only implemented for maximal orders + in imaginary-quadratic fields. + + EXAMPLES:: + + sage: from sage.rings.number_field.order import quadratic_order_approximate_class_number + sage: QuadraticField(-419).class_number() + 9 + sage: quadratic_order_approximate_class_number(-419) # rel tol .01 + 9.01653836091712 + + :: + + sage: from sage.rings.number_field.order import quadratic_order_approximate_class_number + sage: d = 100000000000031 + sage: QuadraticField(-d).class_number(proof=False) + 14414435 + sage: round(quadratic_order_approximate_class_number(-d)) # rel tol .01 + 14407657 + sage: round(quadratic_order_approximate_class_number(-d, bound=10**6)) # rel tol .01 + 14413626 + + Test it against the exact class number computed for the CSIDH-512 prime (source: https://eprint.iacr.org/2019/498.pdf):: + + sage: from sage.rings.number_field.order import quadratic_order_approximate_class_number + sage: p = 4 * prod(primes(3,374)) * 587 - 1 + sage: hreal = 84884147409828091725676728670213067387206838101828807864190286991865870575397 + sage: assert not hreal * BQFClassGroup(-p).random_element() + sage: h = round(quadratic_order_approximate_class_number(-p, bound=10**3)); h # rel tol .01 + 85020334529027955331134025285584937708277525762952749243936367281020276898911 + sage: RR(h / hreal) # abs tol .01 + 1.00160438813790 + sage: h = round(quadratic_order_approximate_class_number(-p, bound=10**4)); h # rel tol .01 + 84396416322932187013685015858232028818143153418602750494380687212063263982976 + sage: RR(h / hreal) # abs tol .01 + 0.994254155790231 + sage: h = round(quadratic_order_approximate_class_number(-p, bound=10**5)); h # rel tol .01 + 84823787383264935642168065590216697209124465727576867303448690101681967969348 + sage: RR(h / hreal) # abs tol .01 + 0.999288912848806 + sage: h = round(quadratic_order_approximate_class_number(-p, bound=10**6)); h # rel tol .01, long time (2s) + 84884627342209070883738394179700676127184729325091471467872632417652836154863 + sage: RR(h / hreal) # abs tol .01, long time (2s) + 1.00000565396951 + + ALGORITHM: Finite approximation of the infinite product given by + the analytic class number formula, using primes up to ``bound``. + """ + disc = ZZ(disc) + if disc >= 0: + raise NotImplementedError('only imaginary-quadratic fields supported') + if not disc.is_fundamental_discriminant(): + raise NotImplementedError('only fundamental discriminants supported') + + from sage.rings.real_mpfr import RealField + from sage.arith.misc import primes, kronecker_symbol + from sage.symbolic.constants import pi + + w = 6 if disc == -3 else 4 if disc == -4 else 2 + RR = RealField(max(53, disc.bit_length())) # wild guess! + + # compute numerator and denominator separately for speed + L1 = L2 = RR(1) + for ell in primes(bound): + L1 *= ell + L2 *= ell - kronecker_symbol(disc, ell) + L = L1 / L2 + + return RR(w * abs(disc).sqrt() * L / (2 * pi)) + + class OrderFactory(UniqueFactory): r""" Abstract base class for factories creating orders, such as @@ -1091,6 +1167,15 @@ def class_number(self, proof=None): r""" Return the class number of this order. + .. NOTE:: + + For some applications (e.g., in algorithms for computing + class groups) it is required to merely *approximate* the + class number. The function + :func:`quadratic_order_approximate_class_number` + can be used to compute such an approximation (currently + restricted to maximal imaginary-quadratic orders). + EXAMPLES:: sage: ZZ[2^(1/3)].class_number() # needs sage.symbolic