lapack: add SchurComp and SchurUpdateComp types and consts, rename EVComp consts
MODIFY
lapack/gonum/dgeev.go
lapack/gonum/dgeev.go
@@ -116,7 +116,7 @@ func (impl Implementation) Dgeev(jobvl lapack.LeftEVJob, jobvr lapack.RightEVJob
maxwrk := 2*n + n*impl.Ilaenv(1, "DGEHRD", " ", n, 1, n, 0)
if wantvl || wantvr {
maxwrk = max(maxwrk, 2*n+(n-1)*impl.Ilaenv(1, "DORGHR", " ", n, 1, n, -1))
- impl.Dhseqr(lapack.EigenvaluesAndSchur, lapack.OriginalEV, n, 0, n-1,
+ impl.Dhseqr(lapack.EigenvaluesAndSchur, lapack.SchurOrig, n, 0, n-1,
nil, 1, nil, nil, nil, 1, work, -1)
maxwrk = max(maxwrk, max(n+1, n+int(work[0])))
side := lapack.EVLeft
@@ -176,7 +176,7 @@ func (impl Implementation) Dgeev(jobvl lapack.LeftEVJob, jobvr lapack.RightEVJob
impl.Dorghr(n, ilo, ihi, vl, ldvl, tau, work[iwrk:], lwork-iwrk)
// Perform QR iteration, accumulating Schur vectors in VL.
iwrk = n
- first = impl.Dhseqr(lapack.EigenvaluesAndSchur, lapack.OriginalEV, n, ilo, ihi,
+ first = impl.Dhseqr(lapack.EigenvaluesAndSchur, lapack.SchurOrig, n, ilo, ihi,
a, lda, wr, wi, vl, ldvl, work[iwrk:], lwork-iwrk)
if wantvr {
// Want left and right eigenvectors.
@@ -192,7 +192,7 @@ func (impl Implementation) Dgeev(jobvl lapack.LeftEVJob, jobvr lapack.RightEVJob
impl.Dorghr(n, ilo, ihi, vr, ldvr, tau, work[iwrk:], lwork-iwrk)
// Perform QR iteration, accumulating Schur vectors in VR.
iwrk = n
- first = impl.Dhseqr(lapack.EigenvaluesAndSchur, lapack.OriginalEV, n, ilo, ihi,
+ first = impl.Dhseqr(lapack.EigenvaluesAndSchur, lapack.SchurOrig, n, ilo, ihi,
a, lda, wr, wi, vr, ldvr, work[iwrk:], lwork-iwrk)
} else {
// Compute eigenvalues only.

MODIFY
lapack/gonum/dhseqr.go
lapack/gonum/dhseqr.go
@@ -27,11 +27,11 @@ import (
// be computed.
// For other values of job Dhseqr will panic.
//
-// If compz == lapack.None, no Schur vectors will be computed and Z will not be
+// If compz == lapack.SchurNone, no Schur vectors will be computed and Z will not be
// referenced.
-// If compz == lapack.HessEV, on return Z will contain the matrix of Schur
+// If compz == lapack.SchurHess, on return Z will contain the matrix of Schur
// vectors of H.
-// If compz == lapack.OriginalEV, on entry z is assumed to contain the orthogonal
+// If compz == lapack.SchurOrig, on entry z is assumed to contain the orthogonal
// matrix Q that is the identity except for the submatrix
// Q[ilo:ihi+1,ilo:ihi+1]. On return z will be updated to the product Q*Z.
//
@@ -96,11 +96,11 @@ import (
// where U is an orthogonal matrix. The final H is upper Hessenberg and
// H[unconverged:ihi+1,unconverged:ihi+1] is upper quasi-triangular.
//
-// If unconverged > 0 and compz == lapack.OriginalEV, then on return
+// If unconverged > 0 and compz == lapack.SchurOrig, then on return
// (final Z) = (initial Z) U,
// where U is the orthogonal matrix in (*) regardless of the value of job.
//
-// If unconverged > 0 and compz == lapack.HessEV, then on return
+// If unconverged > 0 and compz == lapack.SchurHess, then on return
// (final Z) = U,
// where U is the orthogonal matrix in (*) regardless of the value of job.
//
@@ -118,7 +118,7 @@ import (
// URL: http://dx.doi.org/10.1137/S0895479801384585
//
// Dhseqr is an internal routine. It is exported for testing purposes.
-func (impl Implementation) Dhseqr(job lapack.SchurJob, compz lapack.EVComp, n, ilo, ihi int, h []float64, ldh int, wr, wi []float64, z []float64, ldz int, work []float64, lwork int) (unconverged int) {
+func (impl Implementation) Dhseqr(job lapack.SchurJob, compz lapack.SchurComp, n, ilo, ihi int, h []float64, ldh int, wr, wi []float64, z []float64, ldz int, work []float64, lwork int) (unconverged int) {
var wantt bool
switch job {
default:
@@ -130,9 +130,9 @@ func (impl Implementation) Dhseqr(job lapack.SchurJob, compz lapack.EVComp, n, i
var wantz bool
switch compz {
default:
- panic(badEVComp)
- case lapack.None:
- case lapack.HessEV, lapack.OriginalEV:
+ panic(badSchurComp)
+ case lapack.SchurNone:
+ case lapack.SchurHess, lapack.SchurOrig:
wantz = true
}
switch {
@@ -197,7 +197,7 @@ func (impl Implementation) Dhseqr(job lapack.SchurJob, compz lapack.EVComp, n, i
}
// Initialize Z to identity matrix if requested.
- if compz == lapack.HessEV {
+ if compz == lapack.SchurHess {
impl.Dlaset(blas.All, n, n, 0, 1, z, ldz)
}

MODIFY
lapack/gonum/dtrexc.go
lapack/gonum/dtrexc.go
@@ -18,8 +18,8 @@ import "gonum.org/v1/gonum/lapack"
// as Z^T*T*Z, and will be again in Schur canonical form.
//
// If compq is lapack.UpdateSchur, on return the matrix Q of Schur vectors will be
-// updated by postmultiplying it with Z.
-// If compq is lapack.None, the matrix Q is not referenced and will not be
+// updated by post-multiplying it with Z.
+// If compq is lapack.UpdateSchurNone, the matrix Q is not referenced and will not be
// updated.
// For other values of compq Dtrexc will panic.
//
@@ -45,13 +45,13 @@ import "gonum.org/v1/gonum/lapack"
// work must have length at least n, otherwise Dtrexc will panic.
//
// Dtrexc is an internal routine. It is exported for testing purposes.
-func (impl Implementation) Dtrexc(compq lapack.EVComp, n int, t []float64, ldt int, q []float64, ldq int, ifst, ilst int, work []float64) (ifstOut, ilstOut int, ok bool) {
+func (impl Implementation) Dtrexc(compq lapack.UpdateSchurComp, n int, t []float64, ldt int, q []float64, ldq int, ifst, ilst int, work []float64) (ifstOut, ilstOut int, ok bool) {
checkMatrix(n, n, t, ldt)
var wantq bool
switch compq {
default:
panic("lapack: bad value of compq")
- case lapack.None:
+ case lapack.UpdateSchurNone:
// Nothing to do because wantq is already false.
case lapack.UpdateSchur:
wantq = true

MODIFY
lapack/gonum/general.go
lapack/gonum/general.go
@@ -45,6 +45,7 @@ const (
badNorm = "lapack: bad norm"
badPivot = "lapack: bad pivot"
badS = "lapack: s has insufficient length"
+ badSchurComp = "lapack: bad SchurComp"
badSchurJob = "lapack: bad SchurJob"
badShifts = "lapack: bad shifts"
badSide = "lapack: bad side"

MODIFY
lapack/lapack.go
lapack/lapack.go
@@ -127,23 +127,13 @@ const (
GSVDNone GSVDJob = 'N' // Do not compute orthogonal matrix
)
-// EVComp specifies how eigenvectors are computed.
+// EVComp specifies how eigenvectors are computed in Dsteqr.
type EVComp byte
const (
- // OriginalEV specifies to compute the eigenvectors of the original
- // matrix.
- OriginalEV EVComp = 'V'
- // TridiagEV specifies to compute both the eigenvectors of the input
- // tridiagonal matrix.
- TridiagEV EVComp = 'I'
- // HessEV specifies to compute both the eigenvectors of the input upper
- // Hessenberg matrix.
- HessEV EVComp = 'I'
-
- // UpdateSchur specifies that the matrix of Schur vectors will be
- // updated by Dtrexc.
- UpdateSchur EVComp = 'V'
+ EVOrig EVComp = 'V' // Compute eigenvectors of the original symmetric matrix.
+ EVTridiag EVComp = 'I' // Compute eigenvectors of the tridiagonal matrix.
+ EVCompNone EVComp = 'N' // Do not compute eigenvectors.
)
// EVJob specifies whether eigenvectors are computed in Dsyev.
@@ -188,6 +178,23 @@ const (
EigenvaluesAndSchur SchurJob = 'S'
)
+// SchurComp specifies whether and how the Schur vectors are computed in Dhseqr.
+type SchurComp byte
+
+const (
+ SchurNone SchurComp = 'N' // Schur vectors are not computed.
+ SchurHess SchurComp = 'I' // Schur vectors of the upper Hessenberg marix are computed.
+ SchurOrig SchurComp = 'V' // Schur vectors of the original matrix are computed.
+)
+
+// UpdateSchurComp specifies whether the matrix of Schur vectors is updated in Dtrexc.
+type UpdateSchurComp byte
+
+const (
+ UpdateSchur UpdateSchurComp = 'V' // The matrix of Schur vectors is updated.
+ UpdateSchurNone UpdateSchurComp = 'N' // The matrix of Schur vectors is not updated.
+)
+
// EVSide specifies what eigenvectors are computed in Dtrevc3.
type EVSide byte

MODIFY
lapack/testlapack/dhseqr.go
lapack/testlapack/dhseqr.go
@@ -16,7 +16,7 @@ import (
)
type Dhseqrer interface {
- Dhseqr(job lapack.SchurJob, compz lapack.EVComp, n, ilo, ihi int, h []float64, ldh int, wr, wi []float64,
+ Dhseqr(job lapack.SchurJob, compz lapack.SchurComp, n, ilo, ihi int, h []float64, ldh int, wr, wi []float64,
z []float64, ldz int, work []float64, lwork int) int
}
@@ -57,11 +57,11 @@ func testDhseqr(t *testing.T, impl Dhseqrer, i int, test dhseqrTest, job lapack.
copyGeneral(h, blas64.General{Rows: n, Cols: n, Stride: max(1, n), Data: test.h})
hCopy := cloneGeneral(h)
- var compz lapack.EVComp = lapack.None
+ compz := lapack.SchurNone
z := blas64.General{Stride: max(1, n)}
if wantz {
// First, let Dhseqr initialize Z to the identity matrix.
- compz = lapack.HessEV
+ compz = lapack.SchurHess
z = nanGeneral(n, n, n+extra)
}
@@ -70,7 +70,7 @@ func testDhseqr(t *testing.T, impl Dhseqrer, i int, test dhseqrTest, job lapack.
work := nanSlice(max(1, n))
if optwork {
- impl.Dhseqr(job, lapack.HessEV, n, ilo, ihi, nil, h.Stride, nil, nil, nil, z.Stride, work, -1)
+ impl.Dhseqr(job, lapack.SchurHess, n, ilo, ihi, nil, h.Stride, nil, nil, nil, z.Stride, work, -1)
work = nanSlice(int(work[0]))
}
@@ -196,7 +196,7 @@ func testDhseqr(t *testing.T, impl Dhseqrer, i int, test dhseqrTest, job lapack.
copyGeneral(h, hCopy)
// Call Dhseqr again with the identity matrix given explicitly in Q.
q := eye(n, n+extra)
- impl.Dhseqr(job, lapack.OriginalEV, n, ilo, ihi, h.Data, h.Stride, wr, wi, q.Data, q.Stride, work, len(work))
+ impl.Dhseqr(job, lapack.SchurOrig, n, ilo, ihi, h.Data, h.Stride, wr, wi, q.Data, q.Stride, work, len(work))
if !equalApproxGeneral(z, q, 0) {
t.Errorf("%v: Z and Q are not equal", prefix)
}

MODIFY
lapack/testlapack/dsteqr.go
lapack/testlapack/dsteqr.go
@@ -22,7 +22,7 @@ type Dsteqrer interface {
func DsteqrTest(t *testing.T, impl Dsteqrer) {
rnd := rand.New(rand.NewSource(1))
- for _, compz := range []lapack.EVComp{lapack.OriginalEV, lapack.TridiagEV} {
+ for _, compz := range []lapack.EVComp{lapack.EVOrig, lapack.EVTridiag} {
for _, test := range []struct {
n, lda int
}{
@@ -59,7 +59,7 @@ func DsteqrTest(t *testing.T, impl Dsteqrer) {
copy(eCopy, e)
aCopy := make([]float64, len(a))
copy(aCopy, a)
- if compz == lapack.OriginalEV {
+ if compz == lapack.EVOrig {
uplo := blas.Upper
tau := make([]float64, n)
work := make([]float64, 1)
@@ -92,7 +92,7 @@ func DsteqrTest(t *testing.T, impl Dsteqrer) {
copy(dAns, d)
var truth blas64.General
- if compz == lapack.OriginalEV {
+ if compz == lapack.EVOrig {
truth = blas64.General{
Rows: n,
Cols: n,
@@ -130,7 +130,7 @@ func DsteqrTest(t *testing.T, impl Dsteqrer) {
}
if !eigenDecompCorrect(d, truth, V) {
t.Errorf("Eigen reconstruction mismatch. fromFull = %v, n = %v",
- compz == lapack.OriginalEV, n)
+ compz == lapack.EVOrig, n)
}
// Compare eigenvalues when not computing eigenvectors.

MODIFY
lapack/testlapack/dtrexc.go
lapack/testlapack/dtrexc.go
@@ -17,13 +17,13 @@ import (
)
type Dtrexcer interface {
- Dtrexc(compq lapack.EVComp, n int, t []float64, ldt int, q []float64, ldq int, ifst, ilst int, work []float64) (ifstOut, ilstOut int, ok bool)
+ Dtrexc(compq lapack.UpdateSchurComp, n int, t []float64, ldt int, q []float64, ldq int, ifst, ilst int, work []float64) (ifstOut, ilstOut int, ok bool)
}
func DtrexcTest(t *testing.T, impl Dtrexcer) {
rnd := rand.New(rand.NewSource(1))
- for _, compq := range []lapack.EVComp{lapack.None, lapack.UpdateSchur} {
+ for _, compq := range []lapack.UpdateSchurComp{lapack.UpdateSchurNone, lapack.UpdateSchur} {
for _, n := range []int{1, 2, 3, 4, 5, 6, 10, 18, 31, 53} {
for _, extra := range []int{0, 1, 11} {
for cas := 0; cas < 100; cas++ {
@@ -36,7 +36,7 @@ func DtrexcTest(t *testing.T, impl Dtrexcer) {
}
}
- for _, compq := range []lapack.EVComp{lapack.None, lapack.UpdateSchur} {
+ for _, compq := range []lapack.UpdateSchurComp{lapack.UpdateSchurNone, lapack.UpdateSchur} {
for _, extra := range []int{0, 1, 11} {
tmat := randomSchurCanonical(0, extra, rnd)
testDtrexc(t, impl, compq, tmat, 0, 0, extra, rnd)
@@ -44,7 +44,7 @@ func DtrexcTest(t *testing.T, impl Dtrexcer) {
}
}
-func testDtrexc(t *testing.T, impl Dtrexcer, compq lapack.EVComp, tmat blas64.General, ifst, ilst, extra int, rnd *rand.Rand) {
+func testDtrexc(t *testing.T, impl Dtrexcer, compq lapack.UpdateSchurComp, tmat blas64.General, ifst, ilst, extra int, rnd *rand.Rand) {
const tol = 1e-13
n := tmat.Rows

lapack: add SchurComp and SchurUpdateComp types and consts, rename EVComp consts