lapack: rename DecompUpdate and add GenOrtho types and consts

This commit is contained in:
Vladimir Chalupecky
2018-08-23 11:16:14 +02:00
committed by Vladimír Chalupecký
parent 72366fbe54
commit 6d5ac7aa26
8 changed files with 125 additions and 110 deletions

View File

@@ -112,9 +112,9 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
lwork_dorgqr_m := int(work[0]) lwork_dorgqr_m := int(work[0])
impl.Dgebrd(n, n, a, lda, s, nil, nil, nil, work, -1) impl.Dgebrd(n, n, a, lda, s, nil, nil, nil, work, -1)
lwork_dgebrd := int(work[0]) lwork_dgebrd := int(work[0])
impl.Dorgbr(lapack.ApplyP, n, n, n, a, lda, nil, work, -1) impl.Dorgbr(lapack.GeneratePT, n, n, n, a, lda, nil, work, -1)
lwork_dorgbr_p := int(work[0]) lwork_dorgbr_p := int(work[0])
impl.Dorgbr(lapack.ApplyQ, n, n, n, a, lda, nil, work, -1) impl.Dorgbr(lapack.GenerateQ, n, n, n, a, lda, nil, work, -1)
lwork_dorgbr_q := int(work[0]) lwork_dorgbr_q := int(work[0])
if m >= mnthr { if m >= mnthr {
@@ -203,12 +203,12 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
lwork_dgebrd := int(work[0]) lwork_dgebrd := int(work[0])
maxwrk = 3*n + lwork_dgebrd maxwrk = 3*n + lwork_dgebrd
if wantus || wantuo { if wantus || wantuo {
impl.Dorgbr(lapack.ApplyQ, m, n, n, a, lda, nil, work, -1) impl.Dorgbr(lapack.GenerateQ, m, n, n, a, lda, nil, work, -1)
lwork_dorgbr_q = int(work[0]) lwork_dorgbr_q = int(work[0])
maxwrk = max(maxwrk, 3*n+lwork_dorgbr_q) maxwrk = max(maxwrk, 3*n+lwork_dorgbr_q)
} }
if wantua { if wantua {
impl.Dorgbr(lapack.ApplyQ, m, m, n, a, lda, nil, work, -1) impl.Dorgbr(lapack.GenerateQ, m, m, n, a, lda, nil, work, -1)
lwork_dorgbr_q := int(work[0]) lwork_dorgbr_q := int(work[0])
maxwrk = max(maxwrk, 3*n+lwork_dorgbr_q) maxwrk = max(maxwrk, 3*n+lwork_dorgbr_q)
} }
@@ -229,9 +229,9 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
lwork_dorglq_m := int(work[0]) lwork_dorglq_m := int(work[0])
impl.Dgebrd(m, m, a, lda, s, nil, nil, nil, work, -1) impl.Dgebrd(m, m, a, lda, s, nil, nil, nil, work, -1)
lwork_dgebrd := int(work[0]) lwork_dgebrd := int(work[0])
impl.Dorgbr(lapack.ApplyP, m, m, m, a, n, nil, work, -1) impl.Dorgbr(lapack.GeneratePT, m, m, m, a, n, nil, work, -1)
lwork_dorgbr_p := int(work[0]) lwork_dorgbr_p := int(work[0])
impl.Dorgbr(lapack.ApplyQ, m, m, m, a, n, nil, work, -1) impl.Dorgbr(lapack.GenerateQ, m, m, m, a, n, nil, work, -1)
lwork_dorgbr_q := int(work[0]) lwork_dorgbr_q := int(work[0])
if n >= mnthr { if n >= mnthr {
// n >> m // n >> m
@@ -319,12 +319,12 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
lwork_dgebrd = int(work[0]) lwork_dgebrd = int(work[0])
maxwrk = 3*m + lwork_dgebrd maxwrk = 3*m + lwork_dgebrd
if wantvs || wantvo { if wantvs || wantvo {
impl.Dorgbr(lapack.ApplyP, m, n, m, a, n, nil, work, -1) impl.Dorgbr(lapack.GeneratePT, m, n, m, a, n, nil, work, -1)
lwork_dorgbr_p = int(work[0]) lwork_dorgbr_p = int(work[0])
maxwrk = max(maxwrk, 3*m+lwork_dorgbr_p) maxwrk = max(maxwrk, 3*m+lwork_dorgbr_p)
} }
if wantva { if wantva {
impl.Dorgbr(lapack.ApplyP, n, n, m, a, n, nil, work, -1) impl.Dorgbr(lapack.GeneratePT, n, n, m, a, n, nil, work, -1)
lwork_dorgbr_p = int(work[0]) lwork_dorgbr_p = int(work[0])
maxwrk = max(maxwrk, 3*m+lwork_dorgbr_p) maxwrk = max(maxwrk, 3*m+lwork_dorgbr_p)
} }
@@ -398,8 +398,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
ncvt := 0 ncvt := 0
if wantvo || wantvas { if wantvo || wantvas {
// Generate P^T. impl.Dorgbr(lapack.GeneratePT, n, n, n, a, lda, work[itaup:],
impl.Dorgbr(lapack.ApplyP, n, n, n, a, lda, work[itaup:],
work[iwork:], lwork-iwork) work[iwork:], lwork-iwork)
ncvt = n ncvt = n
} }
@@ -453,7 +452,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
work[itauq:], work[itaup:], work[iwork:], lwork-iwork) work[itauq:], work[itaup:], work[iwork:], lwork-iwork)
// Generate left vectors bidiagonalizing R in work[ir:]. // Generate left vectors bidiagonalizing R in work[ir:].
impl.Dorgbr(lapack.ApplyQ, n, n, n, work[ir:], ldworkr, impl.Dorgbr(lapack.GenerateQ, n, n, n, work[ir:], ldworkr,
work[itauq:], work[iwork:], lwork-iwork) work[itauq:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -536,11 +535,11 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
impl.Dlacpy(blas.Upper, n, n, work[iu:], ldworku, vt, ldvt) impl.Dlacpy(blas.Upper, n, n, work[iu:], ldworku, vt, ldvt)
// Generate left bidiagonalizing vectors in work[iu:]. // Generate left bidiagonalizing vectors in work[iu:].
impl.Dorgbr(lapack.ApplyQ, n, n, n, work[iu:], ldworku, impl.Dorgbr(lapack.GenerateQ, n, n, n, work[iu:], ldworku,
work[itauq:], work[iwork:], lwork-iwork) work[itauq:], work[iwork:], lwork-iwork)
// Generate right bidiagonalizing vectors in VT. // Generate right bidiagonalizing vectors in VT.
impl.Dorgbr(lapack.ApplyP, n, n, n, vt, ldvt, impl.Dorgbr(lapack.GeneratePT, n, n, n, vt, ldvt,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -584,7 +583,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
vt, ldvt, work[itauq:], u, ldu, work[iwork:], lwork-iwork) vt, ldvt, work[itauq:], u, ldu, work[iwork:], lwork-iwork)
// Generate right bidiagonalizing vectors in VT. // Generate right bidiagonalizing vectors in VT.
impl.Dorgbr(lapack.ApplyP, n, n, n, vt, ldvt, impl.Dorgbr(lapack.GeneratePT, n, n, n, vt, ldvt,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -630,7 +629,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
work[itauq:], work[itaup:], work[iwork:], lwork-iwork) work[itauq:], work[itaup:], work[iwork:], lwork-iwork)
// Generate left bidiagonalizing vectors in work[ir:]. // Generate left bidiagonalizing vectors in work[ir:].
impl.Dorgbr(lapack.ApplyQ, n, n, n, work[ir:], ldworkr, impl.Dorgbr(lapack.GenerateQ, n, n, n, work[ir:], ldworkr,
work[itauq:], work[iwork:], lwork-iwork) work[itauq:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -718,11 +717,11 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
impl.Dlacpy(blas.Upper, n, n, work[iu:], ldworku, vt, ldvt) impl.Dlacpy(blas.Upper, n, n, work[iu:], ldworku, vt, ldvt)
// Generate left bidiagonalizing vectors in work[iu:]. // Generate left bidiagonalizing vectors in work[iu:].
impl.Dorgbr(lapack.ApplyQ, n, n, n, work[iu:], ldworku, impl.Dorgbr(lapack.GenerateQ, n, n, n, work[iu:], ldworku,
work[itauq:], work[iwork:], lwork-iwork) work[itauq:], work[iwork:], lwork-iwork)
// Generate right bidiagonalizing vectors in VT. // Generate right bidiagonalizing vectors in VT.
impl.Dorgbr(lapack.ApplyP, n, n, n, vt, ldvt, impl.Dorgbr(lapack.GeneratePT, n, n, n, vt, ldvt,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -750,7 +749,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
m, n, n, vt, ldvt, work[itauq:], u, ldu, work[iwork:], lwork-iwork) m, n, n, vt, ldvt, work[itauq:], u, ldu, work[iwork:], lwork-iwork)
// Generate right bidiagonalizing vectors in VT. // Generate right bidiagonalizing vectors in VT.
impl.Dorgbr(lapack.ApplyP, n, n, n, vt, ldvt, impl.Dorgbr(lapack.GeneratePT, n, n, n, vt, ldvt,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -790,7 +789,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
m, n, n, vt, ldvt, work[itauq:], u, ldu, work[iwork:], lwork-iwork) m, n, n, vt, ldvt, work[itauq:], u, ldu, work[iwork:], lwork-iwork)
// Generate right bidiagonizing vectors in VT. // Generate right bidiagonizing vectors in VT.
impl.Dorgbr(lapack.ApplyP, n, n, n, vt, ldvt, impl.Dorgbr(lapack.GeneratePT, n, n, n, vt, ldvt,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + n iwork = ie + n
@@ -824,13 +823,13 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
if wantua { if wantua {
ncu = m ncu = m
} }
impl.Dorgbr(lapack.ApplyQ, m, ncu, n, u, ldu, work[itauq:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GenerateQ, m, ncu, n, u, ldu, work[itauq:], work[iwork:], lwork-iwork)
} }
if wantvas { if wantvas {
// Right singular vectors are desired in VT. Copy result to VT and // Right singular vectors are desired in VT. Copy result to VT and
// generate left biadiagonalizing vectors in VT. // generate left biadiagonalizing vectors in VT.
impl.Dlacpy(blas.Upper, n, n, a, lda, vt, ldvt) impl.Dlacpy(blas.Upper, n, n, a, lda, vt, ldvt)
impl.Dorgbr(lapack.ApplyP, n, n, n, vt, ldvt, work[itaup:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GeneratePT, n, n, n, vt, ldvt, work[itaup:], work[iwork:], lwork-iwork)
} }
if wantuo { if wantuo {
panic(noSVDO) panic(noSVDO)
@@ -886,7 +885,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
impl.Dgebrd(m, m, a, lda, s, work[ie:itauq], impl.Dgebrd(m, m, a, lda, s, work[ie:itauq],
work[itauq:itaup], work[itaup:iwork], work[iwork:], lwork-iwork) work[itauq:itaup], work[itaup:iwork], work[iwork:], lwork-iwork)
if wantuo || wantuas { if wantuo || wantuas {
impl.Dorgbr(lapack.ApplyQ, m, m, m, a, lda, impl.Dorgbr(lapack.GenerateQ, m, m, m, a, lda,
work[itauq:], work[iwork:], lwork-iwork) work[itauq:], work[iwork:], lwork-iwork)
} }
iwork = ie + m iwork = ie + m
@@ -944,7 +943,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
work[itauq:], work[itaup:], work[iwork:], lwork-iwork) work[itauq:], work[itaup:], work[iwork:], lwork-iwork)
// Generate right vectors bidiagonalizing L in work[ir:]. // Generate right vectors bidiagonalizing L in work[ir:].
impl.Dorgbr(lapack.ApplyP, m, m, m, work[ir:], ldworkr, impl.Dorgbr(lapack.GeneratePT, m, m, m, work[ir:], ldworkr,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + m iwork = ie + m
@@ -1029,11 +1028,11 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
impl.Dlacpy(blas.Lower, m, m, work[iu:], ldworku, u, ldu) impl.Dlacpy(blas.Lower, m, m, work[iu:], ldworku, u, ldu)
// Generate right bidiagionalizing vectors in work[iu:]. // Generate right bidiagionalizing vectors in work[iu:].
impl.Dorgbr(lapack.ApplyP, m, m, m, work[iu:], ldworku, impl.Dorgbr(lapack.GeneratePT, m, m, m, work[iu:], ldworku,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
// Generate left bidiagonalizing vectors in U. // Generate left bidiagonalizing vectors in U.
impl.Dorgbr(lapack.ApplyQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GenerateQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork)
iwork = ie + m iwork = ie + m
// Perform bidiagonal QR iteration, computing left singular // Perform bidiagonal QR iteration, computing left singular
@@ -1076,7 +1075,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
u, ldu, work[itaup:], vt, ldvt, work[iwork:], lwork-iwork) u, ldu, work[itaup:], vt, ldvt, work[iwork:], lwork-iwork)
// Generate left bidiagonalizing vectors in U. // Generate left bidiagonalizing vectors in U.
impl.Dorgbr(lapack.ApplyQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GenerateQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork)
iwork = ie + m iwork = ie + m
// Perform bidiagonal QR iteration, computing left singular // Perform bidiagonal QR iteration, computing left singular
@@ -1122,7 +1121,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
work[itauq:], work[itaup:], work[iwork:], lwork-iwork) work[itauq:], work[itaup:], work[iwork:], lwork-iwork)
// Generate right bidiagonalizing vectors in work[ir:]. // Generate right bidiagonalizing vectors in work[ir:].
impl.Dorgbr(lapack.ApplyP, m, m, m, work[ir:], ldworkr, impl.Dorgbr(lapack.GeneratePT, m, m, m, work[ir:], ldworkr,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
iwork = ie + m iwork = ie + m
@@ -1209,11 +1208,11 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
impl.Dlacpy(blas.Lower, m, m, work[iu:], ldworku, u, ldu) impl.Dlacpy(blas.Lower, m, m, work[iu:], ldworku, u, ldu)
// Generate right bidiagonalizing vectors in work[iu:]. // Generate right bidiagonalizing vectors in work[iu:].
impl.Dorgbr(lapack.ApplyP, m, m, m, work[iu:], ldworku, impl.Dorgbr(lapack.GeneratePT, m, m, m, work[iu:], ldworku,
work[itaup:], work[iwork:], lwork-iwork) work[itaup:], work[iwork:], lwork-iwork)
// Generate left bidiagonalizing vectors in U. // Generate left bidiagonalizing vectors in U.
impl.Dorgbr(lapack.ApplyQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GenerateQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork)
iwork = ie + m iwork = ie + m
// Perform bidiagonal QR iteration, computing left singular // Perform bidiagonal QR iteration, computing left singular
@@ -1259,7 +1258,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
u, ldu, work[itaup:], vt, ldvt, work[iwork:], lwork-iwork) u, ldu, work[itaup:], vt, ldvt, work[iwork:], lwork-iwork)
// Generate left bidiagonalizing vectors in U. // Generate left bidiagonalizing vectors in U.
impl.Dorgbr(lapack.ApplyQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GenerateQ, m, m, m, u, ldu, work[itauq:], work[iwork:], lwork-iwork)
iwork = ie + m iwork = ie + m
// Perform bidiagonal QR iteration, computing left singular // Perform bidiagonal QR iteration, computing left singular
@@ -1284,7 +1283,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
// If left singular vectors desired in U, copy result to U and // If left singular vectors desired in U, copy result to U and
// generate left bidiagonalizing vectors in U. // generate left bidiagonalizing vectors in U.
impl.Dlacpy(blas.Lower, m, m, a, lda, u, ldu) impl.Dlacpy(blas.Lower, m, m, a, lda, u, ldu)
impl.Dorgbr(lapack.ApplyQ, m, m, n, u, ldu, work[itauq:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GenerateQ, m, m, n, u, ldu, work[itauq:], work[iwork:], lwork-iwork)
} }
if wantvas { if wantvas {
// If right singular vectors desired in VT, copy result to VT // If right singular vectors desired in VT, copy result to VT
@@ -1296,7 +1295,7 @@ func (impl Implementation) Dgesvd(jobU, jobVT lapack.SVDJob, m, n int, a []float
} else { } else {
nrvt = m nrvt = m
} }
impl.Dorgbr(lapack.ApplyP, nrvt, n, m, vt, ldvt, work[itaup:], work[iwork:], lwork-iwork) impl.Dorgbr(lapack.GeneratePT, nrvt, n, m, vt, ldvt, work[itaup:], work[iwork:], lwork-iwork)
} }
if wantuo { if wantuo {
panic(noSVDO) panic(noSVDO)

View File

@@ -10,18 +10,25 @@ import "gonum.org/v1/gonum/lapack"
// computed from the decomposition Dgebrd. See Dgebd2 for the description of // computed from the decomposition Dgebrd. See Dgebd2 for the description of
// Q and P^T. // Q and P^T.
// //
// If vect == lapack.ApplyQ, then a is assumed to have been an m×k matrix and // If vect == lapack.GenerateQ, then a is assumed to have been an m×k matrix and
// Q is of order m. If m >= k, then Dorgbr returns the first n columns of Q // Q is of order m. If m >= k, then Dorgbr returns the first n columns of Q
// where m >= n >= k. If m < k, then Dorgbr returns Q as an m×m matrix. // where m >= n >= k. If m < k, then Dorgbr returns Q as an m×m matrix.
// //
// If vect == lapack.ApplyP, then A is assumed to have been a k×n matrix, and // If vect == lapack.GeneratePT, then A is assumed to have been a k×n matrix, and
// P^T is of order n. If k < n, then Dorgbr returns the first m rows of P^T, // P^T is of order n. If k < n, then Dorgbr returns the first m rows of P^T,
// where n >= m >= k. If k >= n, then Dorgbr returns P^T as an n×n matrix. // where n >= m >= k. If k >= n, then Dorgbr returns P^T as an n×n matrix.
// //
// Dorgbr is an internal routine. It is exported for testing purposes. // Dorgbr is an internal routine. It is exported for testing purposes.
func (impl Implementation) Dorgbr(vect lapack.DecompUpdate, m, n, k int, a []float64, lda int, tau, work []float64, lwork int) { func (impl Implementation) Dorgbr(vect lapack.GenOrtho, m, n, k int, a []float64, lda int, tau, work []float64, lwork int) {
mn := min(m, n) mn := min(m, n)
wantq := vect == lapack.ApplyQ var wantq bool
switch vect {
case lapack.GenerateQ:
wantq = true
case lapack.GeneratePT:
default:
panic(badGenOrtho)
}
if wantq { if wantq {
if m < n || n < min(m, k) || m < min(m, k) { if m < n || n < min(m, k) || m < min(m, k) {
panic(badDims) panic(badDims)

View File

@@ -44,7 +44,7 @@ import (
// returns it in work[0]. // returns it in work[0].
// //
// Dormbr is an internal routine. It is exported for testing purposes. // Dormbr is an internal routine. It is exported for testing purposes.
func (impl Implementation) Dormbr(vect lapack.DecompUpdate, side blas.Side, trans blas.Transpose, m, n, k int, a []float64, lda int, tau, c []float64, ldc int, work []float64, lwork int) { func (impl Implementation) Dormbr(vect lapack.ApplyOrtho, side blas.Side, trans blas.Transpose, m, n, k int, a []float64, lda int, tau, c []float64, ldc int, work []float64, lwork int) {
if side != blas.Left && side != blas.Right { if side != blas.Left && side != blas.Right {
panic(badSide) panic(badSide)
} }
@@ -52,7 +52,7 @@ func (impl Implementation) Dormbr(vect lapack.DecompUpdate, side blas.Side, tran
panic(badTrans) panic(badTrans)
} }
if vect != lapack.ApplyP && vect != lapack.ApplyQ { if vect != lapack.ApplyP && vect != lapack.ApplyQ {
panic(badDecompUpdate) panic(badApplyOrtho)
} }
nq := n nq := n
nw := m nw := m

View File

@@ -17,61 +17,62 @@ var _ lapack.Float64 = Implementation{}
// This list is duplicated in netlib/lapack/netlib. Keep in sync. // This list is duplicated in netlib/lapack/netlib. Keep in sync.
const ( const (
absIncNotOne = "lapack: increment not one or negative one" absIncNotOne = "lapack: increment not one or negative one"
badAlpha = "lapack: bad alpha length" badAlpha = "lapack: bad alpha length"
badAuxv = "lapack: auxv has insufficient length" badApplyOrtho = "lapack: bad ApplyOrtho"
badBeta = "lapack: bad beta length" badAuxv = "lapack: auxv has insufficient length"
badD = "lapack: d has insufficient length" badBeta = "lapack: bad beta length"
badDecompUpdate = "lapack: bad decomp update" badD = "lapack: d has insufficient length"
badDiag = "lapack: bad diag" badDiag = "lapack: bad diag"
badDims = "lapack: bad input dimensions" badDims = "lapack: bad input dimensions"
badDirect = "lapack: bad direct" badDirect = "lapack: bad direct"
badE = "lapack: e has insufficient length" badE = "lapack: e has insufficient length"
badEVComp = "lapack: bad EVComp" badEVComp = "lapack: bad EVComp"
badEVHowMany = "lapack: bad EVHowMany" badEVHowMany = "lapack: bad EVHowMany"
badEVJob = "lapack: bad EVJob" badEVJob = "lapack: bad EVJob"
badEVSide = "lapack: bad EVSide" badEVSide = "lapack: bad EVSide"
badGSVDJob = "lapack: bad GSVDJob" badGenOrtho = "lapack: bad GenOrtho"
badIlo = "lapack: ilo out of range" badGSVDJob = "lapack: bad GSVDJob"
badIhi = "lapack: ihi out of range" badIlo = "lapack: ilo out of range"
badIpiv = "lapack: bad permutation length" badIhi = "lapack: ihi out of range"
badBalanceJob = "lapack: bad BalanceJob" badIpiv = "lapack: bad permutation length"
badK1 = "lapack: k1 out of range" badBalanceJob = "lapack: bad BalanceJob"
badK2 = "lapack: k2 out of range" badK1 = "lapack: k1 out of range"
badKperm = "lapack: incorrect permutation length" badK2 = "lapack: k2 out of range"
badLdA = "lapack: index of a out of range" badKperm = "lapack: incorrect permutation length"
badNb = "lapack: nb out of range" badLdA = "lapack: index of a out of range"
badNorm = "lapack: bad norm" badNb = "lapack: nb out of range"
badPivot = "lapack: bad pivot" badNorm = "lapack: bad norm"
badS = "lapack: s has insufficient length" badPivot = "lapack: bad pivot"
badSchurJob = "lapack: bad SchurJob" badS = "lapack: s has insufficient length"
badShifts = "lapack: bad shifts" badSchurJob = "lapack: bad SchurJob"
badSide = "lapack: bad side" badShifts = "lapack: bad shifts"
badSlice = "lapack: bad input slice length" badSide = "lapack: bad side"
badSort = "lapack: bad Sort" badSlice = "lapack: bad input slice length"
badStore = "lapack: bad store" badSort = "lapack: bad Sort"
badTau = "lapack: tau has insufficient length" badStore = "lapack: bad store"
badTauQ = "lapack: tauQ has insufficient length" badTau = "lapack: tau has insufficient length"
badTauP = "lapack: tauP has insufficient length" badTauQ = "lapack: tauQ has insufficient length"
badTrans = "lapack: bad trans" badTauP = "lapack: tauP has insufficient length"
badVn1 = "lapack: vn1 has insufficient length" badTrans = "lapack: bad trans"
badVn2 = "lapack: vn2 has insufficient length" badVn1 = "lapack: vn1 has insufficient length"
badUplo = "lapack: illegal triangle" badVn2 = "lapack: vn2 has insufficient length"
badWork = "lapack: insufficient working memory" badUplo = "lapack: illegal triangle"
badZ = "lapack: insufficient z length" badWork = "lapack: insufficient working memory"
kGTM = "lapack: k > m" badZ = "lapack: insufficient z length"
kGTN = "lapack: k > n" kGTM = "lapack: k > m"
kLT0 = "lapack: k < 0" kGTN = "lapack: k > n"
mLT0 = "lapack: m < 0" kLT0 = "lapack: k < 0"
mLTN = "lapack: m < n" mLT0 = "lapack: m < 0"
nanScale = "lapack: NaN scale factor" mLTN = "lapack: m < n"
negDimension = "lapack: negative matrix dimension" nanScale = "lapack: NaN scale factor"
negZ = "lapack: negative z value" negDimension = "lapack: negative matrix dimension"
nLT0 = "lapack: n < 0" negZ = "lapack: negative z value"
nLTM = "lapack: n < m" nLT0 = "lapack: n < 0"
offsetGTM = "lapack: offset > m" nLTM = "lapack: n < m"
shortWork = "lapack: working array shorter than declared" offsetGTM = "lapack: offset > m"
zeroDiv = "lapack: zero divisor" shortWork = "lapack: working array shorter than declared"
zeroDiv = "lapack: zero divisor"
) )
// checkMatrix verifies the parameters of a matrix input. // checkMatrix verifies the parameters of a matrix input.

View File

@@ -90,11 +90,20 @@ const (
Bottom Pivot = 'B' Bottom Pivot = 'B'
) )
type DecompUpdate byte // ApplyOrtho specifies which orthogonal matrix is applied in Dormbr.
type ApplyOrtho byte
const ( const (
ApplyP DecompUpdate = 'P' ApplyP ApplyOrtho = 'P' // Apply P or P^T.
ApplyQ DecompUpdate = 'Q' ApplyQ ApplyOrtho = 'Q' // Apply Q or Q^T.
)
// GenOrtho specifies which orthogonal matrix is generated in Dorgbr.
type GenOrtho byte
const (
GeneratePT GenOrtho = 'P' // Generate P^T.
GenerateQ GenOrtho = 'Q' // Generate Q.
) )
// SVDJob specifies the singular vector computation type for SVD. // SVDJob specifies the singular vector computation type for SVD.

View File

@@ -15,13 +15,13 @@ import (
) )
type Dorgbrer interface { type Dorgbrer interface {
Dorgbr(vect lapack.DecompUpdate, m, n, k int, a []float64, lda int, tau, work []float64, lwork int) Dorgbr(vect lapack.GenOrtho, m, n, k int, a []float64, lda int, tau, work []float64, lwork int)
Dgebrder Dgebrder
} }
func DorgbrTest(t *testing.T, impl Dorgbrer) { func DorgbrTest(t *testing.T, impl Dorgbrer) {
rnd := rand.New(rand.NewSource(1)) rnd := rand.New(rand.NewSource(1))
for _, vect := range []lapack.DecompUpdate{lapack.ApplyQ, lapack.ApplyP} { for _, vect := range []lapack.GenOrtho{lapack.GenerateQ, lapack.GeneratePT} {
for _, test := range []struct { for _, test := range []struct {
m, n, k, lda int m, n, k, lda int
}{ }{
@@ -52,7 +52,7 @@ func DorgbrTest(t *testing.T, impl Dorgbrer) {
k := test.k k := test.k
lda := test.lda lda := test.lda
// Filter out bad tests // Filter out bad tests
if vect == lapack.ApplyQ { if vect == lapack.GenerateQ {
if m < n || n < min(m, k) || m < min(m, k) { if m < n || n < min(m, k) || m < min(m, k) {
continue continue
} }
@@ -63,7 +63,7 @@ func DorgbrTest(t *testing.T, impl Dorgbrer) {
} }
// Sizes for Dorgbr. // Sizes for Dorgbr.
var ma, na int var ma, na int
if vect == lapack.ApplyQ { if vect == lapack.GenerateQ {
if m >= k { if m >= k {
ma = m ma = m
na = k na = k
@@ -83,7 +83,7 @@ func DorgbrTest(t *testing.T, impl Dorgbrer) {
// a eventually needs to store either P or Q, so it must be // a eventually needs to store either P or Q, so it must be
// sufficiently big. // sufficiently big.
var a []float64 var a []float64
if vect == lapack.ApplyQ { if vect == lapack.GenerateQ {
lda = max(m, lda) lda = max(m, lda)
a = make([]float64, m*lda) a = make([]float64, m*lda)
} else { } else {
@@ -110,7 +110,7 @@ func DorgbrTest(t *testing.T, impl Dorgbrer) {
copy(aCopy, a) copy(aCopy, a)
var tau []float64 var tau []float64
if vect == lapack.ApplyQ { if vect == lapack.GenerateQ {
tau = tauQ tau = tauQ
} else { } else {
tau = tauP tau = tauP
@@ -124,20 +124,20 @@ func DorgbrTest(t *testing.T, impl Dorgbrer) {
var ans blas64.General var ans blas64.General
var nRows, nCols int var nRows, nCols int
equal := true equal := true
if vect == lapack.ApplyQ { if vect == lapack.GenerateQ {
nRows = m nRows = m
nCols = m nCols = m
if m >= k { if m >= k {
nCols = n nCols = n
} }
ans = constructQPBidiagonal(vect, ma, na, min(m, k), aCopy, lda, tau) ans = constructQPBidiagonal(lapack.ApplyQ, ma, na, min(m, k), aCopy, lda, tau)
} else { } else {
nRows = n nRows = n
if k < n { if k < n {
nRows = m nRows = m
} }
nCols = n nCols = n
ansTmp := constructQPBidiagonal(vect, ma, na, min(k, n), aCopy, lda, tau) ansTmp := constructQPBidiagonal(lapack.ApplyP, ma, na, min(k, n), aCopy, lda, tau)
// Dorgbr actually computes P^T // Dorgbr actually computes P^T
ans = transposeGeneral(ansTmp) ans = transposeGeneral(ansTmp)
} }
@@ -149,8 +149,7 @@ func DorgbrTest(t *testing.T, impl Dorgbrer) {
} }
} }
if !equal { if !equal {
applyQ := vect == lapack.ApplyQ t.Errorf("Extracted matrix mismatch. gen = %v, m = %v, n = %v, k = %v", string(vect), m, n, k)
t.Errorf("Extracted matrix mismatch. applyQ: %v, m = %v, n = %v, k = %v", applyQ, m, n, k)
} }
} }
} }

View File

@@ -16,14 +16,14 @@ import (
) )
type Dormbrer interface { type Dormbrer interface {
Dormbr(vect lapack.DecompUpdate, side blas.Side, trans blas.Transpose, m, n, k int, a []float64, lda int, tau, c []float64, ldc int, work []float64, lwork int) Dormbr(vect lapack.ApplyOrtho, side blas.Side, trans blas.Transpose, m, n, k int, a []float64, lda int, tau, c []float64, ldc int, work []float64, lwork int)
Dgebrder Dgebrder
} }
func DormbrTest(t *testing.T, impl Dormbrer) { func DormbrTest(t *testing.T, impl Dormbrer) {
rnd := rand.New(rand.NewSource(1)) rnd := rand.New(rand.NewSource(1))
bi := blas64.Implementation() bi := blas64.Implementation()
for _, vect := range []lapack.DecompUpdate{lapack.ApplyQ, lapack.ApplyP} { for _, vect := range []lapack.ApplyOrtho{lapack.ApplyQ, lapack.ApplyP} {
for _, side := range []blas.Side{blas.Left, blas.Right} { for _, side := range []blas.Side{blas.Left, blas.Right} {
for _, trans := range []blas.Transpose{blas.NoTrans, blas.Trans} { for _, trans := range []blas.Transpose{blas.NoTrans, blas.Trans} {
for _, wl := range []worklen{minimumWork, mediumWork, optimumWork} { for _, wl := range []worklen{minimumWork, mediumWork, optimumWork} {

View File

@@ -681,7 +681,7 @@ func checkBidiagonal(t *testing.T, m, n, nb int, a []float64, lda int, d, e, tau
// constructQPBidiagonal constructs Q or P from the Bidiagonal decomposition // constructQPBidiagonal constructs Q or P from the Bidiagonal decomposition
// computed by dlabrd and bgebd2. // computed by dlabrd and bgebd2.
func constructQPBidiagonal(vect lapack.DecompUpdate, m, n, nb int, a []float64, lda int, tau []float64) blas64.General { func constructQPBidiagonal(vect lapack.ApplyOrtho, m, n, nb int, a []float64, lda int, tau []float64) blas64.General {
sz := n sz := n
if vect == lapack.ApplyQ { if vect == lapack.ApplyQ {
sz = m sz = m