remove _left_orth and _right_orth workarounds

Author: lkdvos

src/operators/ortho.jl

@@ -19,7 +19,7 @@ function left_canonicalize!(
     # QR of second column
     if size(W, 1) == 1
         tmp = transpose(C′, ((2, 1), (3,)))
-        Q, R = _left_orth!(tmp; alg, trunc = trscheme)
+        Q, R = left_orth!(tmp; alg, trunc = trscheme)
 
         if dim(R) == 0 # fully truncated
             V = oneunit(S) ⊞ oneunit(S)
@@ -35,7 +35,7 @@ function left_canonicalize!(
         H[i] = JordanMPOTensor(codomain(W) ← physicalspace(W) ⊗ V, Q1, W.B, Q2, W.D)
     else
         tmp = transpose(cat(insertleftunit(C′, 1), W.A; dims = 1), ((3, 1, 2), (4,)))
-        Q, R = _left_orth!(tmp; alg, trunc = trscheme)
+        Q, R = left_orth!(tmp; alg, trunc = trscheme)
 
         if dim(R) == 0 # fully truncated
             V = oneunit(S) ⊞ oneunit(S)
@@ -105,7 +105,7 @@ function right_canonicalize!(
     # LQ of second row
     if size(W, 4) == 1
         tmp = transpose(B′, ((1,), (3, 2)))
-        R, Q = _right_orth!(tmp; alg, trunc = trscheme)
+        R, Q = right_orth!(tmp; alg, trunc = trscheme)
 
         if dim(R) == 0
             V = oneunit(S) ⊞ oneunit(S)
@@ -121,7 +121,7 @@ function right_canonicalize!(
         H[i] = JordanMPOTensor(V ⊗ physicalspace(W) ← domain(W), Q1, Q2, W.C, W.D)
     else
         tmp = transpose(cat(insertleftunit(B′, 4), W.A; dims = 4), ((1,), (3, 4, 2)))
-        R, Q = _right_orth!(tmp; alg, trunc = trscheme)
+        R, Q = right_orth!(tmp; alg, trunc = trscheme)
         if dim(R) == 0
             V = oneunit(S) ⊞ oneunit(S)
             Q1 = typeof(W.A)(undef, SumSpace{S}() ⊗ physicalspace(W) ← physicalspace(W) ⊗ SumSpace{S}())

src/states/abstractmps.jl

@@ -119,11 +119,11 @@ function makefullrank!(A::PeriodicVector{<:GenericMPSTensor}; alg = Defaults.alg
         i = findfirst(!isfullrank, A)
         isnothing(i) && break
         if !isfullrank(A[i]; side = :left)
-            L, Q = _right_orth!(_transpose_tail(A[i]); alg)
+            L, Q = right_orth!(_transpose_tail(A[i]); alg)
             A[i] = _transpose_front(Q)
             A[i - 1] = A[i - 1] * L
         else
-            A[i], R = _left_orth!(A[i]; alg)
+            A[i], R = left_orth!(A[i]; alg)
             A[i + 1] = _transpose_front(R * _transpose_tail(A[i + 1]))
         end
     end

src/states/ortho.jl

@@ -164,8 +164,8 @@ regauge!
 function regauge!(
         AC::GenericMPSTensor, C::MPSBondTensor; alg = Defaults.alg_qr()
     )
-    Q_AC, _ = _left_orth!(AC; alg)
-    Q_C, _ = _left_orth!(C; alg)
+    Q_AC, _ = left_orth!(AC; alg)
+    Q_C, _ = left_orth!(C; alg)
     return mul!(AC, Q_AC, Q_C')
 end
 function regauge!(AC::Vector{<:GenericMPSTensor}, C::Vector{<:MPSBondTensor}; kwargs...)
@@ -178,8 +178,8 @@ function regauge!(
         CL::MPSBondTensor, AC::GenericMPSTensor; alg = Defaults.alg_lq()
     )
     AC_tail = _transpose_tail(AC)
-    _, Q_AC = _right_orth!(AC_tail; alg)
-    _, Q_C = _right_orth!(CL; alg)
+    _, Q_AC = right_orth!(AC_tail; alg)
+    _, Q_C = right_orth!(CL; alg)
     AR_tail = mul!(AC_tail, Q_C', Q_AC)
     return repartition!(AC, AR_tail)
 end
@@ -240,7 +240,7 @@ function gauge_eigsolve_step!(it::IterativeSolver{LeftCanonical}, state)
     if iter ≥ it.eig_miniter
         alg_eigsolve = updatetol(it.alg_eigsolve, 1, ϵ^2)
         _, vec = fixedpoint(flip(TransferMatrix(A, AL)), C[end], :LM, alg_eigsolve)
-        _, C[end] = _left_orth!(vec; alg = it.alg_orth)
+        _, C[end] = left_orth!(vec; alg = it.alg_orth)
     end
     return C[end]
 end
@@ -251,7 +251,7 @@ function gauge_orth_step!(it::IterativeSolver{LeftCanonical}, state)
         # repartition!(A_tail[i], AL[i])
         mul!(CA_tail[i], C[i - 1], A_tail[i])
         repartition!(AL[i], CA_tail[i])
-        AL[i], C[i] = _left_orth!(AL[i]; alg = it.alg_orth)
+        AL[i], C[i] = left_orth!(AL[i]; alg = it.alg_orth)
     end
     normalize!(C[end])
     return C[end]
@@ -298,7 +298,7 @@ function gauge_eigsolve_step!(it::IterativeSolver{RightCanonical}, state)
     if iter ≥ it.eig_miniter
         alg_eigsolve = updatetol(it.alg_eigsolve, 1, ϵ^2)
         _, vec = fixedpoint(TransferMatrix(A, AR), C[end], :LM, alg_eigsolve)
-        C[end], _ = _right_orth!(vec; alg = it.alg_orth)
+        C[end], _ = right_orth!(vec; alg = it.alg_orth)
     end
     return C[end]
 end
@@ -308,7 +308,7 @@ function gauge_orth_step!(it::IterativeSolver{RightCanonical}, state)
     for i in length(AR):-1:1
         AC = mul!(AR[i], A[i], C[i])   # use AR as temporary storage for A * C
         tmp = repartition!(AC_tail[i], AC)
-        C[i - 1], tmp = _right_orth!(tmp; alg = it.alg_orth)
+        C[i - 1], tmp = right_orth!(tmp; alg = it.alg_orth)
         repartition!(AR[i], tmp)       # TODO: avoid doing this every iteration
     end
     normalize!(C[end])

src/utility/utility.jl

@@ -149,31 +149,3 @@ function check_unambiguous_braiding(V::VectorSpace)
     return check_unambiguous_braiding(Bool, V) ||
         throw(ArgumentError("cannot unambiguously braid $V"))
 end
-
-# temporary workaround for the fact that left_orth and right_orth are poorly designed:
-function _left_orth!(t; alg::MatrixAlgebraKit.AbstractAlgorithm, trunc::MatrixAlgebraKit.TruncationStrategy = notrunc())
-    if alg isa LAPACK_HouseholderQR
-        return left_orth!(t; kind = :qr, alg_qr = alg, trunc)
-    elseif alg isa LAPACK_HouseholderLQ
-        return left_orth!(t; kind = :qr, alg_qr = LAPACK_HouseholderQR(; alg.kwargs...), trunc)
-    elseif alg isa PolarViaSVD
-        return left_orth!(t; kind = :polar, alg_polar = alg, trunc)
-    elseif alg isa LAPACK_SVDAlgorithm
-        return left_orth!(t; kind = :svd, alg_svd = alg, trunc)
-    else
-        error(lazy"unkown algorithm $alg")
-    end
-end
-function _right_orth!(t; alg::MatrixAlgebraKit.AbstractAlgorithm, trunc::TruncationStrategy = notrunc())
-    if alg isa LAPACK_HouseholderLQ
-        return right_orth!(t; kind = :lq, alg_lq = alg, trunc)
-    elseif alg isa LAPACK_HouseholderQR
-        return right_orth!(t; kind = :lq, alg_lq = LAPACK_HouseholderLQ(; alg.kwargs...), trunc)
-    elseif alg isa PolarViaSVD
-        return right_orth!(t; kind = :polar, alg_polar = alg, trunc)
-    elseif alg isa LAPACK_SVDAlgorithm
-        return right_orth!(t; kind = :svd, alg_svd = alg, trunc)
-    else
-        error(lazy"unkown algorithm $alg")
-    end
-end