vendor/gonum.org/v1/gonum/mat/vector.go

   1 // Copyright ©2013 The Gonum Authors. All rights reserved.
   2 // Use of this source code is governed by a BSD-style
   3 // license that can be found in the LICENSE file.
   4
   5 package mat
   6
   7 import (
   8         "gonum.org/v1/gonum/blas"
   9         "gonum.org/v1/gonum/blas/blas64"
  10         "gonum.org/v1/gonum/internal/asm/f64"
  11 )
  12
  13 var (
  14         vector *VecDense
  15
  16         _ Matrix  = vector
  17         _ Vector  = vector
  18         _ Reseter = vector
  19 )
  20
  21 // Vector is a vector.
  22 type Vector interface {
  23         Matrix
  24         AtVec(int) float64
  25         Len() int
  26 }
  27
  28 // TransposeVec is a type for performing an implicit transpose of a Vector.
  29 // It implements the Vector interface, returning values from the transpose
  30 // of the vector within.
  31 type TransposeVec struct {
  32         Vector Vector
  33 }
  34
  35 // At returns the value of the element at row i and column j of the transposed
  36 // matrix, that is, row j and column i of the Vector field.
  37 func (t TransposeVec) At(i, j int) float64 {
  38         return t.Vector.At(j, i)
  39 }
  40
  41 // AtVec returns the element at position i. It panics if i is out of bounds.
  42 func (t TransposeVec) AtVec(i int) float64 {
  43         return t.Vector.AtVec(i)
  44 }
  45
  46 // Dims returns the dimensions of the transposed vector.
  47 func (t TransposeVec) Dims() (r, c int) {
  48         c, r = t.Vector.Dims()
  49         return r, c
  50 }
  51
  52 // T performs an implicit transpose by returning the Vector field.
  53 func (t TransposeVec) T() Matrix {
  54         return t.Vector
  55 }
  56
  57 // Len returns the number of columns in the vector.
  58 func (t TransposeVec) Len() int {
  59         return t.Vector.Len()
  60 }
  61
  62 // TVec performs an implicit transpose by returning the Vector field.
  63 func (t TransposeVec) TVec() Vector {
  64         return t.Vector
  65 }
  66
  67 // Untranspose returns the Vector field.
  68 func (t TransposeVec) Untranspose() Matrix {
  69         return t.Vector
  70 }
  71
  72 func (t TransposeVec) UntransposeVec() Vector {
  73         return t.Vector
  74 }
  75
  76 // VecDense represents a column vector.
  77 type VecDense struct {
  78         mat blas64.Vector
  79         n   int
  80         // A BLAS vector can have a negative increment, but allowing this
  81         // in the mat type complicates a lot of code, and doesn't gain anything.
  82         // VecDense must have positive increment in this package.
  83 }
  84
  85 // NewVecDense creates a new VecDense of length n. If data == nil,
  86 // a new slice is allocated for the backing slice. If len(data) == n, data is
  87 // used as the backing slice, and changes to the elements of the returned VecDense
  88 // will be reflected in data. If neither of these is true, NewVecDense will panic.
  89 func NewVecDense(n int, data []float64) *VecDense {
  90         if len(data) != n && data != nil {
  91                 panic(ErrShape)
  92         }
  93         if data == nil {
  94                 data = make([]float64, n)
  95         }
  96         return &VecDense{
  97                 mat: blas64.Vector{
  98                         Inc:  1,
  99                         Data: data,
 100                 },
 101                 n: n,
 102         }
 103 }
 104
 105 // SliceVec returns a new Vector that shares backing data with the receiver.
 106 // The returned matrix starts at i of the receiver and extends k-i elements.
 107 // SliceVec panics with ErrIndexOutOfRange if the slice is outside the capacity
 108 // of the receiver.
 109 func (v *VecDense) SliceVec(i, k int) Vector {
 110         if i < 0 || k <= i || v.Cap() < k {
 111                 panic(ErrIndexOutOfRange)
 112         }
 113         return &VecDense{
 114                 n: k - i,
 115                 mat: blas64.Vector{
 116                         Inc:  v.mat.Inc,
 117                         Data: v.mat.Data[i*v.mat.Inc : (k-1)*v.mat.Inc+1],
 118                 },
 119         }
 120 }
 121
 122 // Dims returns the number of rows and columns in the matrix. Columns is always 1
 123 // for a non-Reset vector.
 124 func (v *VecDense) Dims() (r, c int) {
 125         if v.IsZero() {
 126                 return 0, 0
 127         }
 128         return v.n, 1
 129 }
 130
 131 // Caps returns the number of rows and columns in the backing matrix. Columns is always 1
 132 // for a non-Reset vector.
 133 func (v *VecDense) Caps() (r, c int) {
 134         if v.IsZero() {
 135                 return 0, 0
 136         }
 137         return v.Cap(), 1
 138 }
 139
 140 // Len returns the length of the vector.
 141 func (v *VecDense) Len() int {
 142         return v.n
 143 }
 144
 145 // Cap returns the capacity of the vector.
 146 func (v *VecDense) Cap() int {
 147         if v.IsZero() {
 148                 return 0
 149         }
 150         return (cap(v.mat.Data)-1)/v.mat.Inc + 1
 151 }
 152
 153 // T performs an implicit transpose by returning the receiver inside a Transpose.
 154 func (v *VecDense) T() Matrix {
 155         return Transpose{v}
 156 }
 157
 158 // TVec performs an implicit transpose by returning the receiver inside a TransposeVec.
 159 func (v *VecDense) TVec() Vector {
 160         return TransposeVec{v}
 161 }
 162
 163 // Reset zeros the length of the vector so that it can be reused as the
 164 // receiver of a dimensionally restricted operation.
 165 //
 166 // See the Reseter interface for more information.
 167 func (v *VecDense) Reset() {
 168         // No change of Inc or n to 0 may be
 169         // made unless both are set to 0.
 170         v.mat.Inc = 0
 171         v.n = 0
 172         v.mat.Data = v.mat.Data[:0]
 173 }
 174
 175 // CloneVec makes a copy of a into the receiver, overwriting the previous value
 176 // of the receiver.
 177 func (v *VecDense) CloneVec(a Vector) {
 178         if v == a {
 179                 return
 180         }
 181         v.n = a.Len()
 182         v.mat = blas64.Vector{
 183                 Inc:  1,
 184                 Data: use(v.mat.Data, v.n),
 185         }
 186         if r, ok := a.(RawVectorer); ok {
 187                 blas64.Copy(v.n, r.RawVector(), v.mat)
 188                 return
 189         }
 190         for i := 0; i < a.Len(); i++ {
 191                 v.SetVec(i, a.AtVec(i))
 192         }
 193 }
 194
 195 // VecDenseCopyOf returns a newly allocated copy of the elements of a.
 196 func VecDenseCopyOf(a Vector) *VecDense {
 197         v := &VecDense{}
 198         v.CloneVec(a)
 199         return v
 200 }
 201
 202 func (v *VecDense) RawVector() blas64.Vector {
 203         return v.mat
 204 }
 205
 206 // CopyVec makes a copy of elements of a into the receiver. It is similar to the
 207 // built-in copy; it copies as much as the overlap between the two vectors and
 208 // returns the number of elements it copied.
 209 func (v *VecDense) CopyVec(a Vector) int {
 210         n := min(v.Len(), a.Len())
 211         if v == a {
 212                 return n
 213         }
 214         if r, ok := a.(RawVectorer); ok {
 215                 blas64.Copy(n, r.RawVector(), v.mat)
 216                 return n
 217         }
 218         for i := 0; i < n; i++ {
 219                 v.setVec(i, a.AtVec(i))
 220         }
 221         return n
 222 }
 223
 224 // ScaleVec scales the vector a by alpha, placing the result in the receiver.
 225 func (v *VecDense) ScaleVec(alpha float64, a Vector) {
 226         n := a.Len()
 227
 228         if v == a {
 229                 if v.mat.Inc == 1 {
 230                         f64.ScalUnitary(alpha, v.mat.Data)
 231                         return
 232                 }
 233                 f64.ScalInc(alpha, v.mat.Data, uintptr(n), uintptr(v.mat.Inc))
 234                 return
 235         }
 236
 237         v.reuseAs(n)
 238
 239         if rv, ok := a.(RawVectorer); ok {
 240                 mat := rv.RawVector()
 241                 v.checkOverlap(mat)
 242                 if v.mat.Inc == 1 && mat.Inc == 1 {
 243                         f64.ScalUnitaryTo(v.mat.Data, alpha, mat.Data)
 244                         return
 245                 }
 246                 f64.ScalIncTo(v.mat.Data, uintptr(v.mat.Inc),
 247                         alpha, mat.Data, uintptr(n), uintptr(mat.Inc))
 248                 return
 249         }
 250
 251         for i := 0; i < n; i++ {
 252                 v.setVec(i, alpha*a.AtVec(i))
 253         }
 254 }
 255
 256 // AddScaledVec adds the vectors a and alpha*b, placing the result in the receiver.
 257 func (v *VecDense) AddScaledVec(a Vector, alpha float64, b Vector) {
 258         if alpha == 1 {
 259                 v.AddVec(a, b)
 260                 return
 261         }
 262         if alpha == -1 {
 263                 v.SubVec(a, b)
 264                 return
 265         }
 266
 267         ar := a.Len()
 268         br := b.Len()
 269
 270         if ar != br {
 271                 panic(ErrShape)
 272         }
 273
 274         var amat, bmat blas64.Vector
 275         fast := true
 276         aU, _ := untranspose(a)
 277         if rv, ok := aU.(RawVectorer); ok {
 278                 amat = rv.RawVector()
 279                 if v != a {
 280                         v.checkOverlap(amat)
 281                 }
 282         } else {
 283                 fast = false
 284         }
 285         bU, _ := untranspose(b)
 286         if rv, ok := bU.(RawVectorer); ok {
 287                 bmat = rv.RawVector()
 288                 if v != b {
 289                         v.checkOverlap(bmat)
 290                 }
 291         } else {
 292                 fast = false
 293         }
 294
 295         v.reuseAs(ar)
 296
 297         switch {
 298         case alpha == 0: // v <- a
 299                 if v == a {
 300                         return
 301                 }
 302                 v.CopyVec(a)
 303         case v == a && v == b: // v <- v + alpha * v = (alpha + 1) * v
 304                 blas64.Scal(ar, alpha+1, v.mat)
 305         case !fast: // v <- a + alpha * b without blas64 support.
 306                 for i := 0; i < ar; i++ {
 307                         v.setVec(i, a.AtVec(i)+alpha*b.AtVec(i))
 308                 }
 309         case v == a && v != b: // v <- v + alpha * b
 310                 if v.mat.Inc == 1 && bmat.Inc == 1 {
 311                         // Fast path for a common case.
 312                         f64.AxpyUnitaryTo(v.mat.Data, alpha, bmat.Data, amat.Data)
 313                 } else {
 314                         f64.AxpyInc(alpha, bmat.Data, v.mat.Data,
 315                                 uintptr(ar), uintptr(bmat.Inc), uintptr(v.mat.Inc), 0, 0)
 316                 }
 317         default: // v <- a + alpha * b or v <- a + alpha * v
 318                 if v.mat.Inc == 1 && amat.Inc == 1 && bmat.Inc == 1 {
 319                         // Fast path for a common case.
 320                         f64.AxpyUnitaryTo(v.mat.Data, alpha, bmat.Data, amat.Data)
 321                 } else {
 322                         f64.AxpyIncTo(v.mat.Data, uintptr(v.mat.Inc), 0,
 323                                 alpha, bmat.Data, amat.Data,
 324                                 uintptr(ar), uintptr(bmat.Inc), uintptr(amat.Inc), 0, 0)
 325                 }
 326         }
 327 }
 328
 329 // AddVec adds the vectors a and b, placing the result in the receiver.
 330 func (v *VecDense) AddVec(a, b Vector) {
 331         ar := a.Len()
 332         br := b.Len()
 333
 334         if ar != br {
 335                 panic(ErrShape)
 336         }
 337
 338         v.reuseAs(ar)
 339
 340         aU, _ := untranspose(a)
 341         bU, _ := untranspose(b)
 342
 343         if arv, ok := aU.(RawVectorer); ok {
 344                 if brv, ok := bU.(RawVectorer); ok {
 345                         amat := arv.RawVector()
 346                         bmat := brv.RawVector()
 347
 348                         if v != a {
 349                                 v.checkOverlap(amat)
 350                         }
 351                         if v != b {
 352                                 v.checkOverlap(bmat)
 353                         }
 354
 355                         if v.mat.Inc == 1 && amat.Inc == 1 && bmat.Inc == 1 {
 356                                 // Fast path for a common case.
 357                                 f64.AxpyUnitaryTo(v.mat.Data, 1, bmat.Data, amat.Data)
 358                                 return
 359                         }
 360                         f64.AxpyIncTo(v.mat.Data, uintptr(v.mat.Inc), 0,
 361                                 1, bmat.Data, amat.Data,
 362                                 uintptr(ar), uintptr(bmat.Inc), uintptr(amat.Inc), 0, 0)
 363                         return
 364                 }
 365         }
 366
 367         for i := 0; i < ar; i++ {
 368                 v.setVec(i, a.AtVec(i)+b.AtVec(i))
 369         }
 370 }
 371
 372 // SubVec subtracts the vector b from a, placing the result in the receiver.
 373 func (v *VecDense) SubVec(a, b Vector) {
 374         ar := a.Len()
 375         br := b.Len()
 376
 377         if ar != br {
 378                 panic(ErrShape)
 379         }
 380
 381         v.reuseAs(ar)
 382
 383         aU, _ := untranspose(a)
 384         bU, _ := untranspose(b)
 385
 386         if arv, ok := aU.(RawVectorer); ok {
 387                 if brv, ok := bU.(RawVectorer); ok {
 388                         amat := arv.RawVector()
 389                         bmat := brv.RawVector()
 390
 391                         if v != a {
 392                                 v.checkOverlap(amat)
 393                         }
 394                         if v != b {
 395                                 v.checkOverlap(bmat)
 396                         }
 397
 398                         if v.mat.Inc == 1 && amat.Inc == 1 && bmat.Inc == 1 {
 399                                 // Fast path for a common case.
 400                                 f64.AxpyUnitaryTo(v.mat.Data, -1, bmat.Data, amat.Data)
 401                                 return
 402                         }
 403                         f64.AxpyIncTo(v.mat.Data, uintptr(v.mat.Inc), 0,
 404                                 -1, bmat.Data, amat.Data,
 405                                 uintptr(ar), uintptr(bmat.Inc), uintptr(amat.Inc), 0, 0)
 406                         return
 407                 }
 408         }
 409
 410         for i := 0; i < ar; i++ {
 411                 v.setVec(i, a.AtVec(i)-b.AtVec(i))
 412         }
 413 }
 414
 415 // MulElemVec performs element-wise multiplication of a and b, placing the result
 416 // in the receiver.
 417 func (v *VecDense) MulElemVec(a, b Vector) {
 418         ar := a.Len()
 419         br := b.Len()
 420
 421         if ar != br {
 422                 panic(ErrShape)
 423         }
 424
 425         v.reuseAs(ar)
 426
 427         aU, _ := untranspose(a)
 428         bU, _ := untranspose(b)
 429
 430         if arv, ok := aU.(RawVectorer); ok {
 431                 if brv, ok := bU.(RawVectorer); ok {
 432                         amat := arv.RawVector()
 433                         bmat := brv.RawVector()
 434
 435                         if v != a {
 436                                 v.checkOverlap(amat)
 437                         }
 438                         if v != b {
 439                                 v.checkOverlap(bmat)
 440                         }
 441
 442                         if v.mat.Inc == 1 && amat.Inc == 1 && bmat.Inc == 1 {
 443                                 // Fast path for a common case.
 444                                 for i, a := range amat.Data {
 445                                         v.mat.Data[i] = a * bmat.Data[i]
 446                                 }
 447                                 return
 448                         }
 449                         var ia, ib int
 450                         for i := 0; i < ar; i++ {
 451                                 v.setVec(i, amat.Data[ia]*bmat.Data[ib])
 452                                 ia += amat.Inc
 453                                 ib += bmat.Inc
 454                         }
 455                         return
 456                 }
 457         }
 458
 459         for i := 0; i < ar; i++ {
 460                 v.setVec(i, a.AtVec(i)*b.AtVec(i))
 461         }
 462 }
 463
 464 // DivElemVec performs element-wise division of a by b, placing the result
 465 // in the receiver.
 466 func (v *VecDense) DivElemVec(a, b Vector) {
 467         ar := a.Len()
 468         br := b.Len()
 469
 470         if ar != br {
 471                 panic(ErrShape)
 472         }
 473
 474         v.reuseAs(ar)
 475
 476         aU, _ := untranspose(a)
 477         bU, _ := untranspose(b)
 478
 479         if arv, ok := aU.(RawVectorer); ok {
 480                 if brv, ok := bU.(RawVectorer); ok {
 481                         amat := arv.RawVector()
 482                         bmat := brv.RawVector()
 483
 484                         if v != a {
 485                                 v.checkOverlap(amat)
 486                         }
 487                         if v != b {
 488                                 v.checkOverlap(bmat)
 489                         }
 490
 491                         if v.mat.Inc == 1 && amat.Inc == 1 && bmat.Inc == 1 {
 492                                 // Fast path for a common case.
 493                                 for i, a := range amat.Data {
 494                                         v.setVec(i, a/bmat.Data[i])
 495                                 }
 496                                 return
 497                         }
 498                         var ia, ib int
 499                         for i := 0; i < ar; i++ {
 500                                 v.setVec(i, amat.Data[ia]/bmat.Data[ib])
 501                                 ia += amat.Inc
 502                                 ib += bmat.Inc
 503                         }
 504                 }
 505         }
 506
 507         for i := 0; i < ar; i++ {
 508                 v.setVec(i, a.AtVec(i)/b.AtVec(i))
 509         }
 510 }
 511
 512 // MulVec computes a * b. The result is stored into the receiver.
 513 // MulVec panics if the number of columns in a does not equal the number of rows in b
 514 // or if the number of columns in b does not equal 1.
 515 func (v *VecDense) MulVec(a Matrix, b Vector) {
 516         r, c := a.Dims()
 517         br, bc := b.Dims()
 518         if c != br || bc != 1 {
 519                 panic(ErrShape)
 520         }
 521
 522         aU, trans := untranspose(a)
 523         var bmat blas64.Vector
 524         fast := true
 525         bU, _ := untranspose(b)
 526         if rv, ok := bU.(RawVectorer); ok {
 527                 bmat = rv.RawVector()
 528                 if v != b {
 529                         v.checkOverlap(bmat)
 530                 }
 531         } else {
 532                 fast = false
 533         }
 534
 535         v.reuseAs(r)
 536         var restore func()
 537         if v == aU {
 538                 v, restore = v.isolatedWorkspace(aU.(*VecDense))
 539                 defer restore()
 540         } else if v == b {
 541                 v, restore = v.isolatedWorkspace(b)
 542                 defer restore()
 543         }
 544
 545         // TODO(kortschak): Improve the non-fast paths.
 546         switch aU := aU.(type) {
 547         case Vector:
 548                 if b.Len() == 1 {
 549                         // {n,1} x {1,1}
 550                         v.ScaleVec(b.AtVec(0), aU)
 551                         return
 552                 }
 553
 554                 // {1,n} x {n,1}
 555                 if fast {
 556                         if rv, ok := aU.(RawVectorer); ok {
 557                                 amat := rv.RawVector()
 558                                 if v != aU {
 559                                         v.checkOverlap(amat)
 560                                 }
 561
 562                                 if amat.Inc == 1 && bmat.Inc == 1 {
 563                                         // Fast path for a common case.
 564                                         v.setVec(0, f64.DotUnitary(amat.Data, bmat.Data))
 565                                         return
 566                                 }
 567                                 v.setVec(0, f64.DotInc(amat.Data, bmat.Data,
 568                                         uintptr(c), uintptr(amat.Inc), uintptr(bmat.Inc), 0, 0))
 569                                 return
 570                         }
 571                 }
 572                 var sum float64
 573                 for i := 0; i < c; i++ {
 574                         sum += aU.AtVec(i) * b.AtVec(i)
 575                 }
 576                 v.setVec(0, sum)
 577                 return
 578         case RawSymmetricer:
 579                 if fast {
 580                         amat := aU.RawSymmetric()
 581                         // We don't know that a is a *SymDense, so make
 582                         // a temporary SymDense to check overlap.
 583                         (&SymDense{mat: amat}).checkOverlap(v.asGeneral())
 584                         blas64.Symv(1, amat, bmat, 0, v.mat)
 585                         return
 586                 }
 587         case RawTriangular:
 588                 v.CopyVec(b)
 589                 amat := aU.RawTriangular()
 590                 // We don't know that a is a *TriDense, so make
 591                 // a temporary TriDense to check overlap.
 592                 (&TriDense{mat: amat}).checkOverlap(v.asGeneral())
 593                 ta := blas.NoTrans
 594                 if trans {
 595                         ta = blas.Trans
 596                 }
 597                 blas64.Trmv(ta, amat, v.mat)
 598         case RawMatrixer:
 599                 if fast {
 600                         amat := aU.RawMatrix()
 601                         // We don't know that a is a *Dense, so make
 602                         // a temporary Dense to check overlap.
 603                         (&Dense{mat: amat}).checkOverlap(v.asGeneral())
 604                         t := blas.NoTrans
 605                         if trans {
 606                                 t = blas.Trans
 607                         }
 608                         blas64.Gemv(t, 1, amat, bmat, 0, v.mat)
 609                         return
 610                 }
 611         default:
 612                 if fast {
 613                         for i := 0; i < r; i++ {
 614                                 var f float64
 615                                 for j := 0; j < c; j++ {
 616                                         f += a.At(i, j) * bmat.Data[j*bmat.Inc]
 617                                 }
 618                                 v.setVec(i, f)
 619                         }
 620                         return
 621                 }
 622         }
 623
 624         for i := 0; i < r; i++ {
 625                 var f float64
 626                 for j := 0; j < c; j++ {
 627                         f += a.At(i, j) * b.AtVec(j)
 628                 }
 629                 v.setVec(i, f)
 630         }
 631 }
 632
 633 // reuseAs resizes an empty vector to a r×1 vector,
 634 // or checks that a non-empty matrix is r×1.
 635 func (v *VecDense) reuseAs(r int) {
 636         if v.IsZero() {
 637                 v.mat = blas64.Vector{
 638                         Inc:  1,
 639                         Data: use(v.mat.Data, r),
 640                 }
 641                 v.n = r
 642                 return
 643         }
 644         if r != v.n {
 645                 panic(ErrShape)
 646         }
 647 }
 648
 649 // IsZero returns whether the receiver is zero-sized. Zero-sized vectors can be the
 650 // receiver for size-restricted operations. VecDenses can be zeroed using Reset.
 651 func (v *VecDense) IsZero() bool {
 652         // It must be the case that v.Dims() returns
 653         // zeros in this case. See comment in Reset().
 654         return v.mat.Inc == 0
 655 }
 656
 657 func (v *VecDense) isolatedWorkspace(a Vector) (n *VecDense, restore func()) {
 658         l := a.Len()
 659         n = getWorkspaceVec(l, false)
 660         return n, func() {
 661                 v.CopyVec(n)
 662                 putWorkspaceVec(n)
 663         }
 664 }
 665
 666 // asDense returns a Dense representation of the receiver with the same
 667 // underlying data.
 668 func (v *VecDense) asDense() *Dense {
 669         return &Dense{
 670                 mat:     v.asGeneral(),
 671                 capRows: v.n,
 672                 capCols: 1,
 673         }
 674 }
 675
 676 // asGeneral returns a blas64.General representation of the receiver with the
 677 // same underlying data.
 678 func (v *VecDense) asGeneral() blas64.General {
 679         return blas64.General{
 680                 Rows:   v.n,
 681                 Cols:   1,
 682                 Stride: v.mat.Inc,
 683                 Data:   v.mat.Data,
 684         }
 685 }
 686
 687 // ColViewOf reflects the column j of the RawMatrixer m, into the receiver
 688 // backed by the same underlying data. The length of the receiver must either be
 689 // zero or match the number of rows in m.
 690 func (v *VecDense) ColViewOf(m RawMatrixer, j int) {
 691         rm := m.RawMatrix()
 692
 693         if j >= rm.Cols || j < 0 {
 694                 panic(ErrColAccess)
 695         }
 696         if !v.IsZero() && v.n != rm.Rows {
 697                 panic(ErrShape)
 698         }
 699
 700         v.mat.Inc = rm.Stride
 701         v.mat.Data = rm.Data[j : (rm.Rows-1)*rm.Stride+j+1]
 702         v.n = rm.Rows
 703 }
 704
 705 // RowViewOf reflects the row i of the RawMatrixer m, into the receiver
 706 // backed by the same underlying data. The length of the receiver must either be
 707 // zero or match the number of columns in m.
 708 func (v *VecDense) RowViewOf(m RawMatrixer, i int) {
 709         rm := m.RawMatrix()
 710
 711         if i >= rm.Rows || i < 0 {
 712                 panic(ErrRowAccess)
 713         }
 714         if !v.IsZero() && v.n != rm.Cols {
 715                 panic(ErrShape)
 716         }
 717
 718         v.mat.Inc = 1
 719         v.mat.Data = rm.Data[i*rm.Stride : i*rm.Stride+rm.Cols]
 720         v.n = rm.Cols
 721 }