Module: NMatrix::BLAS

Defined in:
lib/nmatrix/blas.rb

Overview

–

NMatrix

A linear algebra library for scientific computation in Ruby. NMatrix is part of SciRuby.

NMatrix was originally inspired by and derived from NArray, by Masahiro Tanaka: narray.rubyforge.org

Copyright Information

SciRuby is Copyright © 2010 - 2014, Ruby Science Foundation NMatrix is Copyright © 2012 - 2014, John Woods and the Ruby Science Foundation

Please see LICENSE.txt for additional copyright notices.

Contributing

By contributing source code to SciRuby, you agree to be bound by our Contributor Agreement:

blas.rb

This file contains the safer accessors for the BLAS functions supported by NMatrix. ++

Class Method Summary collapse

Class Method Details

.asum(x, incx = 1, n = nil) ⇒ Object

call-seq:

asum(x, incx, n) -> Numeric

Calculate the sum of absolute values of the entries of a vector x of size n

  • Arguments :

    • x -> an NMatrix (will also allow an NMatrix, but will treat it as if it’s a vector )

    • incx -> the skip size (defaults to 1)

    • n -> the size of x (defaults to x.size / incx)

  • Returns :

    • The sum

  • Raises :

    • ArgumentError -> Expected dense NMatrix for arg 0

    • RangeError -> n out of range

Raises:

  • (ArgumentError) 


260
261
262
263
264
265
 
# File 'lib/nmatrix/blas.rb', line 260

def asum(x, incx = 1, n = nil)
  n ||= x.size / incx
  raise(ArgumentError, "Expected dense NMatrix for arg 0") unless x.is_a?(NMatrix)
  raise(RangeError, "n out of range") if n*incx > x.size || n*incx <= 0 || n <= 0
  ::NMatrix::BLAS.cblas_asum(n, x, incx)
end

.cblas_herk(order, uplo, trans, n, k, alpha, a, lda, beta, c, ldc) ⇒ Object

Raises:

  • (NotImplementedError) 


301
302
303
 
# File 'lib/nmatrix/blas.rb', line 301

def cblas_herk(order, uplo, trans, n, k, alpha, a, lda, beta, c, ldc)
  raise(NotImplementedError,"cblas_herk requires either the nmatrix-lapacke or nmatrix-atlas gem")
end

.cblas_syrk(order, uplo, trans, n, k, alpha, a, lda, beta, c, ldc) ⇒ Object

Raises:

  • (NotImplementedError) 


297
298
299
 
# File 'lib/nmatrix/blas.rb', line 297

def cblas_syrk(order, uplo, trans, n, k, alpha, a, lda, beta, c, ldc)
  raise(NotImplementedError,"cblas_syrk requires either the nmatrix-lapacke or nmatrix-atlas gem")
end

.cblas_trmm(order, side, uplo, trans_a, diag, m, n, alpha, a, lda, b, ldb) ⇒ Object

The following are functions that used to be implemented in C, but now require nmatrix-atlas or nmatrix-lapcke to run properly, so we can just implemented their stubs in Ruby.

Raises:

  • (NotImplementedError) 


293
294
295
 
# File 'lib/nmatrix/blas.rb', line 293

def cblas_trmm(order, side, uplo, trans_a, diag, m, n, alpha, a, lda, b, ldb)
  raise(NotImplementedError,"cblas_trmm requires either the nmatrix-lapacke or nmatrix-atlas gem")
end

.gemm(a, b, c = nil, alpha = 1.0, beta = 0.0, transpose_a = false, transpose_b = false, m = nil, n = nil, k = nil, lda = nil, ldb = nil, ldc = nil) ⇒ Object

call-seq:

gemm(a, b) -> NMatrix
gemm(a, b, c) -> NMatrix
gemm(a, b, c, alpha, beta) -> NMatrix

Updates the value of C via the matrix multiplication

C = (alpha * A * B) + (beta * C)

where alpha and beta are scalar values.

  • Arguments :

    • a -> Matrix A.

    • b -> Matrix B.

    • c -> Matrix C.

    • alpha -> A scalar value that multiplies A * B.

    • beta -> A scalar value that multiplies C.

    • transpose_a ->

    • transpose_b ->

    • m ->

    • n ->

    • k ->

    • lda ->

    • ldb ->

    • ldc ->

  • Returns :

    • A NMatrix equal to (alpha * A * B) + (beta * C).

  • Raises :

    • ArgumentError -> a and b must be dense matrices.

    • ArgumentError -> c must be nil or a dense matrix.

    • ArgumentError -> The dtype of the matrices must be equal.

Raises:

  • (ArgumentError) 


71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
 
# File 'lib/nmatrix/blas.rb', line 71

def gemm(a, b, c = nil, alpha = 1.0, beta = 0.0, transpose_a = false, transpose_b = false, m = nil, n = nil, k = nil, lda = nil, ldb = nil, ldc = nil)
  raise(ArgumentError, 'Expected dense NMatrices as first two arguments.') unless a.is_a?(NMatrix) and b.is_a?(NMatrix) and a.stype == :dense and b.stype == :dense
  raise(ArgumentError, 'Expected nil or dense NMatrix as third argument.') unless c.nil? or (c.is_a?(NMatrix) and c.stype == :dense)
  raise(ArgumentError, 'NMatrix dtype mismatch.')              unless a.dtype == b.dtype and (c ? a.dtype == c.dtype : true)

  # First, set m, n, and k, which depend on whether we're taking the
  # transpose of a and b.
  if c
    m ||= c.shape[0]
    n ||= c.shape[1]
    k ||= transpose_a ? a.shape[0] : a.shape[1]

  else
    if transpose_a
      # Either :transpose or :complex_conjugate.
      m ||= a.shape[1]
      k ||= a.shape[0]

    else
      # No transpose.
      m ||= a.shape[0]
      k ||= a.shape[1]
    end

    n ||= transpose_b ? b.shape[0] : b.shape[1]
    c  = NMatrix.new([m, n], dtype: a.dtype)
  end

  # I think these are independent of whether or not a transpose occurs.
  lda ||= a.shape[1]
  ldb ||= b.shape[1]
  ldc ||= c.shape[1]

  # NM_COMPLEX64 and NM_COMPLEX128 both require complex alpha and beta.
  if a.dtype == :complex64 or a.dtype == :complex128
    alpha = Complex(1.0, 0.0) if alpha == 1.0
    beta  = Complex(0.0, 0.0) if beta  == 0.0
  end

  # For argument descriptions, see: http://www.netlib.org/blas/dgemm.f
  ::NMatrix::BLAS.cblas_gemm(:row, transpose_a, transpose_b, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)

  return c
end

.gemv(a, x, y = nil, alpha = 1.0, beta = 0.0, transpose_a = false, m = nil, n = nil, lda = nil, incx = nil, incy = nil) ⇒ Object

call-seq:

gemv(a, x) -> NMatrix
gemv(a, x, y) -> NMatrix
gemv(a, x, y, alpha, beta) -> NMatrix

Implements matrix-vector product via

y = (alpha * A * x) + (beta * y)

where alpha and beta are scalar values.

  • Arguments :

    • a -> Matrix A.

    • x -> Vector x.

    • y -> Vector y.

    • alpha -> A scalar value that multiplies A * x.

    • beta -> A scalar value that multiplies y.

    • transpose_a ->

    • m ->

    • n ->

    • lda ->

    • incx ->

    • incy ->

  • Returns : -

  • Raises :

    • ++ ->

Raises:

  • (ArgumentError) 


143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
 
# File 'lib/nmatrix/blas.rb', line 143

def gemv(a, x, y = nil, alpha = 1.0, beta = 0.0, transpose_a = false, m = nil, n = nil, lda = nil, incx = nil, incy = nil)
  raise(ArgumentError, 'Expected dense NMatrices as first two arguments.') unless a.is_a?(NMatrix) and x.is_a?(NMatrix) and a.stype == :dense and x.stype == :dense
  raise(ArgumentError, 'Expected nil or dense NMatrix as third argument.') unless y.nil? or (y.is_a?(NMatrix) and y.stype == :dense)
  raise(ArgumentError, 'NMatrix dtype mismatch.')              unless a.dtype == x.dtype and (y ? a.dtype == y.dtype : true)

  m ||= transpose_a == :transpose ? a.shape[1] : a.shape[0]
  n ||= transpose_a == :transpose ? a.shape[0] : a.shape[1]
  raise(ArgumentError, "dimensions don't match") unless x.shape[0] == n && x.shape[1] == 1

  if y
    raise(ArgumentError, "dimensions don't match") unless y.shape[0] == m && y.shape[1] == 1
  else
    y = NMatrix.new([m,1], dtype: a.dtype)
  end

  lda  ||= a.shape[1]
  incx ||= 1
  incy ||= 1

  ::NMatrix::BLAS.cblas_gemv(transpose_a, m, n, alpha, a, lda, x, incx, beta, y, incy)

  return y
end

.nrm2(x, incx = 1, n = nil) ⇒ Object

call-seq:

nrm2(x, incx, n)

Calculate the 2-norm of a vector x of size n

  • Arguments :

    • x -> an NMatrix (will also allow an NMatrix, but will treat it as if it’s a vector )

    • incx -> the skip size (defaults to 1)

    • n -> the size of x (defaults to x.size / incx)

  • Returns :

    • The 2-norm

  • Raises :

    • ArgumentError -> Expected dense NMatrix for arg 0

    • RangeError -> n out of range

Raises:

  • (ArgumentError) 


283
284
285
286
287
288
 
# File 'lib/nmatrix/blas.rb', line 283

def nrm2(x, incx = 1, n = nil)
  n ||= x.size / incx
  raise(ArgumentError, "Expected dense NMatrix for arg 0") unless x.is_a?(NMatrix)
  raise(RangeError, "n out of range") if n*incx > x.size || n*incx <= 0 || n <= 0
  ::NMatrix::BLAS.cblas_nrm2(n, x, incx)
end

.rot(x, y, c, s, incx = 1, incy = 1, n = nil, in_place = false) ⇒ Object

call-seq:

rot(x, y, c, s) -> [NMatrix, NMatrix]

Apply plane rotation.

  • Arguments :

    • x -> NMatrix

    • y -> NMatrix

    • c -> cosine of the angle of rotation

    • s -> sine of the angle of rotation

    • incx -> stride of NMatrix x

    • incy -> stride of NMatrix y

    • n -> number of elements to consider in x and y

    • in_place -> true if it’s okay to modify the supplied x and y parameters directly; false if not. Default is false.

  • Returns :

    • Array with the results, in the format [xx, yy]

  • Raises :

    • ArgumentError -> Expected dense NMatrices as first two arguments.

    • ArgumentError -> NMatrix dtype mismatch.

    • ArgumentError -> Need to supply n for non-standard incx, incy values.

Raises:

  • (ArgumentError) 


189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
 
# File 'lib/nmatrix/blas.rb', line 189

def rot(x, y, c, s, incx = 1, incy = 1, n = nil, in_place=false)
  raise(ArgumentError, 'Expected dense NMatrices as first two arguments.')    unless x.is_a?(NMatrix) and y.is_a?(NMatrix) and x.stype == :dense and y.stype == :dense
  raise(ArgumentError, 'NMatrix dtype mismatch.')                 unless x.dtype == y.dtype
  raise(ArgumentError, 'Need to supply n for non-standard incx, incy values') if n.nil? && incx != 1 && incx != -1 && incy != 1 && incy != -1

  n ||= [x.size/incx.abs, y.size/incy.abs].min

  if in_place
    ::NMatrix::BLAS.cblas_rot(n, x, incx, y, incy, c, s)
    return [x,y]
  else
    xx = x.clone
    yy = y.clone

    ::NMatrix::BLAS.cblas_rot(n, xx, incx, yy, incy, c, s)

    return [xx,yy]
  end
end

.rot!(x, y, c, s, incx = 1, incy = 1, n = nil) ⇒ Object

call-seq:

rot!(x, y, c, s) -> [NMatrix, NMatrix]

Apply plane rotation directly to x and y.

See rot for arguments.



217
218
219
 
# File 'lib/nmatrix/blas.rb', line 217

def rot!(x, y, c, s, incx = 1, incy = 1, n = nil)
  rot(x,y,c,s,incx,incy,n,true)
end

.rotg(ab) ⇒ Object

call-seq:

rotg(ab) -> [Numeric, Numeric]

Apply givens plane rotation to the coordinates (a,b), returning the cosine and sine of the angle theta.

Since the givens rotation includes a square root, integers are disallowed.

  • Arguments :

    • ab -> NMatrix with two elements

  • Returns :

    • Array with the results, in the format [cos(theta), sin(theta)]

  • Raises :

    • ArgumentError -> Expected dense NMatrix of size 2

Raises:

  • (ArgumentError) 


237
238
239
240
241
 
# File 'lib/nmatrix/blas.rb', line 237

def rotg(ab)
  raise(ArgumentError, "Expected dense NMatrix of shape [2,1] or [1,2]") unless ab.is_a?(NMatrix) && ab.stype == :dense && ab.size == 2

  ::NMatrix::BLAS.cblas_rotg(ab)
end