Attempt to improve some codes in `QuadForm` #1254

StevellM · 2023-10-20T12:45:20Z

Mainly trying to avoid creation of lists when it is not necessary and use more the the ! functions whenever possible.

I have tried to improve some type stability by giving new names to variables which were changing types during the sub-processes and inferring more type for lists. Hopefully this with help for compilation time (?).

I have added some semi-colons and spaces in some places where I noticed then missing, for code readability: this was not my first goal so the code is not perfectly reformatted.

src/QuadForm/Herm/Lattices.jl

thofma · 2023-10-20T13:49:11Z

src/QuadForm/Herm/Lattices.jl

@@ -251,7 +251,7 @@ function scale(L::HermLat)
  for i in 1:d
    push!(to_sum, involution(L)(to_sum[i]))
  end
-  s = sum(to_sum, init = zero(base_field(L))*base_ring(L))
+  s = sum(to_sum; init = 0*base_ring(L))


thofma · 2023-10-20T13:49:31Z

src/QuadForm/Herm/Lattices.jl

@@ -359,7 +359,7 @@ function jordan_decomposition(L::HermLat, p)

  k = 1
  while k <= n
-    G = S * F * transpose(_map(S, aut))
+    G = S * F * transpose(map(aut, S))


What did _map do?

The same, but I run several examples using both functions and the underscore one seems to be slightly faster, allocations are the same... I will revert those changed. I first thought that this was a work-around before map was fully implemented.

thofma · 2023-10-20T13:50:49Z

src/QuadForm/Lattices.jl

@@ -305,8 +305,8 @@ function generators(L::AbstractLat; minimal::Bool = false)
    d = ncols(St)
    for i in 1:nrows(St)
      if base_ring(L) isa NfOrd
-        I = numerator(St.coeffs[i])
-        den = denominator(St.coeffs[i])
+        I = numerator(coefficient_ideals(St)[i])


Please check that this does not make a copy

No copy, access directly St.coeffs

thofma · 2023-10-20T13:54:17Z

src/QuadForm/Quad/GenusRep.jl

@@ -155,7 +155,7 @@ function spinor_genera_in_genus(L, mod_out)
  # The smaller the element, the better
  for d in diagonal(Gr)
    if (iszero(spinornorm) && !iszero(d)) || (!iszero(d) && abs(norm(d)) < abs(norm(spinornorm)))
-      spinornorm = d
+      add!(spinornorm, spinornorm, d)


Why add? Wasn't there a = before?

spinornorm is already initialized to be zero, I thought that doing so would lighten the allocation process (?)
add! will just modify the content of spinornorm

No. Nothing is gained from this. It is actually worse with the change. Before it was just an assignment, which is a noop (no operation; does not cost anything). Now there is an addition, which is not for free.

I see, I will correct it.

thofma · 2023-10-20T13:54:32Z

src/QuadForm/Quad/GenusRep.jl

@@ -167,7 +167,7 @@ function spinor_genera_in_genus(L, mod_out)
      end
    end
    @assert !iszero(Gr[1, i])
-    spinornorm = 2 * Gr[1, i]
+    add!(spinornorm, spinornorm, 2 * Gr[1, i])


thofma · 2023-10-20T13:55:28Z

src/QuadForm/Quad/GenusRep.jl

@@ -511,7 +511,7 @@ function maximal_norm_splitting(L, p)
    end
  end

-  @assert all(let G = G; k -> nrows(G[k]) in [1,2]; end, 1:length(G))
+  @assert all(k -> nrows(G[k]) in [1,2], 1:length(G))


Did you check with code_warntype? This is added to work around a julia bug.

I ran the code_warntype and notice wrong appears.. The assertion works also (like if I change some values it returns an error, and does nothing otherwise). So it seems to be fixed (?)

thofma · 2023-10-20T13:55:48Z

src/QuadForm/Quad/GenusRep.jl

+  @hassert :GenRep 1 all(k -> isone(quo(R, factors[k])[2](FacElem(s * the_idele[k]))), 1:length(the_idele))
+  @hassert :GenRep 1 all(j -> sign(s * the_idele_inf[j], IP[j]) == 1, 1:length(IP))


thofma · 2023-10-20T13:56:11Z

src/QuadForm/Quad/Lattices.jl

@@ -232,8 +232,7 @@ function scale(L::QuadLat)
  end
  G = gram_matrix_of_rational_span(L)
  C = coefficient_ideals(L)
-  to_sum = [ G[i, j] * C[i] * involution(L)(C[j]) for j in 1:length(C) for i in 1:j]
-  s = sum(to_sum)
+  s = sum(G[i, j] * C[i] * involution(L)(C[j]) for j in 1:length(C) for i in 1:j)


Needs an init

thofma · 2023-10-20T13:57:41Z

src/QuadForm/Quad/Spaces.jl

@@ -613,13 +612,13 @@ end
 function _quadratic_form_with_invariants(dim::Int, det::nf_elem, finite::Vector, negative::Dict{<:InfPlc, Int})
  @hassert :Lattice 1 dim >= 1
  @hassert :Lattice 1 !iszero(det)
-  K::AnticNumberField = parent(det)
+  K = parent(det)


This was added for type stability. Have you checked that everything is still fine?

I will have to re-put it. There are very odd things happening in this function, which I barely understand on the code_warntype. The variable finite0 is defined and used only at the very end of the function in a hassert. If we remove the hassert and infer the type for finite0, for instance

local finite0::Vector{<: NfOrdIdl} finite0 = copy(finite)

or

local finite0::typeof(finite) finite0 = copy(finite)

or

finite0 = eltype(finite)[p for p in finite]

in all cases, if we remove the only line calling this variable, then code_warntype will still show finite0::Any in red...

This is all related to JuliaLang/julia#15276. Every time one uses a p -> begin; x[1] + y[1] + z + p; end, every "captured" variable which is used (like x, y, z), might become unstable in the surrounding function. This also effects generators like (p + z for p in 1:10). The let is a workaround, also mentioned in the issue above.

I see, I completely forgot about this. So maybe the code_warntype might not capture the issue. Since it was a recent fix I will put this back too until we make sure the problem is solved (to avoid any performance issue)

thofma · 2023-10-20T13:58:29Z

src/QuadForm/Quad/ZGenus.jl

@@ -135,20 +134,19 @@ function _p_adic_symbol(A::MatElem, p, val)
  q = p^m0
  n = nrows(A)
  A = divexact(A, q)
-  Fp = Native.GF(p)
+  Fp = GF(p)


No, should be Native.GF

I will revert them. What is the reason for this ?

thofma · 2023-10-20T13:58:50Z

src/QuadForm/Quad/ZGenus.jl

  q = ZZ(2)^m0
  A = divexact(A, q)
-  A_2 = change_base_ring(Native.GF(2), A)
+  A_2 = change_base_ring(GF(2), A)


thofma · 2023-10-20T14:01:55Z

src/QuadForm/Quad/ZGenus.jl

@@ -1693,9 +1695,9 @@ function _gram_from_jordan_block(p::ZZRingElem, block, discr_form=false)
      q = zero_matrix(QQ, 0, 0)
    end
    if discr_form
-        q = q * (1//2)^level
+      map_entries!(x -> x * (1//2)^level, q, q)


Does that work and yield the correct result?

It works, one has to take care of the type but if everything is made of QQFieldElem, then this overwrites the content of q; the tests show that nothing goes wrong.

thofma · 2023-10-20T14:02:58Z

src/QuadForm/Quad/ZGenus.jl

@@ -2148,11 +2150,11 @@ function _M_p(species, p)
  p = QQ(p)


_p etc while you are at it

thofma · 2023-10-20T14:04:36Z

src/QuadForm/Quad/ZGenus.jl

-      push!(primes, p )
-    end
-  end
+  primes = union!(Hecke.primes(G1), Hecke.primes(G2))


Does primes(G1) make a copy? If not this is wrong

thofma · 2023-10-20T14:05:13Z

src/QuadForm/Quad/ZLattices.jl

@@ -437,8 +439,8 @@ function assert_has_automorphisms(L::ZZLat; redo::Bool = false,
  end

  # Now gens are with respect to the basis of L
-  @hassert :Lattice 1 all(let gens = gens; i -> change_base_ring(FlintQQ, gens[i]) * GL *


thofma · 2023-10-20T14:05:33Z

src/QuadForm/Quad/ZLattices.jl

@@ -1198,7 +1201,7 @@ function is_maximal_even(L::ZZLat, p)
    return true, L
  end
  G = change_base_ring(ZZ, gram_matrix(L))
-  k = Native.GF(p)
+  k = GF(p)


thofma · 2023-10-20T14:07:14Z

src/QuadForm/Torsion.jl

@@ -1501,7 +1500,7 @@ Return the radical `\{x \in T | b(x,T) = 0 and q(x)=0\}` of the quadratic form
 function radical_quadratic(T::TorQuadModule)
  Kb, ib = radical_bilinear(T)
  G = gram_matrix_quadratic(Kb)*1//modulus_bilinear_form(Kb)
-  F = Native.GF(2, cached=false)
+  F = GF(2; cached=false)


thofma

Most of this looks innocent. Most of the let statements came from painstakingly working through issues with @code_warntype. So please make sure that code stability is not degraded. Same is true for sum(iterator) vs sum(Bla[x for x in iterator]). There once were some compiler bugs. Not sure if they are around anymore.

StevellM · 2023-10-20T14:16:13Z

I will check all that, yes, thanks for all the comment.

codecov · 2023-10-23T09:22:45Z

Codecov Report

Merging #1254 (5973970) into master (d8aa190) will decrease coverage by 0.14%.
Report is 8 commits behind head on master.
The diff coverage is 91.51%.

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #1254      +/-   ##
==========================================
- Coverage   74.55%   74.42%   -0.14%     
==========================================
  Files         346      346              
  Lines      110993   111276     +283     
==========================================
+ Hits        82756    82813      +57     
- Misses      28237    28463     +226

Files	Coverage Δ
src/AlgAss/AlgAss.jl	`86.66% <100.00%> (-0.59%)`	⬇️
src/AlgAss/Elem.jl	`89.89% <ø> (ø)`
src/EllCrv/EllCrv.jl	`95.57% <100.00%> (ø)`
src/EllCrv/Finite.jl	`94.28% <100.00%> (-0.17%)`	⬇️
src/EllCrv/Heights.jl	`92.43% <100.00%> (ø)`
src/EllCrv/Isogeny.jl	`90.68% <100.00%> (ø)`
src/EllCrv/Isomorphisms.jl	`90.85% <100.00%> (ø)`
src/FunField/Factor.jl	`43.79% <100.00%> (ø)`
src/GenOrd/Types.jl	`88.78% <100.00%> (ø)`
src/HypellCrv/HypellCrv.jl	`82.70% <100.00%> (ø)`
... and 54 more

... and 44 files with indirect coverage changes

thofma · 2023-10-26T08:14:25Z

Not sure if you were finished with the revision?

StevellM · 2023-10-26T08:53:47Z

Sorry, I forgot to mention it - I am done with the revision.

thofma · 2023-10-27T05:44:01Z

src/QuadForm/Herm/Genus.jl

        end
      end

-      for d in cartesian_product_iterator(dets, inplace = true)# Iterators.product(dets...)


inplace = true is the default for the function so I do not see why keeping it

Sorry, I forgot about that, you are right.

thofma · 2023-10-27T05:44:20Z

src/QuadForm/Herm/Genus.jl

@@ -1748,7 +1748,7 @@ function hermitian_local_genera(E, p, rank::Int, det_val::Int, min_scale::Int, m
      end
    end

-    for dn in cartesian_product_iterator(det_norms, inplace = true)


Same here, why?

thofma · 2023-10-27T05:46:44Z

src/QuadForm/Herm/Mass.jl

@@ -292,7 +292,7 @@ function _standard_mass(L::HermLat, prec::Int = 10)
  local relzeta::arb

  while true
-    relzeta = prod(arb[_L_function(E, 1 - i, wprec) for i in 1:2:m])
+    relzeta = prod(_L_function(E, 1 - i, wprec) for i in 1:2:m; init = 1)


Suggested change

relzeta = prod(_L_function(E, 1 - i, wprec) for i in 1:2:m; init = 1)

relzeta = prod(_L_function(E, 1 - i, wprec) for i in 1:2:m; init = one(ArbField(wprec, cached = false)))

Ah thanks, that was the thing I was trying to get but couldn't make it work.

thofma · 2023-10-27T05:52:40Z

I added a few things. There are still some map instead of _map. Do you want to change it?

StevellM · 2023-10-27T07:02:59Z

Yes I guess so, I just missed them in my last commit. Since we have to make a new commit for your suggestion, I will put everything together.

thofma · 2023-10-27T19:24:01Z

71 min for the long tests is pretty quick. Not sure it is due to the changes here.

StevellM added 3 commits October 20, 2023 14:25

hopefully some improvements

8b4a071

missing bit

286e078

fix docs

efd98dd

thofma reviewed Oct 20, 2023

View reviewed changes

src/QuadForm/Herm/Lattices.jl Outdated Show resolved Hide resolved

thofma reviewed Oct 20, 2023

View reviewed changes

apply suggestions: revert some changes, add more init

aeef2ff

thofma reviewed Oct 27, 2023

View reviewed changes

apply further suggestions

5973970

thofma merged commit 787d345 into thofma:master Oct 27, 2023
15 of 18 checks passed

This pull request was closed.

		@hassert :GenRep 1 all(k -> isone(quo(R, factors[k])[2](FacElem(s * the_idele[k]))), 1:length(the_idele))
		@hassert :GenRep 1 all(j -> sign(s * the_idele_inf[j], IP[j]) == 1, 1:length(IP))

	relzeta = prod(_L_function(E, 1 - i, wprec) for i in 1:2:m; init = 1)
	relzeta = prod(_L_function(E, 1 - i, wprec) for i in 1:2:m; init = one(ArbField(wprec, cached = false)))

Attempt to improve some codes in QuadForm #1254

Attempt to improve some codes in QuadForm #1254

Conversation

StevellM commented Oct 20, 2023

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

thofma Oct 23, 2023 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

thofma left a comment

Choose a reason for hiding this comment

StevellM commented Oct 20, 2023

codecov bot commented Oct 23, 2023 • edited Loading

Codecov Report

thofma commented Oct 26, 2023

StevellM commented Oct 26, 2023

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

thofma commented Oct 27, 2023

StevellM commented Oct 27, 2023

thofma commented Oct 27, 2023

Attempt to improve some codes in `QuadForm` #1254

Attempt to improve some codes in `QuadForm` #1254

thofma Oct 23, 2023 •

edited

Loading

codecov bot commented Oct 23, 2023 •

edited

Loading