Skip to content

Commit

Permalink
Switch to the Solve.* methods (#1671)
Browse files Browse the repository at this point in the history
  • Loading branch information
joschmitt authored Feb 12, 2024
1 parent abfd06b commit 2098cf8
Show file tree
Hide file tree
Showing 30 changed files with 429 additions and 461 deletions.
2 changes: 1 addition & 1 deletion Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ RandomExtensions = "fb686558-2515-59ef-acaa-46db3789a887"
SHA = "ea8e919c-243c-51af-8825-aaa63cd721ce"

[compat]
AbstractAlgebra = "0.37.6"
AbstractAlgebra = "0.38.1"
Antic_jll = "~0.201.500"
Arb_jll = "~200.2300.000"
Calcium_jll = "~0.401.100"
Expand Down
2 changes: 1 addition & 1 deletion benchmarks/solve_polynomials.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ function benchmark_solve_poly()
end
end

tt = @elapsed solve_rational(M, b)
tt = @elapsed _solve_rational(M, b)
println("$tt")
end

23 changes: 0 additions & 23 deletions docs/src/matrix.md
Original file line number Diff line number Diff line change
Expand Up @@ -128,29 +128,6 @@ c = det_divisor(A)
d = det_given_divisor(A, c)
```

### Linear solving

```@docs
cansolve(::ZZMatrix, ::ZZMatrix)
```

```@docs
solve_dixon(::ZZMatrix, ::ZZMatrix)
solve_dixon(::QQMatrix, ::QQMatrix)
```

**Examples**

```julia
S = matrix_space(ZZ, 3, 3)
T = matrix_space(ZZ, 3, 1)

A = S([ZZ(2) 3 5; 1 4 7; 9 2 2])
B = T([ZZ(4), 5, 7])

X, m = solve_dixon(A, B)
```

### Pseudo inverse

```@docs
Expand Down
9 changes: 0 additions & 9 deletions src/Exports.jl
Original file line number Diff line number Diff line change
Expand Up @@ -65,8 +65,6 @@ export CalciumField
export CalciumFieldElem
export CalciumQQBar
export canonical_unit
export cansolve
export cansolve_with_nullspace
export cdiv
export cdivpow2
export cdivrem
Expand Down Expand Up @@ -587,13 +585,6 @@ export SkewDiagram
export snf
export snf_diagonal
export solve
export solve_cholesky_precomp
export solve_cholesky_precomp!
export solve_dixon
export solve_lu_precomp
export solve_lu_precomp!
export solve_rational
export solve!
export sort_terms!
export SparsePolynomialRing
export sqr_classical
Expand Down
20 changes: 10 additions & 10 deletions src/HeckeMiscMatrix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -181,24 +181,24 @@ end
#
################################################################################

function right_kernel(x::fpMatrix)
function _right_kernel(x::fpMatrix)
z = zero_matrix(base_ring(x), ncols(x), max(nrows(x), ncols(x)))
n = ccall((:nmod_mat_nullspace, libflint), Int, (Ref{fpMatrix}, Ref{fpMatrix}), z, x)
return n, z
end

function left_kernel(x::fpMatrix)
n, M = right_kernel(transpose(x))
function _left_kernel(x::fpMatrix)
n, M = _right_kernel(transpose(x))
return n, transpose(M)
end

@doc raw"""
left_kernel(a::ZZMatrix) -> Int, ZZMatrix
_left_kernel(a::ZZMatrix) -> Int, ZZMatrix
It returns a tuple $(n, M)$ where $M$ is a matrix whose rows generate
the kernel of $a$ and $n$ is the rank of the kernel.
"""
function left_kernel(x::ZZMatrix)
function _left_kernel(x::ZZMatrix)
if nrows(x) == 0
return 0, zero(x, 0, 0)
end
Expand All @@ -217,12 +217,12 @@ function left_kernel(x::ZZMatrix)
end
end

function right_kernel(x::ZZMatrix)
n, M = left_kernel(transpose(x))
function _right_kernel(x::ZZMatrix)
n, M = _left_kernel(transpose(x))
return n, transpose(M)
end

function right_kernel(M::zzModMatrix)
function _right_kernel(M::zzModMatrix)
R = base_ring(M)
if is_prime(modulus(R))
k = zero_matrix(R, ncols(M), ncols(M))
Expand Down Expand Up @@ -250,7 +250,7 @@ function right_kernel(M::zzModMatrix)
return 0, zero_matrix(R, nrows(M), 0)
end

function right_kernel(M::ZZModMatrix)
function _right_kernel(M::ZZModMatrix)
R = base_ring(M)
N = hcat(transpose(M), identity_matrix(R, ncols(M)))
if nrows(N) < ncols(N)
Expand Down Expand Up @@ -784,7 +784,7 @@ function eigenspace(M::MatElem{T}, lambda::T; side::Symbol = :left) where T <: F
for i = 1:ncols(N)
N[i, i] -= lambda
end
return kernel(N, side = side)[2]
return kernel(N, side = side)
end

@doc raw"""
Expand Down
14 changes: 7 additions & 7 deletions src/arb/ComplexMat.jl
Original file line number Diff line number Diff line change
Expand Up @@ -566,34 +566,34 @@ function lu!(P::Generic.Perm, x::ComplexMat)
return min(nrows(x), ncols(x))
end

function solve!(z::ComplexMat, x::ComplexMat, y::ComplexMat)
function _solve!(z::ComplexMat, x::ComplexMat, y::ComplexMat)
r = ccall((:acb_mat_solve, libarb), Cint,
(Ref{ComplexMat}, Ref{ComplexMat}, Ref{ComplexMat}, Int),
z, x, y, precision(Balls))
r == 0 && error("Matrix cannot be inverted numerically")
nothing
end

function solve(x::ComplexMat, y::ComplexMat)
function _solve(x::ComplexMat, y::ComplexMat)
ncols(x) != nrows(x) && error("First argument must be square")
ncols(x) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve!(z, x, y)
_solve!(z, x, y)
return z
end

function solve_lu_precomp!(z::ComplexMat, P::Generic.Perm, LU::ComplexMat, y::ComplexMat)
function _solve_lu_precomp!(z::ComplexMat, P::Generic.Perm, LU::ComplexMat, y::ComplexMat)
Q = inv(P)
ccall((:acb_mat_solve_lu_precomp, libarb), Nothing,
(Ref{ComplexMat}, Ptr{Int}, Ref{ComplexMat}, Ref{ComplexMat}, Int),
z, Q.d .- 1, LU, y, precision(Balls))
nothing
end

function solve_lu_precomp(P::Generic.Perm, LU::ComplexMat, y::ComplexMat)
function _solve_lu_precomp(P::Generic.Perm, LU::ComplexMat, y::ComplexMat)
ncols(LU) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve_lu_precomp!(z, P, LU, y)
_solve_lu_precomp!(z, P, LU, y)
return z
end

Expand Down Expand Up @@ -622,7 +622,7 @@ end
#
################################################################################

function AbstractAlgebra.Solve.solve_init(A::ComplexMat)
function solve_init(A::ComplexMat)
return AbstractAlgebra.Solve.SolveCtx{ComplexFieldElem, ComplexMat, ComplexMat}(A)
end

Expand Down
14 changes: 7 additions & 7 deletions src/arb/RealMat.jl
Original file line number Diff line number Diff line change
Expand Up @@ -508,34 +508,34 @@ function lu!(P::Generic.Perm, x::RealMat)
return min(nrows(x), ncols(x))
end

function solve!(z::RealMat, x::RealMat, y::RealMat)
function _solve!(z::RealMat, x::RealMat, y::RealMat)
r = ccall((:arb_mat_solve, libarb), Cint,
(Ref{RealMat}, Ref{RealMat}, Ref{RealMat}, Int),
z, x, y, precision(Balls))
r == 0 && error("Matrix cannot be inverted numerically")
nothing
end

function solve(x::RealMat, y::RealMat)
function _solve(x::RealMat, y::RealMat)
ncols(x) != nrows(x) && error("First argument must be square")
ncols(x) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve!(z, x, y)
_solve!(z, x, y)
return z
end

function solve_lu_precomp!(z::RealMat, P::Generic.Perm, LU::RealMat, y::RealMat)
function _solve_lu_precomp!(z::RealMat, P::Generic.Perm, LU::RealMat, y::RealMat)
Q = inv(P)
ccall((:arb_mat_solve_lu_precomp, libarb), Nothing,
(Ref{RealMat}, Ptr{Int}, Ref{RealMat}, Ref{RealMat}, Int),
z, Q.d .- 1, LU, y, precision(Balls))
nothing
end

function solve_lu_precomp(P::Generic.Perm, LU::RealMat, y::RealMat)
function _solve_lu_precomp(P::Generic.Perm, LU::RealMat, y::RealMat)
ncols(LU) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve_lu_precomp!(z, P, LU, y)
_solve_lu_precomp!(z, P, LU, y)
return z
end

Expand Down Expand Up @@ -564,7 +564,7 @@ end
#
################################################################################

function AbstractAlgebra.Solve.solve_init(A::RealMat)
function solve_init(A::RealMat)
return AbstractAlgebra.Solve.SolveCtx{RealFieldElem, RealMat, RealMat}(A)
end

Expand Down
14 changes: 7 additions & 7 deletions src/arb/acb_mat.jl
Original file line number Diff line number Diff line change
Expand Up @@ -569,34 +569,34 @@ function lu!(P::Generic.Perm, x::AcbMatrix)
return nrows(x)
end

function solve!(z::AcbMatrix, x::AcbMatrix, y::AcbMatrix)
function _solve!(z::AcbMatrix, x::AcbMatrix, y::AcbMatrix)
r = ccall((:acb_mat_solve, libarb), Cint,
(Ref{AcbMatrix}, Ref{AcbMatrix}, Ref{AcbMatrix}, Int),
z, x, y, precision(base_ring(x)))
r == 0 && error("Matrix cannot be inverted numerically")
nothing
end

function solve(x::AcbMatrix, y::AcbMatrix)
function _solve(x::AcbMatrix, y::AcbMatrix)
ncols(x) != nrows(x) && error("First argument must be square")
ncols(x) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve!(z, x, y)
_solve!(z, x, y)
return z
end

function solve_lu_precomp!(z::AcbMatrix, P::Generic.Perm, LU::AcbMatrix, y::AcbMatrix)
function _solve_lu_precomp!(z::AcbMatrix, P::Generic.Perm, LU::AcbMatrix, y::AcbMatrix)
Q = inv(P)
ccall((:acb_mat_solve_lu_precomp, libarb), Nothing,
(Ref{AcbMatrix}, Ptr{Int}, Ref{AcbMatrix}, Ref{AcbMatrix}, Int),
z, Q.d .- 1, LU, y, precision(base_ring(LU)))
nothing
end

function solve_lu_precomp(P::Generic.Perm, LU::AcbMatrix, y::AcbMatrix)
function _solve_lu_precomp(P::Generic.Perm, LU::AcbMatrix, y::AcbMatrix)
ncols(LU) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve_lu_precomp!(z, P, LU, y)
_solve_lu_precomp!(z, P, LU, y)
return z
end

Expand Down Expand Up @@ -625,7 +625,7 @@ end
#
################################################################################

function AbstractAlgebra.Solve.solve_init(A::AcbMatrix)
function solve_init(A::AcbMatrix)
return AbstractAlgebra.Solve.SolveCtx{AcbFieldElem, AcbMatrix, AcbMatrix}(A)
end

Expand Down
20 changes: 10 additions & 10 deletions src/arb/arb_mat.jl
Original file line number Diff line number Diff line change
Expand Up @@ -521,48 +521,48 @@ function lu!(P::Generic.Perm, x::ArbMatrix)
return nrows(x)
end

function solve!(z::ArbMatrix, x::ArbMatrix, y::ArbMatrix)
function _solve!(z::ArbMatrix, x::ArbMatrix, y::ArbMatrix)
r = ccall((:arb_mat_solve, libarb), Cint,
(Ref{ArbMatrix}, Ref{ArbMatrix}, Ref{ArbMatrix}, Int),
z, x, y, precision(base_ring(x)))
r == 0 && error("Matrix cannot be inverted numerically")
nothing
end

function solve(x::ArbMatrix, y::ArbMatrix)
function _solve(x::ArbMatrix, y::ArbMatrix)
ncols(x) != nrows(x) && error("First argument must be square")
ncols(x) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve!(z, x, y)
_solve!(z, x, y)
return z
end

function solve_lu_precomp!(z::ArbMatrix, P::Generic.Perm, LU::ArbMatrix, y::ArbMatrix)
function _solve_lu_precomp!(z::ArbMatrix, P::Generic.Perm, LU::ArbMatrix, y::ArbMatrix)
Q = inv(P)
ccall((:arb_mat_solve_lu_precomp, libarb), Nothing,
(Ref{ArbMatrix}, Ptr{Int}, Ref{ArbMatrix}, Ref{ArbMatrix}, Int),
z, Q.d .- 1, LU, y, precision(base_ring(LU)))
nothing
end

function solve_lu_precomp(P::Generic.Perm, LU::ArbMatrix, y::ArbMatrix)
function _solve_lu_precomp(P::Generic.Perm, LU::ArbMatrix, y::ArbMatrix)
ncols(LU) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve_lu_precomp!(z, P, LU, y)
_solve_lu_precomp!(z, P, LU, y)
return z
end

function solve_cholesky_precomp!(z::ArbMatrix, cho::ArbMatrix, y::ArbMatrix)
function _solve_cholesky_precomp!(z::ArbMatrix, cho::ArbMatrix, y::ArbMatrix)
ccall((:arb_mat_solve_cho_precomp, libarb), Nothing,
(Ref{ArbMatrix}, Ref{ArbMatrix}, Ref{ArbMatrix}, Int),
z, cho, y, precision(base_ring(cho)))
nothing
end

function solve_cholesky_precomp(cho::ArbMatrix, y::ArbMatrix)
function _solve_cholesky_precomp(cho::ArbMatrix, y::ArbMatrix)
ncols(cho) != nrows(y) && error("Matrix dimensions are wrong")
z = similar(y)
solve_cholesky_precomp!(z, cho, y)
_solve_cholesky_precomp!(z, cho, y)
return z
end

Expand Down Expand Up @@ -591,7 +591,7 @@ end
#
################################################################################

function AbstractAlgebra.Solve.solve_init(A::ArbMatrix)
function solve_init(A::ArbMatrix)
return AbstractAlgebra.Solve.SolveCtx{ArbFieldElem, ArbMatrix, ArbMatrix}(A)
end

Expand Down
Loading

0 comments on commit 2098cf8

Please sign in to comment.