Hardcode bands in __bibimul! instead of calling _diag

Author: jishnub

src/bidiag.jl

@@ -698,9 +698,9 @@ function __bibimul!(C, A, B::Bidiagonal, _add)
     Al = _diag(A, -1)
     Ad = _diag(A, 0)
     Au = _diag(A, 1)
-    Bd = _diag(B, 0)
+    Bd = B.dv
     if B.uplo == 'U'
-        Bu = _diag(B, 1)
+        Bu = B.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋2j₋1 = Au[j-2]
@@ -719,7 +719,7 @@ function __bibimul!(C, A, B::Bidiagonal, _add)
             end
         end
     else # B.uplo == 'L'
-        Bl = _diag(B, -1)
+        Bl = B.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋1j   = Au[j-1]
@@ -745,9 +745,9 @@ function __bibimul!(C, A::Bidiagonal, B, _add)
     Bl = _diag(B, -1)
     Bd = _diag(B, 0)
     Bu = _diag(B, 1)
-    Ad = _diag(A, 0)
+    Ad = A.dv
     if A.uplo == 'U'
-        Au = _diag(A, 1)
+        Au = A.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋2j₋1 = Au[j-2]
@@ -767,7 +767,7 @@ function __bibimul!(C, A::Bidiagonal, B, _add)
             end
         end
     else # A.uplo == 'L'
-        Al = _diag(A, -1)
+        Al = A.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋1j₋1 = Ad[j-1]
@@ -791,11 +791,11 @@ function __bibimul!(C, A::Bidiagonal, B, _add)
 end
 function __bibimul!(C, A::Bidiagonal, B::Bidiagonal, _add)
     n = size(A,1)
-    Ad = _diag(A, 0)
-    Bd = _diag(B, 0)
+    Ad = A.dv
+    Bd = B.dv
     if A.uplo == 'U' && B.uplo == 'U'
-        Au = _diag(A, 1)
-        Bu = _diag(B, 1)
+        Au = A.ev
+        Bu = B.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋2j₋1 = Au[j-2]
@@ -811,8 +811,8 @@ function __bibimul!(C, A::Bidiagonal, B::Bidiagonal, _add)
             end
         end
     elseif A.uplo == 'U' && B.uplo == 'L'
-        Au = _diag(A, 1)
-        Bl = _diag(B, -1)
+        Au = A.ev
+        Bl = B.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋1j   = Au[j-1]
@@ -828,8 +828,8 @@ function __bibimul!(C, A::Bidiagonal, B::Bidiagonal, _add)
             end
         end
     elseif A.uplo == 'L' && B.uplo == 'U'
-        Al = _diag(A, -1)
-        Bu = _diag(B, 1)
+        Al = A.ev
+        Bu = B.ev
         @inbounds begin
             for j in 3:n-2
                 Aj₋1j₋1 = Ad[j-1]
@@ -845,8 +845,8 @@ function __bibimul!(C, A::Bidiagonal, B::Bidiagonal, _add)
             end
         end
     else # A.uplo == 'L' && B.uplo == 'L'
-        Al = _diag(A, -1)
-        Bl = _diag(B, -1)
+        Al = A.ev
+        Bl = B.ev
         @inbounds begin
             for j in 3:n-2
                 Ajj     = Ad[j]