vendor/gonum.org/v1/gonum/blas/blas32/blas32.go

   1 // Copyright ©2015 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 blas32
   6
   7 import (
   8         "gonum.org/v1/gonum/blas"
   9         "gonum.org/v1/gonum/blas/gonum"
  10 )
  11
  12 var blas32 blas.Float32 = gonum.Implementation{}
  13
  14 // Use sets the BLAS float32 implementation to be used by subsequent BLAS calls.
  15 // The default implementation is native.Implementation.
  16 func Use(b blas.Float32) {
  17         blas32 = b
  18 }
  19
  20 // Implementation returns the current BLAS float32 implementation.
  21 //
  22 // Implementation allows direct calls to the current the BLAS float32 implementation
  23 // giving finer control of parameters.
  24 func Implementation() blas.Float32 {
  25         return blas32
  26 }
  27
  28 // Vector represents a vector with an associated element increment.
  29 type Vector struct {
  30         Inc  int
  31         Data []float32
  32 }
  33
  34 // General represents a matrix using the conventional storage scheme.
  35 type General struct {
  36         Rows, Cols int
  37         Stride     int
  38         Data       []float32
  39 }
  40
  41 // Band represents a band matrix using the band storage scheme.
  42 type Band struct {
  43         Rows, Cols int
  44         KL, KU     int
  45         Stride     int
  46         Data       []float32
  47 }
  48
  49 // Triangular represents a triangular matrix using the conventional storage scheme.
  50 type Triangular struct {
  51         N      int
  52         Stride int
  53         Data   []float32
  54         Uplo   blas.Uplo
  55         Diag   blas.Diag
  56 }
  57
  58 // TriangularBand represents a triangular matrix using the band storage scheme.
  59 type TriangularBand struct {
  60         N, K   int
  61         Stride int
  62         Data   []float32
  63         Uplo   blas.Uplo
  64         Diag   blas.Diag
  65 }
  66
  67 // TriangularPacked represents a triangular matrix using the packed storage scheme.
  68 type TriangularPacked struct {
  69         N    int
  70         Data []float32
  71         Uplo blas.Uplo
  72         Diag blas.Diag
  73 }
  74
  75 // Symmetric represents a symmetric matrix using the conventional storage scheme.
  76 type Symmetric struct {
  77         N      int
  78         Stride int
  79         Data   []float32
  80         Uplo   blas.Uplo
  81 }
  82
  83 // SymmetricBand represents a symmetric matrix using the band storage scheme.
  84 type SymmetricBand struct {
  85         N, K   int
  86         Stride int
  87         Data   []float32
  88         Uplo   blas.Uplo
  89 }
  90
  91 // SymmetricPacked represents a symmetric matrix using the packed storage scheme.
  92 type SymmetricPacked struct {
  93         N    int
  94         Data []float32
  95         Uplo blas.Uplo
  96 }
  97
  98 // Level 1
  99
 100 const negInc = "blas32: negative vector increment"
 101
 102 // Dot computes the dot product of the two vectors:
 103 //  \sum_i x[i]*y[i].
 104 func Dot(n int, x, y Vector) float32 {
 105         return blas32.Sdot(n, x.Data, x.Inc, y.Data, y.Inc)
 106 }
 107
 108 // DDot computes the dot product of the two vectors:
 109 //  \sum_i x[i]*y[i].
 110 func DDot(n int, x, y Vector) float64 {
 111         return blas32.Dsdot(n, x.Data, x.Inc, y.Data, y.Inc)
 112 }
 113
 114 // SDDot computes the dot product of the two vectors adding a constant:
 115 //  alpha + \sum_i x[i]*y[i].
 116 func SDDot(n int, alpha float32, x, y Vector) float32 {
 117         return blas32.Sdsdot(n, alpha, x.Data, x.Inc, y.Data, y.Inc)
 118 }
 119
 120 // Nrm2 computes the Euclidean norm of the vector x:
 121 //  sqrt(\sum_i x[i]*x[i]).
 122 //
 123 // Nrm2 will panic if the vector increment is negative.
 124 func Nrm2(n int, x Vector) float32 {
 125         if x.Inc < 0 {
 126                 panic(negInc)
 127         }
 128         return blas32.Snrm2(n, x.Data, x.Inc)
 129 }
 130
 131 // Asum computes the sum of the absolute values of the elements of x:
 132 //  \sum_i |x[i]|.
 133 //
 134 // Asum will panic if the vector increment is negative.
 135 func Asum(n int, x Vector) float32 {
 136         if x.Inc < 0 {
 137                 panic(negInc)
 138         }
 139         return blas32.Sasum(n, x.Data, x.Inc)
 140 }
 141
 142 // Iamax returns the index of an element of x with the largest absolute value.
 143 // If there are multiple such indices the earliest is returned.
 144 // Iamax returns -1 if n == 0.
 145 //
 146 // Iamax will panic if the vector increment is negative.
 147 func Iamax(n int, x Vector) int {
 148         if x.Inc < 0 {
 149                 panic(negInc)
 150         }
 151         return blas32.Isamax(n, x.Data, x.Inc)
 152 }
 153
 154 // Swap exchanges the elements of the two vectors:
 155 //  x[i], y[i] = y[i], x[i] for all i.
 156 func Swap(n int, x, y Vector) {
 157         blas32.Sswap(n, x.Data, x.Inc, y.Data, y.Inc)
 158 }
 159
 160 // Copy copies the elements of x into the elements of y:
 161 //  y[i] = x[i] for all i.
 162 func Copy(n int, x, y Vector) {
 163         blas32.Scopy(n, x.Data, x.Inc, y.Data, y.Inc)
 164 }
 165
 166 // Axpy adds x scaled by alpha to y:
 167 //  y[i] += alpha*x[i] for all i.
 168 func Axpy(n int, alpha float32, x, y Vector) {
 169         blas32.Saxpy(n, alpha, x.Data, x.Inc, y.Data, y.Inc)
 170 }
 171
 172 // Rotg computes the parameters of a Givens plane rotation so that
 173 //  ⎡ c s⎤   ⎡a⎤   ⎡r⎤
 174 //  ⎣-s c⎦ * ⎣b⎦ = ⎣0⎦
 175 // where a and b are the Cartesian coordinates of a given point.
 176 // c, s, and r are defined as
 177 //  r = ±Sqrt(a^2 + b^2),
 178 //  c = a/r, the cosine of the rotation angle,
 179 //  s = a/r, the sine of the rotation angle,
 180 // and z is defined such that
 181 //  if |a| > |b|,        z = s,
 182 //  otherwise if c != 0, z = 1/c,
 183 //  otherwise            z = 1.
 184 func Rotg(a, b float32) (c, s, r, z float32) {
 185         return blas32.Srotg(a, b)
 186 }
 187
 188 // Rotmg computes the modified Givens rotation. See
 189 // http://www.netlib.org/lapack/explore-html/df/deb/drotmg_8f.html
 190 // for more details.
 191 func Rotmg(d1, d2, b1, b2 float32) (p blas.SrotmParams, rd1, rd2, rb1 float32) {
 192         return blas32.Srotmg(d1, d2, b1, b2)
 193 }
 194
 195 // Rot applies a plane transformation to n points represented by the vectors x
 196 // and y:
 197 //  x[i] =  c*x[i] + s*y[i],
 198 //  y[i] = -s*x[i] + c*y[i], for all i.
 199 func Rot(n int, x, y Vector, c, s float32) {
 200         blas32.Srot(n, x.Data, x.Inc, y.Data, y.Inc, c, s)
 201 }
 202
 203 // Rotm applies the modified Givens rotation to n points represented by the
 204 // vectors x and y.
 205 func Rotm(n int, x, y Vector, p blas.SrotmParams) {
 206         blas32.Srotm(n, x.Data, x.Inc, y.Data, y.Inc, p)
 207 }
 208
 209 // Scal scales the vector x by alpha:
 210 //  x[i] *= alpha for all i.
 211 //
 212 // Scal will panic if the vector increment is negative.
 213 func Scal(n int, alpha float32, x Vector) {
 214         if x.Inc < 0 {
 215                 panic(negInc)
 216         }
 217         blas32.Sscal(n, alpha, x.Data, x.Inc)
 218 }
 219
 220 // Level 2
 221
 222 // Gemv computes
 223 //  y = alpha * A * x + beta * y,   if t == blas.NoTrans,
 224 //  y = alpha * A^T * x + beta * y, if t == blas.Trans or blas.ConjTrans,
 225 // where A is an m×n dense matrix, x and y are vectors, and alpha and beta are scalars.
 226 func Gemv(t blas.Transpose, alpha float32, a General, x Vector, beta float32, y Vector) {
 227         blas32.Sgemv(t, a.Rows, a.Cols, alpha, a.Data, a.Stride, x.Data, x.Inc, beta, y.Data, y.Inc)
 228 }
 229
 230 // Gbmv computes
 231 //  y = alpha * A * x + beta * y,   if t == blas.NoTrans,
 232 //  y = alpha * A^T * x + beta * y, if t == blas.Trans or blas.ConjTrans,
 233 // where A is an m×n band matrix, x and y are vectors, and alpha and beta are scalars.
 234 func Gbmv(t blas.Transpose, alpha float32, a Band, x Vector, beta float32, y Vector) {
 235         blas32.Sgbmv(t, a.Rows, a.Cols, a.KL, a.KU, alpha, a.Data, a.Stride, x.Data, x.Inc, beta, y.Data, y.Inc)
 236 }
 237
 238 // Trmv computes
 239 //  x = A * x,   if t == blas.NoTrans,
 240 //  x = A^T * x, if t == blas.Trans or blas.ConjTrans,
 241 // where A is an n×n triangular matrix, and x is a vector.
 242 func Trmv(t blas.Transpose, a Triangular, x Vector) {
 243         blas32.Strmv(a.Uplo, t, a.Diag, a.N, a.Data, a.Stride, x.Data, x.Inc)
 244 }
 245
 246 // Tbmv computes
 247 //  x = A * x,   if t == blas.NoTrans,
 248 //  x = A^T * x, if t == blas.Trans or blas.ConjTrans,
 249 // where A is an n×n triangular band matrix, and x is a vector.
 250 func Tbmv(t blas.Transpose, a TriangularBand, x Vector) {
 251         blas32.Stbmv(a.Uplo, t, a.Diag, a.N, a.K, a.Data, a.Stride, x.Data, x.Inc)
 252 }
 253
 254 // Tpmv computes
 255 //  x = A * x,   if t == blas.NoTrans,
 256 //  x = A^T * x, if t == blas.Trans or blas.ConjTrans,
 257 // where A is an n×n triangular matrix in packed format, and x is a vector.
 258 func Tpmv(t blas.Transpose, a TriangularPacked, x Vector) {
 259         blas32.Stpmv(a.Uplo, t, a.Diag, a.N, a.Data, x.Data, x.Inc)
 260 }
 261
 262 // Trsv solves
 263 //  A * x = b,   if t == blas.NoTrans,
 264 //  A^T * x = b, if t == blas.Trans or blas.ConjTrans,
 265 // where A is an n×n triangular matrix, and x and b are vectors.
 266 //
 267 // At entry to the function, x contains the values of b, and the result is
 268 // stored in-place into x.
 269 //
 270 // No test for singularity or near-singularity is included in this
 271 // routine. Such tests must be performed before calling this routine.
 272 func Trsv(t blas.Transpose, a Triangular, x Vector) {
 273         blas32.Strsv(a.Uplo, t, a.Diag, a.N, a.Data, a.Stride, x.Data, x.Inc)
 274 }
 275
 276 // Tbsv solves
 277 //  A * x = b,   if t == blas.NoTrans,
 278 //  A^T * x = b, if t == blas.Trans or blas.ConjTrans,
 279 // where A is an n×n triangular band matrix, and x and b are vectors.
 280 //
 281 // At entry to the function, x contains the values of b, and the result is
 282 // stored in place into x.
 283 //
 284 // No test for singularity or near-singularity is included in this
 285 // routine. Such tests must be performed before calling this routine.
 286 func Tbsv(t blas.Transpose, a TriangularBand, x Vector) {
 287         blas32.Stbsv(a.Uplo, t, a.Diag, a.N, a.K, a.Data, a.Stride, x.Data, x.Inc)
 288 }
 289
 290 // Tpsv solves
 291 //  A * x = b,   if t == blas.NoTrans,
 292 //  A^T * x = b, if t == blas.Trans or blas.ConjTrans,
 293 // where A is an n×n triangular matrix in packed format, and x and b are
 294 // vectors.
 295 //
 296 // At entry to the function, x contains the values of b, and the result is
 297 // stored in place into x.
 298 //
 299 // No test for singularity or near-singularity is included in this
 300 // routine. Such tests must be performed before calling this routine.
 301 func Tpsv(t blas.Transpose, a TriangularPacked, x Vector) {
 302         blas32.Stpsv(a.Uplo, t, a.Diag, a.N, a.Data, x.Data, x.Inc)
 303 }
 304
 305 // Symv computes
 306 //    y = alpha * A * x + beta * y,
 307 // where A is an n×n symmetric matrix, x and y are vectors, and alpha and
 308 // beta are scalars.
 309 func Symv(alpha float32, a Symmetric, x Vector, beta float32, y Vector) {
 310         blas32.Ssymv(a.Uplo, a.N, alpha, a.Data, a.Stride, x.Data, x.Inc, beta, y.Data, y.Inc)
 311 }
 312
 313 // Sbmv performs
 314 //  y = alpha * A * x + beta * y,
 315 // where A is an n×n symmetric band matrix, x and y are vectors, and alpha
 316 // and beta are scalars.
 317 func Sbmv(alpha float32, a SymmetricBand, x Vector, beta float32, y Vector) {
 318         blas32.Ssbmv(a.Uplo, a.N, a.K, alpha, a.Data, a.Stride, x.Data, x.Inc, beta, y.Data, y.Inc)
 319 }
 320
 321 // Spmv performs
 322 //    y = alpha * A * x + beta * y,
 323 // where A is an n×n symmetric matrix in packed format, x and y are vectors,
 324 // and alpha and beta are scalars.
 325 func Spmv(alpha float32, a SymmetricPacked, x Vector, beta float32, y Vector) {
 326         blas32.Sspmv(a.Uplo, a.N, alpha, a.Data, x.Data, x.Inc, beta, y.Data, y.Inc)
 327 }
 328
 329 // Ger performs a rank-1 update
 330 //  A += alpha * x * y^T,
 331 // where A is an m×n dense matrix, x and y are vectors, and alpha is a scalar.
 332 func Ger(alpha float32, x, y Vector, a General) {
 333         blas32.Sger(a.Rows, a.Cols, alpha, x.Data, x.Inc, y.Data, y.Inc, a.Data, a.Stride)
 334 }
 335
 336 // Syr performs a rank-1 update
 337 //  A += alpha * x * x^T,
 338 // where A is an n×n symmetric matrix, x is a vector, and alpha is a scalar.
 339 func Syr(alpha float32, x Vector, a Symmetric) {
 340         blas32.Ssyr(a.Uplo, a.N, alpha, x.Data, x.Inc, a.Data, a.Stride)
 341 }
 342
 343 // Spr performs the rank-1 update
 344 //  A += alpha * x * x^T,
 345 // where A is an n×n symmetric matrix in packed format, x is a vector, and
 346 // alpha is a scalar.
 347 func Spr(alpha float32, x Vector, a SymmetricPacked) {
 348         blas32.Sspr(a.Uplo, a.N, alpha, x.Data, x.Inc, a.Data)
 349 }
 350
 351 // Syr2 performs a rank-2 update
 352 //  A += alpha * x * y^T + alpha * y * x^T,
 353 // where A is a symmetric n×n matrix, x and y are vectors, and alpha is a scalar.
 354 func Syr2(alpha float32, x, y Vector, a Symmetric) {
 355         blas32.Ssyr2(a.Uplo, a.N, alpha, x.Data, x.Inc, y.Data, y.Inc, a.Data, a.Stride)
 356 }
 357
 358 // Spr2 performs a rank-2 update
 359 //  A += alpha * x * y^T + alpha * y * x^T,
 360 // where A is an n×n symmetric matrix in packed format, x and y are vectors,
 361 // and alpha is a scalar.
 362 func Spr2(alpha float32, x, y Vector, a SymmetricPacked) {
 363         blas32.Sspr2(a.Uplo, a.N, alpha, x.Data, x.Inc, y.Data, y.Inc, a.Data)
 364 }
 365
 366 // Level 3
 367
 368 // Gemm computes
 369 //  C = alpha * A * B + beta * C,
 370 // where A, B, and C are dense matrices, and alpha and beta are scalars.
 371 // tA and tB specify whether A or B are transposed.
 372 func Gemm(tA, tB blas.Transpose, alpha float32, a, b General, beta float32, c General) {
 373         var m, n, k int
 374         if tA == blas.NoTrans {
 375                 m, k = a.Rows, a.Cols
 376         } else {
 377                 m, k = a.Cols, a.Rows
 378         }
 379         if tB == blas.NoTrans {
 380                 n = b.Cols
 381         } else {
 382                 n = b.Rows
 383         }
 384         blas32.Sgemm(tA, tB, m, n, k, alpha, a.Data, a.Stride, b.Data, b.Stride, beta, c.Data, c.Stride)
 385 }
 386
 387 // Symm performs
 388 //  C = alpha * A * B + beta * C, if s == blas.Left,
 389 //  C = alpha * B * A + beta * C, if s == blas.Right,
 390 // where A is an n×n or m×m symmetric matrix, B and C are m×n matrices, and
 391 // alpha is a scalar.
 392 func Symm(s blas.Side, alpha float32, a Symmetric, b General, beta float32, c General) {
 393         var m, n int
 394         if s == blas.Left {
 395                 m, n = a.N, b.Cols
 396         } else {
 397                 m, n = b.Rows, a.N
 398         }
 399         blas32.Ssymm(s, a.Uplo, m, n, alpha, a.Data, a.Stride, b.Data, b.Stride, beta, c.Data, c.Stride)
 400 }
 401
 402 // Syrk performs a symmetric rank-k update
 403 //  C = alpha * A * A^T + beta * C, if t == blas.NoTrans,
 404 //  C = alpha * A^T * A + beta * C, if t == blas.Trans or blas.ConjTrans,
 405 // where C is an n×n symmetric matrix, A is an n×k matrix if t == blas.NoTrans and
 406 // a k×n matrix otherwise, and alpha and beta are scalars.
 407 func Syrk(t blas.Transpose, alpha float32, a General, beta float32, c Symmetric) {
 408         var n, k int
 409         if t == blas.NoTrans {
 410                 n, k = a.Rows, a.Cols
 411         } else {
 412                 n, k = a.Cols, a.Rows
 413         }
 414         blas32.Ssyrk(c.Uplo, t, n, k, alpha, a.Data, a.Stride, beta, c.Data, c.Stride)
 415 }
 416
 417 // Syr2k performs a symmetric rank-2k update
 418 //  C = alpha * A * B^T + alpha * B * A^T + beta * C, if t == blas.NoTrans,
 419 //  C = alpha * A^T * B + alpha * B^T * A + beta * C, if t == blas.Trans or blas.ConjTrans,
 420 // where C is an n×n symmetric matrix, A and B are n×k matrices if t == NoTrans
 421 // and k×n matrices otherwise, and alpha and beta are scalars.
 422 func Syr2k(t blas.Transpose, alpha float32, a, b General, beta float32, c Symmetric) {
 423         var n, k int
 424         if t == blas.NoTrans {
 425                 n, k = a.Rows, a.Cols
 426         } else {
 427                 n, k = a.Cols, a.Rows
 428         }
 429         blas32.Ssyr2k(c.Uplo, t, n, k, alpha, a.Data, a.Stride, b.Data, b.Stride, beta, c.Data, c.Stride)
 430 }
 431
 432 // Trmm performs
 433 //  B = alpha * A * B,   if tA == blas.NoTrans and s == blas.Left,
 434 //  B = alpha * A^T * B, if tA == blas.Trans or blas.ConjTrans, and s == blas.Left,
 435 //  B = alpha * B * A,   if tA == blas.NoTrans and s == blas.Right,
 436 //  B = alpha * B * A^T, if tA == blas.Trans or blas.ConjTrans, and s == blas.Right,
 437 // where A is an n×n or m×m triangular matrix, B is an m×n matrix, and alpha is
 438 // a scalar.
 439 func Trmm(s blas.Side, tA blas.Transpose, alpha float32, a Triangular, b General) {
 440         blas32.Strmm(s, a.Uplo, tA, a.Diag, b.Rows, b.Cols, alpha, a.Data, a.Stride, b.Data, b.Stride)
 441 }
 442
 443 // Trsm solves
 444 //  A * X = alpha * B,   if tA == blas.NoTrans and s == blas.Left,
 445 //  A^T * X = alpha * B, if tA == blas.Trans or blas.ConjTrans, and s == blas.Left,
 446 //  X * A = alpha * B,   if tA == blas.NoTrans and s == blas.Right,
 447 //  X * A^T = alpha * B, if tA == blas.Trans or blas.ConjTrans, and s == blas.Right,
 448 // where A is an n×n or m×m triangular matrix, X and B are m×n matrices, and
 449 // alpha is a scalar.
 450 //
 451 // At entry to the function, X contains the values of B, and the result is
 452 // stored in-place into X.
 453 //
 454 // No check is made that A is invertible.
 455 func Trsm(s blas.Side, tA blas.Transpose, alpha float32, a Triangular, b General) {
 456         blas32.Strsm(s, a.Uplo, tA, a.Diag, b.Rows, b.Cols, alpha, a.Data, a.Stride, b.Data, b.Stride)
 457 }