Remove outdated linear solving code (#1413)

examples/Jordan.jl

@@ -32,7 +32,7 @@ function rational_normal_form(M::T) where {T <: MatElem{S} where {S <: FieldElem
     m = (g^v)(M)
     this_g = zero_matrix(base_ring(M), nrows(M), 0)
     while true
-      d, km = nullspace(m)
+      km = kernel(m, side = :right)
       km = quo_gens(km, this_g)
       if ncols(km) == 0
         break

examples/MultDep.jl

@@ -234,16 +234,16 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN
           end
         end
       end
-      @vtime :qAdic 1 k = Hecke._left_kernel_basis(lv)
-      @assert length(k) < 2
-      if length(k) == 0
+      @vtime :qAdic 1 k = kernel(lv, side = :left)
+      @assert nrows(k) < 2
+      if nrows(k) == 0
         println("new ")
         push!(u, a)
         lu = vcat(lu, la)
         @assert length(u) <= r
       else # length == 1 extend the module
         s = QQFieldElem[]
-        for x = k[1]
+        for x in k[1, :]
           @vtime :qAdic 1 y = lift_reco(FlintQQ, x, reco = true)
           if y === nothing
             prec *= 2
@@ -253,7 +253,7 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN
           end
           push!(s, y)
         end
-        if length(s) < length(k[1])
+        if length(s) < ncols(k)
           continue
         end
         d = reduce(lcm, map(denominator, s))

examples/RelNeq.jl

@@ -123,7 +123,8 @@ function norm_1_subgroup(A::RelNeq)
       end
       f = [ZZRingElem(div(degree(Q[1]), degree(P[1]))) for Q = lP]
       m = matrix(FlintZZ, 1, length(f), f)
-      r, n = nullspace(m)
+      n = kernel(m, side = :right)
+      r = ncols(n)
 
       decom = [mq(mr\Q[1]) for Q = lP]
       for i=1:r
@@ -387,7 +388,8 @@ function n1group(A::RelNeq, B::Int)
       lq = Hecke.prime_decomposition_nonindex(A.m_k_K, P)
     end
     f = matrix(FlintZZ, 1, length(lq), ZZRingElem[div(degree(Q[1]), degree(P)) for Q = lq])
-    r, n = nullspace(f)
+    n = kernel(f, side = :right)
+    r = ncols(n)
     res = false
     for i = 1:r
       I = evaluate(FacElem(Dict((lq[j][1], n[j,i]) for j = 1:length(lq))), coprime = true)

src/AlgAss/AbsAlgAss.jl

@@ -174,8 +174,8 @@ function kernel_of_frobenius(A::AbstractAssociativeAlgebra)
     end
   end
 
-  V = _right_kernel_basis(B)
-  return [ A(v) for v in V ]
+  V = kernel(B, side = :right)
+  return [ A(V[:, i]) for i in 1:ncols(V) ]
 end
 
 @doc raw"""

src/AlgAss/AlgAss.jl

@@ -899,7 +899,6 @@ function _build_subalgebra_mult_table!(A::StructureConstantAlgebra{T}, B::MatEle
   r = rref!(B)
   if r == 0
     if return_LU == Val{true}
-      #return Array{elem_type(K), 3}(undef, 0, 0, 0),  SymmetricGroup(ncols(B))(), zero_matrix(K, 0, 0), zero_matrix(K, 0, 0), LinearSolveCtx{typeof(B)}
       return Array{elem_type(K), 3}(undef, 0, 0, 0), solve_init(B)
     else
       return Array{elem_type(K), 3}(undef, 0, 0, 0)
@@ -1132,7 +1131,8 @@ function center(A::StructureConstantAlgebra{T}) where {T}
   M = zero_matrix(base_ring(A), n^2, n)
   # I concatenate the difference between the right and left representation matrices.
   _rep_for_center!(M, A)
-  k, B = nullspace(M)
+  B = kernel(M, side = :right)
+  k = ncols(B)
   res = Vector{elem_type(A)}(undef, k)
   for i=1:k
     res[i]= A(T[B[j,i] for j=1:n])
@@ -1338,7 +1338,7 @@ function _matrix_basis(A::StructureConstantAlgebra{T}, idempotents::Vector{S}) w
     M4 = representation_matrix(bb - one(eAe), :right)
 
     M = hcat(M1, M2, M3, M4)
-    xx = eAe(_left_kernel_basis(M)[1])
+    xx = eAe(kernel(M, side = :left)[1, :])
     x = m1(m2(xx))
 
     new_basis[1 + (i - 1)*k] = x # this is e_1i

src/AlgAss/AlgMat.jl

@@ -144,19 +144,17 @@ function basis_matrix_rref(A::MatAlgebra{S, T}) where {S, T}
   return A.basis_matrix_rref::Tuple{dense_matrix_type(S), dense_matrix_type(S), Vector{Int}}
 end
 
-function assure_has_basis_matrix_transpose_rref(A::MatAlgebra)
-  if isdefined(A, :basis_matrix_transpose_rref)
+function assure_has_basis_matrix_solve_context(A::MatAlgebra)
+  if isdefined(A, :basis_matrix_solve_ctx)
     return nothing
   end
-  s, R, U = _rref_with_trans(transpose(basis_matrix(A, copy = false)))
-  pivots = _get_pivots_ut(R)
-  A.basis_matrix_transpose_rref = (R, U, pivots)
+  A.basis_matrix_solve_ctx = solve_init(basis_matrix(A))
   return nothing
 end
 
-function basis_matrix_transpose_rref(A::MatAlgebra{S, T}) where {S, T}
-  assure_has_basis_matrix_transpose_rref(A)
-  return A.basis_matrix_transpose_rref::Tuple{dense_matrix_type(S), dense_matrix_type(S), Vector{Int}}
+function basis_matrix_solve_context(A::MatAlgebra{S, T}) where {S, T}
+  assure_has_basis_matrix_solve_context(A)
+  return A.basis_matrix_solve_ctx
 end
 
 ################################################################################
@@ -580,13 +578,8 @@ function _check_matrix_in_algebra(M::S, A::MatAlgebra{T, S}, short::Type{Val{U}}
       end
     end
   end
-  R, UU, piv = basis_matrix_transpose_rref(A)
-  #@show B, tt
-  #@show UU
-  #@show UU * B == R
-  b, N = _can_solve_given_rref(R, UU, piv, t)
-  #b, N = can_solve_with_solution(B, tt, side = :left)
-  #@assert b == bb && NN == [N[1, i] for i in 1:length(N)]
+  C = basis_matrix_solve_context(A)
+  b, N = can_solve_with_solution(C, t, side = :left)
   if short == Val{true}
     return b
   end

src/AlgAss/Types.jl

@@ -282,7 +282,7 @@ end
   basis
   basis_matrix # matrix over the base_ring
   basis_matrix_rref
-  basis_matrix_transpose_rref
+  basis_matrix_solve_ctx
   dim::Int
   degree::Int
   is_simple::Int

src/AlgAss/radical.jl

@@ -54,7 +54,8 @@ end
 
 function _radical_zero(A::AbstractAssociativeAlgebra{T}) where { T <: Union{ QQFieldElem, NumFieldElem } }
   M = trace_matrix(A)
-  n, N = nullspace(M)
+  N = kernel(M; side = :right)
+  n = ncols(N)
   b = Vector{elem_type(A)}(undef, n)
   t = zeros(base_ring(A), dim(A))
   # the construct A(t) will make a copy (hopefully :))
@@ -150,7 +151,8 @@ function _radical_finite_prime_field(A::MatAlgebra{<:Union{fpFieldElem, FpFieldE
   # We have to take the matrix trace : M_n(K) -> K
   # and not the "algebra" trace
   MF = trace_matrix(A, x -> tr(matrix(x, copy = false)))
-  d, B = nullspace(MF)
+  B = kernel(MF; side = :right)
+  d = ncols(B)
   if d == 0
     return elem_type(A)[]
   end
@@ -182,7 +184,8 @@ function _radical_finite_prime_field(A::MatAlgebra{<:Union{fpFieldElem, FpFieldE
           M[j, i] = F(divexact(t, pl))
         end
       end
-      d, B = nullspace(M)
+      B = kernel(M; side = :right)
+      d = ncols(B)
       if d == 0
         return elem_type(A)[]
       end
@@ -200,7 +203,8 @@ function _radical_finite_prime_field(A::AbstractAssociativeAlgebra)
   k = flog(ZZRingElem(dim(A)), p)
 
   MF = trace_matrix(A)
-  d, B = nullspace(MF)
+  B = kernel(MF; side = :right)
+  d = ncols(B)
   if d == 0
     return elem_type(A)[]
   end
@@ -230,8 +234,8 @@ function _radical_finite_prime_field(A::AbstractAssociativeAlgebra)
         M[j, i] = F(divexact(t, pl))
       end
     end
-    d, B = nullspace(M)
-    if d == 0
+    B = kernel(M; side = :right)
+    if ncols(B) == 0
       return elem_type(A)[]
     end
     C = transpose(B)*C
@@ -282,8 +286,8 @@ function _radical_finite_generic(A::AbstractAssociativeAlgebra{T}) where {T <: U
   K, a = number_field(f, "a")
 
   MF = trace_matrix(A2)
-  d, B = nullspace(MF)
-  if d == 0
+  B = kernel(MF; side = :right)
+  if ncols(B) == 0
     return elem_type(A)[]
   end
 
@@ -309,8 +313,8 @@ function _radical_finite_generic(A::AbstractAssociativeAlgebra{T}) where {T <: U
         end
       end
     end
-    d, B = nullspace(M)
-    if d == 0
+    B = kernel(M; side = :right)
+    if ncols(B) == 0
       return elem_type(A)[]
     end
     C = transpose(B)*C

src/AlgAssAbsOrd/ICM.jl

@@ -259,7 +259,7 @@ function matrix_to_ideal(O::AlgAssAbsOrd, M::ZZMatrix)
   result = zeros(A, dim(A))
   for (K, AtoK) in fields_and_maps
     MK = change_base_ring(K, M) - gen(K)*identity_matrix(K, dim(A))
-    _, B = nullspace(MK)
+    B = kernel(MK, side = :right)
     for j = 1:ncols(B)
       for i = 1:dim(A)
         result[i] += AtoK\B[i, j]
@@ -275,7 +275,7 @@ function matrix_to_ideal(O::AbsNumFieldOrder, M::ZZMatrix)
   @assert K.pol == parent(K.pol)(f)
   result = zeros(K, degree(K))
   MK = change_base_ring(K, M) - gen(K)*identity_matrix(K, degree(K))
-  _, B = nullspace(MK)
+  B = kernel(MK, side = :right)
   for j = 1:ncols(B)
     for i = 1:degree(K)
       result[i] += B[i, j]

src/AlgAssAbsOrd/Ideal.jl

@@ -1178,7 +1178,7 @@ function pradical_meataxe(O::AlgAssAbsOrd, p::Int)
     return J
   end
   M1 = view(M1, 1:r, 1:degree(O))
-  dM = transpose(nullspace(M1)[2])
+  dM = transpose(kernel(M1, side = :right))
   g = Vector{elem_type(algebra(O))}(undef, nrows(dM) + 1)
   m = zero_matrix(FlintZZ, degree(O), degree(O))
   for i = 1:nrows(dM)
@@ -1210,7 +1210,8 @@ function pradical(O::AlgAssAbsOrd, p::IntegerUnion)
   F = GF(p, cached = false)
 
   I = change_base_ring(F, trred_matrix(O))
-  k, B = nullspace(I)
+  B = kernel(I, side = :right)
+  k = ncols(B)
   # The columns of B give the coordinates of the elements in the order.
   if k == 0
     J = ideal(algebra(O), O, p*basis_matrix(FakeFmpqMat, O, copy = false); side=:twosided)

src/GrpAb/Map.jl

@@ -143,7 +143,8 @@ function hom(A::Vector{FinGenAbGroupElem}, B::Vector{FinGenAbGroupElem}; check::
   if (check)
     m = vcat(ZZMatrix[x.coeff for x in A])
     m = vcat(m, rels(parent(A[1])))
-    i, T = nullspace(transpose(m))
+    T = kernel(transpose(m), side = :right)
+    i = ncols(T)
     T = transpose(T)
     T = sub(T, 1:nrows(T), 1:length(A))
     n = vcat(ZZMatrix[x.coeff for x in B])

src/GrpAb/stable_sub.jl

@@ -152,12 +152,12 @@ end
 
 function _dualize_1(M::zzModMatrix, snf_struct::Vector{ZZRingElem})
 
-  A=nullspace(transpose(M))
+  A=kernel(transpose(M), side = :right)
   B=vcat(transpose(A),zero_matrix(M[1,1].parent, ncols(A),ncols(A)))
   for j=1:ncols(A)
     B[nrows(A)+j,j]=snf_struct[j]
   end
-  S=nullspace(B)
+  S=kernel(B, side = :right)
   C=vcat(transpose(A),zero_matrix(M[1,1].parent, ncols(A),ncols(A)))
   return S*C
 

src/LocalField/map.jl

@@ -10,8 +10,7 @@ mutable struct LocalFieldMor{S, T, U, V, W} <: Map{S, T, HeckeMap, LocalFieldMor
   inverse_data::V
   absolute_basis::Vector{W}
   absolute_basis_matrix_image::Generic.MatSpaceElem{PadicFieldElem}
-  rref::Tuple{Generic.MatSpaceElem{PadicFieldElem}, Generic.MatSpaceElem{PadicFieldElem}}
-  pivots_of_rref::Vector{Int}
+  solve_context::Solve.SolveCtx{PadicFieldElem, Generic.MatSpaceElem{PadicFieldElem}, Solve.LazyTransposeMatElem{PadicFieldElem, Generic.MatSpaceElem{PadicFieldElem}}}
 
   function LocalFieldMor{S, T, U, V}() where {S, T, U, V}
     z = new{S, T, U, V, elem_type(S)}()
@@ -349,13 +348,9 @@ function _assert_has_preimage_data(f::LocalFieldMor)
     end
   end
 
-  r, R, U =  _rref_with_trans(M) #terribly unstable numerically
-  pivots = _get_pivots_ut(R)
-
+  f.solve_context = solve_init(M)
   f.absolute_basis_matrix_image = M
   f.absolute_basis = b
-  f.pivots_of_rref = pivots
-  f.rref = R, U
 
   return nothing
 end
@@ -369,7 +364,7 @@ function has_preimage_with_preimage(f::LocalFieldMor, g::Union{QadicFieldElem, L
   d = absolute_degree(K)
   cc = absolute_coordinates(g)
   _assert_has_preimage_data(f)
-  fl, s = _can_solve_given_rref(f.rref[1], f.rref[2], f.pivots_of_rref, cc)
+  fl, s = can_solve_with_solution(f.solve_context, cc, side = :right)
   if !fl
     return false, zero(K)
   else

src/Map/NumField.jl

@@ -109,8 +109,7 @@ mutable struct NumFieldHom{S, T, U, V, W} <: Map{S, T, HeckeMap, NumFieldHom}
   inverse_data::V
   absolute_basis::Vector{W}
   absolute_basis_matrix_image::QQMatrix
-  rref::Tuple{QQMatrix, QQMatrix}
-  pivots_of_rref::Vector{Int}
+  solve_context::Solve.SolveCtx{QQFieldElem, QQMatrix, Solve.LazyTransposeMatElem{QQFieldElem, QQMatrix}}
 
   function NumFieldHom{S, T, U, V}() where {S, T, U, V}
     z = new{S, T, U, V, elem_type(S)}()
@@ -753,13 +752,9 @@ function _assert_has_preimage_data(f::NumFieldHom)
     end
   end
 
-  r, R, U =  _rref_with_trans(M)
-  pivots = _get_pivots_ut(R)
-
+  f.solve_context = solve_init(M)
   f.absolute_basis_matrix_image = M
   f.absolute_basis = b
-  f.pivots_of_rref = pivots
-  f.rref = R, U
 
   return nothing
 end
@@ -772,7 +767,7 @@ function has_preimage_with_preimage(f::NumFieldHom, g::NumFieldElem)
   cc = absolute_coordinates(g)
   K = domain(f)
   _assert_has_preimage_data(f)
-  fl, s = _can_solve_given_rref(f.rref[1], f.rref[2], f.pivots_of_rref, cc)
+  fl, s = can_solve_with_solution(f.solve_context, cc, side = :right)
   if !fl
     return false, zero(K)
   else