fix: make GRH=false more robust

src/AlgAss/AlgQuat.jl

@@ -363,9 +363,9 @@ function unit_group_modulo_scalars(O::AlgAssAbsOrd)
   return enumerate(O, 1)
 end
 
-function _unit_group_generators_quaternion(O::Union{AlgAssRelOrd, AlgAssAbsOrd})
+function _unit_group_generators_quaternion(O::Union{AlgAssRelOrd, AlgAssAbsOrd}; GRH::Bool = true)
   gens1 = unit_group_modulo_scalars(O)
-  u, mu = unit_group(base_ring(O))
+  u, mu = unit_group(base_ring(O); GRH = GRH)
   A = algebra(O)
   gens2 = [ O(A(elem_in_nf(mu(u[i])))) for i in 1:ngens(u) ]
   return append!(gens1, gens2)

src/AlgAssAbsOrd/LocallyFreeClassGroup.jl

@@ -16,7 +16,7 @@
 Given an order $O$ in a semisimple algebra over $\mathbb Q$, this function
 returns the locally free class group of $O$.
 """
-function locally_free_class_group(O::AlgAssAbsOrd, cond::Symbol = :center, ::Val{return_disc_log_data} = Val(false)) where return_disc_log_data
+function locally_free_class_group(O::AlgAssAbsOrd, cond::Symbol = :center, ::Val{return_disc_log_data} = Val(false); GRH::Bool = true) where return_disc_log_data
   A = algebra(O)
   OA = maximal_order(O)
   Z, ZtoA = center(A)
@@ -70,7 +70,7 @@
   end
 
   @vprintln :LocallyFreeClassGroup "Compute ray class group"
-  R, mR = ray_class_group(FinZ, inf_plc)
+  R, mR = ray_class_group(FinZ, inf_plc; GRH = GRH)
 
   # Compute K_1(O/F) and the subgroup of R generated by nr(a)*OZ for a in k1 where
   # nr is the reduced norm and OZ the maximal order in Z
@@ -86,7 +86,7 @@
     while !t
       k += 1
       r = rand(F, 1)
-      x += r
+      x += parent(x)(r)
       t = is_invertible(elem_in_algebra(x, copy = false))[1]
       if k > 100
         error("Something wrong")

src/AlgAssAbsOrd/PIP/bley_johnston.jl

@@ -5,7 +5,7 @@
 ################################################################################
 
 # Compute a set S, such that pi(S) = im(pi), where pi : M^* -> (M/F)^*.
-function _unit_reps(M, F)
+function _unit_reps(M, F; GRH::Bool = true)
   D = get_attribute!(M, :unit_reps) do
     return Dict{typeof(F), Vector{Vector{elem_type(algebra(M))}}}()
   end::Dict{typeof(F), Vector{Vector{elem_type(algebra(M))}}}
@@ -14,20 +14,20 @@ function _unit_reps(M, F)
     @vprintln :PIP "Unit representatives cached for this conductor ideal"
     return D[F]
   else
-    u = __unit_reps(M, F)
+    u = __unit_reps(M, F; GRH = GRH)
     D[F] = u
     return u
   end
 end
 
-function __unit_reps_simple(M, F)
+function __unit_reps_simple(M, F; GRH::Bool = true)
   B = algebra(M)
   @vprintln :PIP _describe(B)
   @vprintln :PIP "Computing generators of the maximal order"
   if dim(B) == 0
     return [zero(B)]
   end
-  UB = _unit_group_generators_maximal_simple(M)
+  UB = _unit_group_generators_maximal_simple(M; GRH = GRH)
   Q, MtoQ = quo(M, F)
   for u in UB
     @assert u in M && inv(u) in M
@@ -137,7 +137,7 @@ function __unit_reps_estimates(M, F)
   #return unit_reps
 end
 
-function __unit_reps(M, F)
+function __unit_reps(M, F; GRH::Bool = true)
   #_assert_has_refined_wedderburn_decomposition(algebra(M))
   A = algebra(M)
   dec = decompose(algebra(M))
@@ -147,7 +147,7 @@ function __unit_reps(M, F)
     FinB = ideal_from_lattice_gens(B, elem_type(B)[(mB\(b)) for b in absolute_basis(F)])
     @assert Hecke._test_ideal_sidedness(FinB, MinB, :right)
     FinB.order = MinB
-    _unit_reps =  __unit_reps_simple(MinB, FinB)
+    _unit_reps =  __unit_reps_simple(MinB, FinB; GRH = GRH)
     @vprintln :PIP "Mapping back once more"
     to_return = Vector{elem_type(A)}(undef, length(_unit_reps))
     Threads.@threads for i in 1:length(_unit_reps)
@@ -158,7 +158,7 @@ function __unit_reps(M, F)
   return unit_reps
 end
 
-function _is_principal_with_data_bj(I, O; side = :right, _alpha = nothing, local_freeness::Bool = false)
+function _is_principal_with_data_bj(I, O; side = :right, _alpha = nothing, local_freeness::Bool = false, GRH::Bool = true)
   # local_freeness needs to be accepted since the generic interface uses it
   A = algebra(O)
   if _alpha === nothing
@@ -241,7 +241,7 @@ function _is_principal_with_data_bj(I, O; side = :right, _alpha = nothing, local
     #@show FinB
   end
 
-  unit_reps = _unit_reps(M, F)
+  unit_reps = _unit_reps(M, F; GRH = GRH)
 
   decc = copy(dec)
   p = sortperm(unit_reps, by = x -> length(x))

src/AlgAssAbsOrd/PIP/unit_group_generators.jl

@@ -11,15 +11,15 @@ function _unit_group_generators(O)
   return Y
 end
 
-function _unit_group_generators_maximal(M)
+function _unit_group_generators_maximal(M; GRH::Bool = true)
   res = decompose(algebra(M))
   Mbas = basis(M)
   idems = [mB(one(B)) for (B, mB) in res]
   gens = elem_type(algebra(M))[]
   for i in 1:length(res)
     B, mB = res[i]
     MinB = Order(B, [(mB\(mB(one(B)) * elem_in_algebra(b))) for b in Mbas])
-    UB = _unit_group_generators_maximal_simple(MinB)
+    UB = _unit_group_generators_maximal_simple(MinB; GRH = GRH)
     e = sum(idems[j] for j in 1:length(res) if j != i; init = zero(algebra(M)))
     @assert isone(e + mB(one(B)))
     for u in UB
@@ -31,19 +31,21 @@ function _unit_group_generators_maximal(M)
   return gens
 end
 
-function _unit_group_generators_maximal_simple(M)
+function _unit_group_generators_maximal_simple(M; GRH::Bool = true)
   A = algebra(M)
   if dim(A) == 0
     return [zero(A)]
   end
   ZA, ZAtoA = _as_algebra_over_center(A)
   if dim(ZA) == 1
+    !GRH && error("Not implemented (yet)")
     # this is a field
     K, AtoK = _as_field_with_isomorphism(A)
     OK = maximal_order(K)
     u, mu = unit_group(OK)
     return [preimage(AtoK, elem_in_nf(mu(u[i]))) for i in 1:ngens(u)]
   elseif isdefined(A, :isomorphic_full_matrix_algebra)
+    !GRH && error("Not implemented (yet)")
     B, AtoB = A.isomorphic_full_matrix_algebra
     OB = _get_order_from_gens(B, [AtoB(elem_in_algebra(b)) for b in absolute_basis(M)])
     N, S = nice_order(OB)
@@ -58,7 +60,7 @@ function _unit_group_generators_maximal_simple(M)
   elseif dim(ZA) == 4 && !is_split(ZA) && !isdefined(A, :isomorphic_full_matrix_algebra)
     Q, QtoZA = is_quaternion_algebra(ZA)
     MQ = _get_order_from_gens(Q, [QtoZA\(ZAtoA\(elem_in_algebra(b))) for b in absolute_basis(M)])
-    _gens =  _unit_group_generators_quaternion(MQ)
+    _gens =  _unit_group_generators_quaternion(MQ; GRH = GRH)
     gens_in_M = [ ZAtoA(QtoZA(elem_in_algebra(u))) for u in _gens]
     @assert all(b in M for b in gens_in_M)
     return gens_in_M

src/AlgAssAbsOrd/PicardGroup.jl

@@ -545,7 +545,7 @@ end
 ################################################################################
 
 # inf_pcl[i] may contain places for the field A.maps_to_numberfields[i][1]
-function ray_class_group(m::AlgAssAbsOrdIdl, inf_plc::Vector{Vector{T}} = Vector{Vector{InfPlc{AbsSimpleNumField, AbsSimpleNumFieldEmbedding}}}()) where {T}
+function ray_class_group(m::AlgAssAbsOrdIdl, inf_plc::Vector{Vector{T}} = Vector{Vector{InfPlc{AbsSimpleNumField, AbsSimpleNumFieldEmbedding}}}(); GRH::Bool = true) where {T}
   O = order(m)
   A = algebra(O)
   fields_and_maps = as_number_fields(A)
@@ -555,9 +555,9 @@ function ray_class_group(m::AlgAssAbsOrdIdl, inf_plc::Vector{Vector{T}} = Vector
   for i = 1:length(fields_and_maps)
     mi = _as_ideal_of_number_field(m, fields_and_maps[i][2])
     if length(inf_plc) != 0
-      r, mr = ray_class_group(mi, inf_plc[i])
+      r, mr = ray_class_group(mi, inf_plc[i]; GRH = GRH)
     else
-      r, mr = ray_class_group(mi)
+      r, mr = ray_class_group(mi; GRH = GRH)
     end
     push!(groups, _make_disc_exp_deterministic(mr))
     #push!(groups, ray_class_group(mi))

src/AlgAssAbsOrd/UnitGroup.jl

@@ -4,14 +4,14 @@
 #
 ################################################################################
 
-function unit_group(O::AlgAssAbsOrd)
+function unit_group(O::AlgAssAbsOrd; GRH::Bool = true)
   @assert is_commutative(O)
   mU = get_attribute!(O, :unit_group) do
     if is_maximal(O)
-      U, mU = _unit_group_maximal(O)
+      U, mU = _unit_group_maximal(O; GRH = GRH)
     else
       OK = maximal_order(O)
-      UU, mUU = unit_group(OK)
+      UU, mUU = unit_group(OK; GRH = GRH)
       U, mU = _unit_group_non_maximal(O, OK, mUU)
     end
     return mU
@@ -101,10 +101,10 @@ function _unit_group_maximal_fac_elem(O::AlgAssAbsOrd)
   return S, StoO
 end
 
-function _unit_group_maximal(O::AlgAssAbsOrd)
+function _unit_group_maximal(O::AlgAssAbsOrd; GRH::Bool = true)
   A = algebra(O)
   fields_and_maps = as_number_fields(A)
-  unit_groups = Tuple{FinGenAbGroup, MapUnitGrp{AbsSimpleNumFieldOrder}}[ unit_group(maximal_order(field)) for (field, map) in fields_and_maps ]
+  unit_groups = Tuple{FinGenAbGroup, MapUnitGrp{AbsSimpleNumFieldOrder}}[ unit_group(maximal_order(field); GRH = GRH) for (field, map) in fields_and_maps ]
   G = unit_groups[1][1]
   for i = 2:length(unit_groups)
     G = direct_product(G, unit_groups[i][1], task = :none)::FinGenAbGroup

src/NumFieldOrd/NfOrd/Clgp.jl

@@ -345,11 +345,32 @@
 
   if !GRH
     class_group_proof(c, ZZRingElem(2), factor_base_bound_minkowski(O))
+    sat = ZZRingElem[]
     for (p, _) in factor(c.h)
+      push!(sat, p)
       while saturate!(c, U, Int(p), 3.5)
       end
     end
     c.GRH = false
+
+    if unit_group_rank(O) > 0
+      # need to make sure that the unit group is correct
+      tent_reg = tentative_regulator(U)
+      low_reg = lower_regulator_bound(Hecke.nf(O))
+      fl, regindexbound = unique_integer(floor(tent_reg/low_reg))
+      @vprintln :ClassGroupProof "Saturating unit group up to $(regindexbound))"
+      for p in PrimesSet(1, Int(regindexbound))
+        if p in sat
+          continue
+        end
+        while saturate!(c, U, Int(p), 3.5)
+        end
+      end
+    end
+  end
+
+  if degree(O) == 1
+    @assert c.h == 1
   end
 
   return c, U, _validate_class_unit_group(c, U)[1]
@@ -470,7 +491,7 @@
 function unit_group_fac_elem(O::AbsSimpleNumFieldOrder; method::Int = 3, unit_method::Int = 1, use_aut::Bool = false, GRH::Bool = true, redo::Bool = false)
   if !is_maximal(O)
     OK = maximal_order(nf(O))
-    UUU, mUUU = unit_group_fac_elem(OK)::Tuple{FinGenAbGroup, MapUnitGrp{FacElemMon{AbsSimpleNumField}}}
+    UUU, mUUU = unit_group_fac_elem(OK; GRH = GRH)::Tuple{FinGenAbGroup, MapUnitGrp{FacElemMon{AbsSimpleNumField}}}
     return _unit_group_non_maximal(O, OK, mUUU)::Tuple{FinGenAbGroup, MapUnitGrp{FacElemMon{AbsSimpleNumField}}}
   end
 

src/NumFieldOrd/NfOrd/Unit/Regulator.jl

@@ -53,5 +53,9 @@ function regulator(x::Vector{T}, abs_tol::Int = 64) where T
 end
 
 function lower_regulator_bound(K::AbsSimpleNumField)
-  return ArbField(64, cached = false)("0.054")
+  R = ArbField(64, cached = false)
+  r1, r2 = signature(K)
+  w = torsion_units_order(K)
+  zimmert = w * 4//100 * exp(46//100 * R(r1) + 1/100 * R(r2))
+  return max(zimmert, ArbField(64, cached = false)("0.054"))
 end