diff --git a/src/arb/ComplexMat.jl b/src/arb/ComplexMat.jl index 4522d529a..f0b9864ef 100644 --- a/src/arb/ComplexMat.jl +++ b/src/arb/ComplexMat.jl @@ -568,30 +568,7 @@ function lu!(P::Generic.Perm, x::ComplexMat) r == 0 && error("Could not find $(nrows(x)) invertible pivot elements") P.d .+= 1 inv!(P) - return nrows(x) -end - -function lu(P::Generic.Perm, x::ComplexMat) - ncols(x) != nrows(x) && error("Matrix must be square") - parent(P).n != nrows(x) && error("Permutation does not match matrix") - R = base_ring(x) - L = similar(x) - U = deepcopy(x) - n = ncols(x) - lu!(P, U) - for i = 1:n - for j = 1:n - if i > j - L[i, j] = U[i, j] - U[i, j] = R() - elseif i == j - L[i, j] = one(R) - else - L[i, j] = R() - end - end - end - return L, U + return min(nrows(x), ncols(x)) end function solve!(z::ComplexMat, x::ComplexMat, y::ComplexMat) diff --git a/src/arb/RealMat.jl b/src/arb/RealMat.jl index ea08f7a5b..0be07cddf 100644 --- a/src/arb/RealMat.jl +++ b/src/arb/RealMat.jl @@ -502,7 +502,6 @@ end ############################################################################### function lu!(P::Generic.Perm, x::RealMat) - ncols(x) != nrows(x) && error("Matrix must be square") parent(P).n != nrows(x) && error("Permutation does not match matrix") P.d .-= 1 r = ccall((:arb_mat_lu, libarb), Cint, @@ -511,29 +510,7 @@ function lu!(P::Generic.Perm, x::RealMat) r == 0 && error("Could not find $(nrows(x)) invertible pivot elements") P.d .+= 1 inv!(P) - return nrows(x) -end - -function lu(x::RealMat, P = SymmetricGroup(nrows(x))) - p = one(P) - R = base_ring(x) - L = similar(x) - U = deepcopy(x) - n = ncols(x) - r = lu!(p, U) - for i = 1:n - for j = 1:n - if i > j - L[i, j] = U[i, j] - U[i, j] = R() - elseif i == j - L[i, j] = one(R) - else - L[i, j] = R() - end - end - end - return r, p, L, U + return min(nrows(x), ncols(x)) end function solve!(z::RealMat, x::RealMat, y::RealMat) diff --git a/src/arb/acb_mat.jl b/src/arb/acb_mat.jl index d0557baff..0e65402bf 100644 --- a/src/arb/acb_mat.jl +++ b/src/arb/acb_mat.jl @@ -574,29 +574,6 @@ function lu!(P::Generic.Perm, x::acb_mat) return nrows(x) end -function lu(P::Generic.Perm, x::acb_mat) - ncols(x) != nrows(x) && error("Matrix must be square") - parent(P).n != nrows(x) && error("Permutation does not match matrix") - R = base_ring(x) - L = similar(x) - U = deepcopy(x) - n = ncols(x) - lu!(P, U) - for i = 1:n - for j = 1:n - if i > j - L[i, j] = U[i, j] - U[i, j] = R() - elseif i == j - L[i, j] = one(R) - else - L[i, j] = R() - end - end - end - return L, U -end - function solve!(z::acb_mat, x::acb_mat, y::acb_mat) r = ccall((:acb_mat_solve, libarb), Cint, (Ref{acb_mat}, Ref{acb_mat}, Ref{acb_mat}, Int), diff --git a/src/arb/arb_mat.jl b/src/arb/arb_mat.jl index eab79a27b..f9efbd083 100644 --- a/src/arb/arb_mat.jl +++ b/src/arb/arb_mat.jl @@ -506,7 +506,6 @@ end ############################################################################### function lu!(P::Generic.Perm, x::arb_mat) - ncols(x) != nrows(x) && error("Matrix must be square") parent(P).n != nrows(x) && error("Permutation does not match matrix") P.d .-= 1 r = ccall((:arb_mat_lu, libarb), Cint, @@ -518,28 +517,6 @@ function lu!(P::Generic.Perm, x::arb_mat) return nrows(x) end -function lu(x::arb_mat, P = SymmetricGroup(nrows(x))) - p = one(P) - R = base_ring(x) - L = similar(x) - U = deepcopy(x) - n = ncols(x) - r = lu!(p, U) - for i = 1:n - for j = 1:n - if i > j - L[i, j] = U[i, j] - U[i, j] = R() - elseif i == j - L[i, j] = one(R) - else - L[i, j] = R() - end - end - end - return r, p, L, U -end - function solve!(z::arb_mat, x::arb_mat, y::arb_mat) r = ccall((:arb_mat_solve, libarb), Cint, (Ref{arb_mat}, Ref{arb_mat}, Ref{arb_mat}, Int), diff --git a/test/arb/ComplexMat-test.jl b/test/arb/ComplexMat-test.jl index 1a854cfe8..72f18219a 100644 --- a/test/arb/ComplexMat-test.jl +++ b/test/arb/ComplexMat-test.jl @@ -403,6 +403,18 @@ end @test overlaps(B, C) end +@testset "ComplexMat.lu_nonsquare" begin + S = matrix_space(CC, 2, 3) + + A = S(["1.0 +/- 0.01" "1.0 +/- 0.01" "1.0 +/- 0.01"; + "1.0 +/- 0.01" "0.0 +/- 0.01" "-1.0 +/- 0.01"]) + + r, p, L, U = lu(A) + + @test overlaps(L*U, p*A) + @test r == 2 +end + @testset "ComplexMat.linear_solving" begin S = matrix_space(CC, 3, 3) T = matrix_space(ZZ, 3, 3) diff --git a/test/arb/RealMat-test.jl b/test/arb/RealMat-test.jl index b07eb779b..c5598ecff 100644 --- a/test/arb/RealMat-test.jl +++ b/test/arb/RealMat-test.jl @@ -369,6 +369,18 @@ end @test overlaps(B, C) end +@testset "RealMat.lu_nonsquare" begin + S = matrix_space(RR, 2, 3) + + A = S(["1.0 +/- 0.01" "1.0 +/- 0.01" "1.0 +/- 0.01"; + "1.0 +/- 0.01" "0.0 +/- 0.01" "-1.0 +/- 0.01"]) + + r, p, L, U = lu(A) + + @test overlaps(L*U, p*A) + @test r == 2 +end + @testset "RealMat.linear_solving" begin S = matrix_space(RR, 3, 3) T = matrix_space(ZZ, 3, 3) diff --git a/test/arb/acb_mat-test.jl b/test/arb/acb_mat-test.jl index e7fcad142..150187558 100644 --- a/test/arb/acb_mat-test.jl +++ b/test/arb/acb_mat-test.jl @@ -403,6 +403,18 @@ end @test overlaps(B, C) end +@testset "acb_mat.lu_nonsquare" begin + S = matrix_space(CC, 2, 3) + + A = S(["1.0 +/- 0.01" "1.0 +/- 0.01" "1.0 +/- 0.01"; + "1.0 +/- 0.01" "0.0 +/- 0.01" "-1.0 +/- 0.01"]) + + r, p, L, U = lu(A) + + @test overlaps(L*U, p*A) + @test r == 2 +end + @testset "acb_mat.linear_solving" begin S = matrix_space(CC, 3, 3) T = matrix_space(ZZ, 3, 3) diff --git a/test/arb/arb_mat-test.jl b/test/arb/arb_mat-test.jl index 60d882f79..8048c08ca 100644 --- a/test/arb/arb_mat-test.jl +++ b/test/arb/arb_mat-test.jl @@ -369,6 +369,18 @@ end @test overlaps(B, C) end +@testset "arb_mat.lu_nonsquare" begin + S = matrix_space(RR, 2, 3) + + A = S(["1.0 +/- 0.01" "1.0 +/- 0.01" "1.0 +/- 0.01"; + "1.0 +/- 0.01" "0.0 +/- 0.01" "-1.0 +/- 0.01"]) + + r, p, L, U = lu(A) + + @test overlaps(L*U, p*A) + @test r == 2 +end + @testset "arb_mat.linear_solving" begin S = matrix_space(RR, 3, 3) T = matrix_space(ZZ, 3, 3)