mat: Add DiagView methods to types
* mat: Add DiagView methods to types
Fixes
MODIFY
mat/band.go
mat/band.go
@@ -188,6 +188,18 @@ func (b *BandDense) RawBand() blas64.Band {
return b.mat
}
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (b *BandDense) DiagView() Diagonal {
+ n := min(b.mat.Rows, b.mat.Cols)
+ return &DiagDense{
+ mat: blas64.Vector{
+ Inc: b.mat.Stride,
+ Data: b.mat.Data[b.mat.KL : (n-1)*b.mat.Stride+b.mat.KL+1],
+ },
+ n: n,
+ }
+}
+
// DoNonZero calls the function fn for each of the non-zero elements of b. The function fn
// takes a row/column index and the element value of b at (i, j).
func (b *BandDense) DoNonZero(fn func(i, j int, v float64)) {

MODIFY
mat/band_test.go
mat/band_test.go
@@ -248,6 +248,30 @@ func TestNewDiagonalRect(t *testing.T) {
}
}
+func TestBandDiagView(t *testing.T) {
+ for cas, test := range []*BandDense{
+ NewBandDense(1, 1, 0, 0, []float64{1}),
+ NewBandDense(6, 6, 1, 2, []float64{
+ -1, 2, 3, 4,
+ 5, 6, 7, 8,
+ 9, 10, 11, 12,
+ 13, 14, 15, 16,
+ 17, 18, 19, -1,
+ 21, 22, -1, -1,
+ }),
+ NewBandDense(6, 6, 2, 1, []float64{
+ -1, -1, 1, 2,
+ -1, 3, 4, 5,
+ 6, 7, 8, 9,
+ 10, 11, 12, 13,
+ 14, 15, 16, 17,
+ 18, 19, 20, -1,
+ }),
+ } {
+ testDiagView(t, cas, test)
+ }
+}
+
func TestBandAtSet(t *testing.T) {
// 2 3 4 0 0 0
// 5 6 7 8 0 0

MODIFY
mat/dense.go
mat/dense.go
@@ -278,6 +278,18 @@ func (m *Dense) rawRowView(i int) []float64 {
return m.mat.Data[i*m.mat.Stride : i*m.mat.Stride+m.mat.Cols]
}
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (m *Dense) DiagView() Diagonal {
+ n := min(m.mat.Rows, m.mat.Cols)
+ return &DiagDense{
+ mat: blas64.Vector{
+ Inc: m.mat.Stride + 1,
+ Data: m.mat.Data[:(n-1)*m.mat.Stride+n],
+ },
+ n: n,
+ }
+}
+
// Slice returns a new Matrix that shares backing data with the receiver.
// The returned matrix starts at {i,j} of the receiver and extends k-i rows
// and l-j columns. The final row in the resulting matrix is k-1 and the

MODIFY
mat/dense_test.go
mat/dense_test.go
@@ -360,6 +360,26 @@ func TestRowColView(t *testing.T) {
}
}
+func TestDiagView(t *testing.T) {
+ for cas, test := range []*Dense{
+ NewDense(1, 1, []float64{1}),
+ NewDense(2, 2, []float64{1, 2, 3, 4}),
+ NewDense(3, 4, []float64{
+ 1, 2, 3, 4,
+ 5, 6, 7, 8,
+ 9, 10, 11, 12,
+ }),
+ NewDense(4, 3, []float64{
+ 1, 2, 3,
+ 4, 5, 6,
+ 7, 8, 9,
+ 10, 11, 12,
+ }),
+ } {
+ testDiagView(t, cas, test)
+ }
+}
+
func TestGrow(t *testing.T) {
m := &Dense{}
m = m.Grow(10, 10).(*Dense)

MODIFY
mat/diagonal.go
mat/diagonal.go
@@ -170,8 +170,13 @@ func (d *DiagDense) Reset() {
d.mat.Data = d.mat.Data[:0]
}
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (d *DiagDense) DiagView() Diagonal {
+ return d
+}
+
// DiagFrom copies the diagonal of m into the receiver. The receiver must
-// be min(r, c) long or zero. Otherwise DiagOf will panic.
+// be min(r, c) long or zero. Otherwise DiagFrom will panic.
func (d *DiagDense) DiagFrom(m Matrix) {
n := min(m.Dims())
d.reuseAs(n)
@@ -204,7 +209,16 @@ func (d *DiagDense) DiagFrom(m Matrix) {
Inc: mat.Stride + 1,
Data: mat.Data[:(n-1)*mat.Stride+n],
}
- // TODO(kortschak): Add banded triangular handling when the type exists.
+ case RawTriBander:
+ mat := r.RawTriBand()
+ data := mat.Data
+ if mat.Uplo == blas.Lower {
+ data = data[mat.K:]
+ }
+ vec = blas64.Vector{
+ Inc: mat.Stride,
+ Data: data[:(n-1)*mat.Stride+1],
+ }
case RawTriangular:
mat := r.RawTriangular()
if mat.Diag == blas.Unit {

MODIFY
mat/diagonal_test.go
mat/diagonal_test.go
@@ -244,6 +244,42 @@ func TestDiagFrom(t *testing.T) {
0, 0, 0, 0, 0, 6,
}),
},
+ {
+ mat: NewTriBandDense(6, 2, Upper, []float64{
+ 1, math.NaN(), math.NaN(),
+ 2, math.NaN(), math.NaN(),
+ 3, math.NaN(), math.NaN(),
+ 4, math.NaN(), math.NaN(),
+ 5, math.NaN(), math.NaN(),
+ 6, math.NaN(), math.NaN(),
+ }),
+ want: NewDense(6, 6, []float64{
+ 1, 0, 0, 0, 0, 0,
+ 0, 2, 0, 0, 0, 0,
+ 0, 0, 3, 0, 0, 0,
+ 0, 0, 0, 4, 0, 0,
+ 0, 0, 0, 0, 5, 0,
+ 0, 0, 0, 0, 0, 6,
+ }),
+ },
+ {
+ mat: NewTriBandDense(6, 2, Lower, []float64{
+ math.NaN(), math.NaN(), 1,
+ math.NaN(), math.NaN(), 2,
+ math.NaN(), math.NaN(), 3,
+ math.NaN(), math.NaN(), 4,
+ math.NaN(), math.NaN(), 5,
+ math.NaN(), math.NaN(), 6,
+ }),
+ want: NewDense(6, 6, []float64{
+ 1, 0, 0, 0, 0, 0,
+ 0, 2, 0, 0, 0, 0,
+ 0, 0, 3, 0, 0, 0,
+ 0, 0, 0, 4, 0, 0,
+ 0, 0, 0, 0, 5, 0,
+ 0, 0, 0, 0, 0, 6,
+ }),
+ },
{
mat: NewTriDense(6, Upper, []float64{
1, math.NaN(), math.NaN(), math.NaN(), math.NaN(), math.NaN(),
@@ -322,6 +358,42 @@ func TestDiagFrom(t *testing.T) {
}
}
+// diagDenseViewer takes the view of the Diagonal with the underlying Diagonal
+// as the DiagDense type.
+type diagDenseViewer interface {
+ Matrix
+ DiagView() Diagonal
+}
+
+func testDiagView(t *testing.T, cas int, test diagDenseViewer) {
+ // Check the DiagView matches the Diagonal.
+ r, c := test.Dims()
+ diagView := test.DiagView()
+ for i := 0; i < min(r, c); i++ {
+ if diagView.At(i, i) != test.At(i, i) {
+ t.Errorf("Diag mismatch case %d, element %d", cas, i)
+ }
+ }
+
+ // Check that changes to the diagonal are reflected.
+ offset := 10.0
+ diag := diagView.(*DiagDense)
+ for i := 0; i < min(r, c); i++ {
+ v := test.At(i, i)
+ diag.SetDiag(i, v+offset)
+ if test.At(i, i) != v+offset {
+ t.Errorf("Diag set mismatch case %d, element %d", cas, i)
+ }
+ }
+
+ // Check that DiagView and DiagFrom match.
+ var diag2 DiagDense
+ diag2.DiagFrom(test)
+ if !Equal(diag, &diag2) {
+ t.Errorf("Cas %d: DiagView and DiagFrom mismatch", cas)
+ }
+}
+
func TestDiagonalAtSet(t *testing.T) {
for _, n := range []int{1, 3, 8} {
for _, nilstart := range []bool{true, false} {

MODIFY
mat/symband.go
mat/symband.go
@@ -145,6 +145,18 @@ func (s *SymBandDense) RawSymBand() blas64.SymmetricBand {
return s.mat
}
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (s *SymBandDense) DiagView() Diagonal {
+ n := s.mat.N
+ return &DiagDense{
+ mat: blas64.Vector{
+ Inc: s.mat.Stride,
+ Data: s.mat.Data[:(n-1)*s.mat.Stride+1],
+ },
+ n: n,
+ }
+}
+
// DoNonZero calls the function fn for each of the non-zero elements of s. The function fn
// takes a row/column index and the element value of s at (i, j).
func (s *SymBandDense) DoNonZero(fn func(i, j int, v float64)) {

MODIFY
mat/symband_test.go
mat/symband_test.go
@@ -179,3 +179,19 @@ func TestSymBandAtSet(t *testing.T) {
}
}
}
+
+func TestSymBandDiagView(t *testing.T) {
+ for cas, test := range []*SymBandDense{
+ NewSymBandDense(1, 0, []float64{1}),
+ NewSymBandDense(6, 2, []float64{
+ 1, 2, 3,
+ 4, 5, 6,
+ 7, 8, 9,
+ 10, 11, 12,
+ 13, 14, -1,
+ 16, -1, -1,
+ }),
+ } {
+ testDiagView(t, cas, test)
+ }
+}

MODIFY
mat/symmetric.go
mat/symmetric.go
@@ -178,6 +178,18 @@ func (s *SymDense) isolatedWorkspace(a Symmetric) (w *SymDense, restore func())
}
}
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (s *SymDense) DiagView() Diagonal {
+ n := s.mat.N
+ return &DiagDense{
+ mat: blas64.Vector{
+ Inc: s.mat.Stride + 1,
+ Data: s.mat.Data[:(n-1)*s.mat.Stride+n],
+ },
+ n: n,
+ }
+}
+
func (s *SymDense) AddSym(a, b Symmetric) {
n := a.Symmetric()
if n != b.Symmetric() {

MODIFY
mat/symmetric_test.go
mat/symmetric_test.go
@@ -129,6 +129,16 @@ func TestSymAtSet(t *testing.T) {
}
}
+func TestSymDiagView(t *testing.T) {
+ for cas, test := range []*SymDense{
+ NewSymDense(1, []float64{1}),
+ NewSymDense(2, []float64{1, 2, 2, 3}),
+ NewSymDense(3, []float64{1, 2, 3, 2, 4, 5, 3, 5, 6}),
+ } {
+ testDiagView(t, cas, test)
+ }
+}
+
func TestSymAdd(t *testing.T) {
for _, test := range []struct {
n int

MODIFY
mat/triangular.go
mat/triangular.go
@@ -296,6 +296,21 @@ func (t *TriDense) isolatedWorkspace(a Triangular) (w *TriDense, restore func())
}
}
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (t *TriDense) DiagView() Diagonal {
+ if t.mat.Diag == blas.Unit {
+ panic("mat: cannot take view of Unit diagonal")
+ }
+ n := t.mat.N
+ return &DiagDense{
+ mat: blas64.Vector{
+ Inc: t.mat.Stride + 1,
+ Data: t.mat.Data[:(n-1)*t.mat.Stride+n],
+ },
+ n: n,
+ }
+}
+
// Copy makes a copy of elements of a into the receiver. It is similar to the
// built-in copy; it copies as much as the overlap between the two matrices and
// returns the number of rows and columns it copied. Only elements within the

MODIFY
mat/triangular_test.go
mat/triangular_test.go
@@ -139,6 +139,19 @@ func TestTriAtSet(t *testing.T) {
}
}
+func TestTriDiagView(t *testing.T) {
+ for cas, test := range []*TriDense{
+ NewTriDense(1, Upper, []float64{1}),
+ NewTriDense(2, Upper, []float64{1, 2, 0, 3}),
+ NewTriDense(3, Upper, []float64{1, 2, 3, 0, 4, 5, 0, 0, 6}),
+ NewTriDense(1, Lower, []float64{1}),
+ NewTriDense(2, Lower, []float64{1, 2, 2, 3}),
+ NewTriDense(3, Lower, []float64{1, 0, 0, 2, 3, 0, 4, 5, 6}),
+ } {
+ testDiagView(t, cas, test)
+ }
+}
+
func TestTriDenseCopy(t *testing.T) {
for i := 0; i < 100; i++ {
size := rand.Intn(100)

MODIFY
mat/triband.go
mat/triband.go
@@ -305,3 +305,22 @@ func (t *TriBandDense) TTriBand() TriBanded {
func (t *TriBandDense) RawTriBand() blas64.TriangularBand {
return t.mat
}
+
+// DiagView returns the diagonal as a matrix backed by the original data.
+func (t *TriBandDense) DiagView() Diagonal {
+ if t.mat.Diag == blas.Unit {
+ panic("mat: cannot take view of Unit diagonal")
+ }
+ n := t.mat.N
+ data := t.mat.Data
+ if !t.isUpper() {
+ data = data[t.mat.K:]
+ }
+ return &DiagDense{
+ mat: blas64.Vector{
+ Inc: t.mat.Stride,
+ Data: data[:(n-1)*t.mat.Stride+1],
+ },
+ n: n,
+ }
+}

MODIFY
mat/triband_test.go
mat/triband_test.go
@@ -332,3 +332,30 @@ func TestTriBandAtSetUpper(t *testing.T) {
}
}
}
+
+func TestTriBandDiagView(t *testing.T) {
+ for cas, test := range []*TriBandDense{
+ NewTriBandDense(1, 0, Upper, []float64{1}),
+ NewTriBandDense(4, 0, Upper, []float64{1, 2, 3, 4}),
+ NewTriBandDense(6, 2, Upper, []float64{
+ 1, 2, 3,
+ 4, 5, 6,
+ 7, 8, 9,
+ 10, 11, 12,
+ 13, 14, -1,
+ 15, -1, -1,
+ }),
+ NewTriBandDense(1, 0, Lower, []float64{1}),
+ NewTriBandDense(4, 0, Lower, []float64{1, 2, 3, 4}),
+ NewTriBandDense(6, 2, Lower, []float64{
+ -1, -1, 1,
+ -1, 2, 3,
+ 4, 5, 6,
+ 7, 8, 9,
+ 10, 11, 12,
+ 13, 14, 15,
+ }),
+ } {
+ testDiagView(t, cas, test)
+ }
+}

mat: Add DiagView methods to types (#711)

* mat: Add DiagView methods to types

Fixes #213