Class: NMatrix

Inherits:
Object show all
Includes:
Enumerable
Defined in:
lib/nmatrix/math.rb,
lib/nmatrix/nmatrix.rb,
lib/nmatrix/version.rb,
lib/nmatrix/enumerate.rb,
lib/nmatrix/shortcuts.rb,
lib/nmatrix/homogeneous.rb,
lib/nmatrix/lapack_core.rb,
lib/nmatrix/io/fortran_format.rb,
lib/nmatrix/io/harwell_boeing.rb,
ext/nmatrix/ruby_nmatrix.c

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:

io/matlab/fortran_format.rb

A parser for making sense of FORTRAN formats.

> Only handles R (real), F (float) and E (exponential) format codes.

++

Defined Under Namespace

Modules: BLAS, FactorizeLUMethods, IO, LAPACK, NMMath, VERSION, YaleFunctions

Class Method Summary collapse

Instance Method Summary collapse

Constructor Details

#initializeObject

Forward Declarations



34
# File 'ext/nmatrix/ruby_nmatrix.c', line 34

static VALUE nm_init(int argc, VALUE* argv, VALUE nm);

Dynamic Method Handling

This class handles dynamic methods through the method_missing method

#method_missing(name, *args, &block) ⇒ Object

:nodoc:



932
933
934
935
936
937
938
939
940
# File 'lib/nmatrix/nmatrix.rb', line 932

def method_missing name, *args, &block #:nodoc:
  if name.to_s =~ /^__list_elementwise_.*__$/
    raise NotImplementedError, "requested undefined list matrix element-wise operation"
  elsif name.to_s =~ /^__yale_scalar_.*__$/
    raise NotImplementedError, "requested undefined yale scalar element-wise operation"
  else
    super(name, *args, &block)
  end
end

Class Method Details

.[](*params) ⇒ Object

call-seq:

NMatrix[Numeric, ..., Numeric, dtype: Symbol] -> NMatrix
NMatrix[Array, dtype: Symbol] -> NMatrix

The default value for dtype is guessed from the first parameter. For example:

NMatrix[1.0, 2.0].dtype # => :float64

But this is just a guess. If the other values can’t be converted to this dtype, a TypeError will be raised.

You can use the N constant in this way:

N = NMatrix
N[1, 2, 3]

NMatrix needs to have a succinct way to create a matrix by specifying the components directly. This is very useful for using it as an advanced calculator, it is useful for learning how to use, for testing language features and for developing algorithms.

The NMatrix::[] method provides a way to create a matrix in a way that is compact and natural. The components are specified using Ruby array syntax. Optionally, one can specify a dtype as the last parameter (default is :float64).

Examples:

a = N[ 1,2,3,4 ]          =>  1  2  3  4

a = N[ 1,2,3,4, :int32 ]  =>  1  2  3  4

a = N[ [1,2,3], [3,4,5] ] =>  1.0  2.0  3.0
                              3.0  4.0  5.0

a = N[ 3,6,9 ].transpose => 3
                            6
                            9

SYNTAX COMPARISON:

MATLAB:  a = [ [1 2 3] ; [4 5 6] ]   or  [ 1 2 3 ; 4 5 6 ]
IDL:   a = [ [1,2,3] , [4,5,6] ]
NumPy:  a = array( [1,2,3], [4,5,6] )

SciRuby:      a = NMatrix[ [1,2,3], [4,5,6] ]
Ruby array:   a =  [ [1,2,3], [4,5,6] ]


100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
# File 'lib/nmatrix/shortcuts.rb', line 100

def [](*params)
  options = params.last.is_a?(Hash) ? params.pop : {}

  # First find the dimensions of the array.
  i = 0
  shape = []
  row = params
  while row.is_a?(Array)
    shape[i] = row.length
    row = row[0]
    i += 1
  end

  # A row vector should be stored as 1xN, not N
  #shape.unshift(1) if shape.size == 1

  # Then flatten the array.
  NMatrix.new(shape, params.flatten, options)
end

.block_diagonal(*params) ⇒ Object Also known as: block_diag

Generate a block-diagonal NMatrix from the supplied 2D square matrices.

  • Arguments

    • *params -> An array that collects all arguments passed to the method. The method

      can receive any number of arguments. Optionally, the last entry of +params+ is 
      a hash of options from NMatrix#initialize. All other entries of +params+ are 
      the blocks of the desired block-diagonal matrix. Each such matrix block can be 
      supplied as a square 2D NMatrix object, or alternatively as an array of arrays 
      (with dimensions corresponding to a square matrix), or alternatively as a number.
      
  • Returns

    • NMatrix of block-diagonal form filled with specified matrices as the blocks along the diagonal.

  • Example

a = NMatrix.new([2,2], [1,2,3,4])
b = NMatrix.new([1,1], [123], dtype: :float64)
c = Array.new(2) { [[10,10], [10,10]] }
d = Array[[1,2,3], [4,5,6], [7,8,9]]
m = NMatrix.block_diagonal(a, b, *c, d, 10.0, 11, dtype: :int64, stype: :yale)
      => 
      [
        [1, 2,   0,  0,  0,  0,  0, 0, 0, 0,  0,  0]
        [3, 4,   0,  0,  0,  0,  0, 0, 0, 0,  0,  0]
        [0, 0, 123,  0,  0,  0,  0, 0, 0, 0,  0,  0]
        [0, 0,   0, 10, 10,  0,  0, 0, 0, 0,  0,  0]
        [0, 0,   0, 10, 10,  0,  0, 0, 0, 0,  0,  0]
        [0, 0,   0,  0,  0, 10, 10, 0, 0, 0,  0,  0]
        [0, 0,   0,  0,  0, 10, 10, 0, 0, 0,  0,  0]
        [0, 0,   0,  0,  0,  0,  0, 1, 2, 3,  0,  0]
        [0, 0,   0,  0,  0,  0,  0, 4, 5, 6,  0,  0]
        [0, 0,   0,  0,  0,  0,  0, 7, 8, 9,  0,  0]
        [0, 0,   0,  0,  0,  0,  0, 0, 0, 0, 10,  0]
        [0, 0,   0,  0,  0,  0,  0, 0, 0, 0,  0, 11]
      ]


313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
# File 'lib/nmatrix/shortcuts.rb', line 313

def block_diagonal(*params)
  options = params.last.is_a?(Hash) ? params.pop : {}

  params.each_index do |i|
    params[i] = params[i].to_nm if params[i].is_a?(Array) # Convert Array to NMatrix
    params[i] = NMatrix.new([1,1], [params[i]]) if params[i].is_a?(Numeric) # Convert number to NMatrix
  end

  block_sizes = [] #holds the size of each matrix block
  params.each do |b|
    unless b.is_a?(NMatrix)
      raise(ArgumentError, "Only NMatrix or appropriate Array objects or single numbers allowed")
    end
    raise(ArgumentError, "Only 2D matrices or 2D arrays allowed") unless b.shape.size == 2
    raise(ArgumentError, "Only square-shaped blocks allowed") unless b.shape[0] == b.shape[1]
    block_sizes << b.shape[0]
  end

  block_diag_mat = NMatrix.zeros(block_sizes.sum, options)
  (0...params.length).each do |n|
    # First determine the size and position of the n'th block in the block-diagonal matrix
    block_size = block_sizes[n]
    block_pos = block_sizes[0...n].sum
    # populate the n'th block in the block-diagonal matrix
    (0...block_size).each do |i|
      (0...block_size).each do |j|
        block_diag_mat[block_pos+i,block_pos+j] = params[n][i,j]
      end
    end
  end

  return block_diag_mat
end

.diagonal(entries, opts = {}) ⇒ Object Also known as: diag, diagonals

call-seq:

diagonals(array) -> NMatrix
diagonals(array, dtype: dtype, stype: stype) -> NMatrix

Creates a matrix filled with specified diagonals.

  • Arguments :

    • entries -> Array containing input values for diagonal matrix

    • options -> (optional) Hash with options for NMatrix#initialize

  • Returns :

    • NMatrix filled with specified diagonal values.

Examples:

NMatrix.diagonal([1.0,2,3,4]) # => 1.0 0.0 0.0 0.0
                                   0.0 2.0 0.0 0.0
                                   0.0 0.0 3.0 0.0
                                   0.0 0.0 0.0 4.0

NMatrix.diagonal([1,2,3,4], dtype: :int32) # => 1 0 0 0
                                                0 2 0 0
                                                0 0 3 0
                                                0 0 0 4


265
266
267
268
269
270
271
272
273
# File 'lib/nmatrix/shortcuts.rb', line 265

def diagonal(entries, opts={})
  m = NMatrix.zeros(entries.size,
                    {:dtype => guess_dtype(entries[0]), :capacity => entries.size + 1}.merge(opts)
                   )
  entries.each_with_index do |n, i|
    m[i,i] = n
  end
  m
end

.eye(shape, opts = {}) ⇒ Object Also known as: identity

call-seq:

eye(shape) -> NMatrix
eye(shape, dtype: dtype) -> NMatrix
eye(shape, stype: stype, dtype: dtype) -> NMatrix

Creates an identity matrix (square matrix rank 2).

  • Arguments :

    • size -> Array (or integer for square matrix) specifying the dimensions.

    • dtype -> (optional) Default is :float64

    • stype -> (optional) Default is :dense.

  • Returns :

    • An identity matrix.

Examples:

NMatrix.eye(3) # =>   1.0   0.0   0.0
                      0.0   1.0   0.0
                      0.0   0.0   1.0

NMatrix.eye(3, dtype: :int32) # =>   1   0   0
                                     0   1   0
                                     0   0   1

NMatrix.eye(2, dtype: :int32, stype: :yale) # =>   1   0
                                                   0   1


229
230
231
232
233
234
235
236
237
# File 'lib/nmatrix/shortcuts.rb', line 229

def eye(shape, opts={})
  # Fill the diagonal with 1's.
  m = NMatrix.zeros(shape, {:dtype => :float64}.merge(opts))
  (0...m.shape[0]).each do |i|
    m[i, i] = 1
  end

  m
end

.guess_dtypeObject



65
# File 'ext/nmatrix/ruby_nmatrix.c', line 65

static VALUE nm_guess_dtype(VALUE self, VALUE v);

.load_matlab_file(file_path) ⇒ Object

call-seq:

load_matlab_file(path) -> Mat5Reader
  • Arguments :

    • file_path -> The path to a version 5 .mat file.

  • Returns :

    • A Mat5Reader object.



88
89
90
# File 'lib/nmatrix/nmatrix.rb', line 88

def load_matlab_file(file_path)
  NMatrix::IO::Mat5Reader.new(File.open(file_path, 'rb')).to_ruby
end

.load_pcd_file(file_path) ⇒ Object

call-seq:

load_pcd_file(path) -> PointCloudReader::MetaReader
  • Arguments :

    • file_path -> The path to a PCL PCD file.

  • Returns :

    • A PointCloudReader::MetaReader object with the matrix stored in its matrix property



99
100
101
# File 'lib/nmatrix/nmatrix.rb', line 99

def load_pcd_file(file_path)
  NMatrix::IO::PointCloudReader::MetaReader.new(file_path)
end

.meshgrid(vectors, options = {}) ⇒ Object

Make N-D coordinate arrays for vectorized evaluations of N-D scalar/vector fields over N-D grids, given N coordinate arrays arrs. N > 1.

call-seq:

meshgrid(arrs) -> Array of NMatrix
meshgrid(arrs, options) -> Array of NMatrix
  • Arguments :

    • vectors -> Array of N coordinate arrays (Array or NMatrix), if any have more than one dimension they will be flatten

    • options -> Hash with options (:sparse Boolean, false by default; :indexing Symbol, may be :ij or :xy, :xy by default)

  • Returns :

    • Array of N N-D NMatrixes

  • Examples :

    x, y = NMatrix::meshgrid([[1, [2, 3]], [4, 5]])
    x.to_a #<= [[1, 2, 3], [1, 2, 3]]
    y.to_a #<= [[4, 4, 4], [5, 5, 5]]
    
  • Using options :

    x, y = NMatrix::meshgrid([[[1, 2], 3], [4, 5]], sparse: true)
    x.to_a #<= [[1, 2, 3]]
    y.to_a #<= [[4], [5]]
    
    x, y = NMatrix::meshgrid([[1, 2, 3], [[4], 5]], indexing: :ij)
    x.to_a #<= [[1, 1], [2, 2], [3, 3]]
    y.to_a #<= [[4, 5], [4, 5], [4, 5]]
    

Raises:

  • (ArgumentError)


136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
# File 'lib/nmatrix/nmatrix.rb', line 136

def meshgrid(vectors, options = {})
  raise(ArgumentError, 'Expected at least 2 arrays.') if vectors.size < 2
  options[:indexing] ||= :xy
  raise(ArgumentError, 'Indexing must be :xy of :ij') unless [:ij, :xy].include? options[:indexing]
  mats = vectors.map { |arr| arr.respond_to?(:flatten) ? arr.flatten : arr.to_flat_array }
  mats[0], mats[1] = mats[1], mats[0] if options[:indexing] == :xy
  new_dim = mats.size
  lengths = mats.map(&:size)
  result = mats.map.with_index do |matrix, axis|
    if options[:sparse]
      new_shape = Array.new(new_dim, 1)
      new_shape[axis] = lengths[axis]
      new_elements = matrix
    else
      before_axis = lengths[0...axis].reduce(:*)
      after_axis = lengths[(axis+1)..-1].reduce(:*)
      new_shape = lengths
      new_elements = after_axis ? matrix.map{ |el| [el] * after_axis }.flatten : matrix
      new_elements *= before_axis if before_axis
    end
    NMatrix.new(new_shape, new_elements)
  end
  result[0], result[1] = result[1], result[0] if options[:indexing] == :xy
  result
end

.min_dtypeObject



66
# File 'ext/nmatrix/ruby_nmatrix.c', line 66

static VALUE nm_min_dtype(VALUE self, VALUE v);

.ones(shape, opts = {}) ⇒ Object

call-seq:

ones(shape) -> NMatrix
ones(shape, dtype: dtype, stype: stype) -> NMatrix

Creates a matrix filled with ones.

  • Arguments :

    • shape -> Array (or integer for square matrix) specifying the shape.

    • opts -> (optional) Hash of options from NMatrix#initialize

  • Returns :

    • NMatrix filled with ones.

Examples:

NMatrix.ones([1, 3]) # =>  1.0   1.0   1.0

NMatrix.ones([2, 3], dtype: :int32) # =>  1  1  1
                                          1  1  1


171
172
173
# File 'lib/nmatrix/shortcuts.rb', line 171

def ones(shape, opts={})
  NMatrix.new(shape, 1, {:dtype => :float64, :default => 1}.merge(opts))
end

.ones_like(nm) ⇒ NMatrix

call-seq:

ones_like(nm) -> NMatrix

Creates a new matrix of ones with the same dtype and shape as the provided matrix.

Parameters:

  • nm (NMatrix)

    the nmatrix whose dtype and shape will be used

Returns:

  • (NMatrix)

    a new nmatrix filled with ones.



184
185
186
# File 'lib/nmatrix/shortcuts.rb', line 184

def ones_like(nm)
  NMatrix.ones(nm.shape, dtype: nm.dtype, stype: nm.stype, capacity: nm.capacity, default: 1)
end

.random(shape, opts = {}) ⇒ Object Also known as: rand

call-seq:

random(shape) -> NMatrix

Creates a :dense NMatrix with random numbers between 0 and 1 generated by Random::rand. The parameter is the dimension of the matrix.

If you use an integer dtype, make sure to specify :scale as a parameter, or you’ll only get a matrix of 0s.

  • Arguments :

    • shape -> Array (or integer for square matrix) specifying the dimensions.

  • Returns :

    • NMatrix filled with random values.

Examples:

NMatrix.random([2, 2]) # => 0.4859439730644226   0.1783195585012436
                            0.23193766176700592  0.4503345191478729

NMatrix.random([2, 2], :dtype => :byte, :scale => 255) # => [ [252, 108] [44, 12] ]


370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
# File 'lib/nmatrix/shortcuts.rb', line 370

def random(shape, opts={})
  scale = opts.delete(:scale) || 1.0

  rng = Random.new

  random_values = []


  # Construct the values of the final matrix based on the dimension.
  if opts[:dtype] == :complex64 || opts[:dtype] == :complex128
    NMatrix.size(shape).times { |i| random_values << Complex(rng.rand(scale), rng.rand(scale)) }
  else
    NMatrix.size(shape).times { |i| random_values << rng.rand(scale) }
  end

  NMatrix.new(shape, random_values, {:dtype => :float64, :stype => :dense}.merge(opts))
end

.readObject



37
# File 'ext/nmatrix/ruby_nmatrix.c', line 37

static VALUE nm_read(int argc, VALUE* argv, VALUE self);

.seq(shape, options = {}) ⇒ Object

call-seq:

seq(shape) -> NMatrix
seq(shape, options) -> NMatrix
bindgen(shape) -> NMatrix of :byte
indgen(shape) -> NMatrix of :int64
findgen(shape) -> NMatrix of :float32
dindgen(shape) -> NMatrix of :float64
cindgen(shape) -> NMatrix of :complex64
zindgen(shape) -> NMatrix of :complex128
rbindgen(shape) -> NMatrix of :object

Creates a matrix filled with a sequence of integers starting at zero.

  • Arguments :

    • shape -> Array (or integer for square matrix) specifying the dimensions.

    • options -> (optional) Options permissible for NMatrix#initialize

  • Returns :

    • NMatrix filled with values 0 through size.

Examples:

NMatrix.seq(2) # =>   0   1
              2   3

NMatrix.seq([3, 3], dtype: :float32) # =>  0.0  1.0  2.0
                                    3.0  4.0  5.0
                                    6.0  7.0  8.0


418
419
420
421
422
423
424
425
# File 'lib/nmatrix/shortcuts.rb', line 418

def seq(shape, options={})

  # Construct the values of the final matrix based on the dimension.
  values = (0 ... NMatrix.size(shape)).to_a

  # It'll produce :int32, except if a dtype is provided.
  NMatrix.new(shape, values, {:stype => :dense}.merge(options))
end

.size(shape) ⇒ Object

Calculate the size of an NMatrix of a given shape.



104
105
106
107
# File 'lib/nmatrix/nmatrix.rb', line 104

def size(shape)
  shape = [shape,shape] unless shape.is_a?(Array)
  (0...shape.size).inject(1) { |x,i| x * shape[i] }
end

.translation(*args) ⇒ Object

call-seq:

translation(x, y, z) -> NMatrix
translation([x,y,z]) -> NMatrix
translation(translation_matrix) -> NMatrix
translation(translation_matrix) -> NMatrix
translation(translation, dtype: dtype) -> NMatrix
translation(x, y, z, dtype: dtype) -> NMatrix

Generate a 4x4 homogeneous transformation matrix representing a translation.

  • Returns :

    • A homogeneous transformation matrix consisting of a translation.

Examples:

NMatrix.translation(4.0,5.0,6.0) # =>
                                      1.0   0.0   0.0   4.0
                                      0.0   1.0   0.0   5.0
                                      0.0   0.0   1.0   6.0
                                      0.0   0.0   0.0   1.0

NMatrix.translation(4.0,5.0,6.0, dtype: :int64) # =>
                                                     1  0  0  4
                                                     0  1  0  5
                                                     0  0  1  6
                                                     0  0  0  1
NMatrix.translation(4,5,6) # =>
                                 1  0  0  4
                                 0  1  0  5
                                 0  0  1  6
                                 0  0  0  1


127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
# File 'lib/nmatrix/homogeneous.rb', line 127

def translation *args
  xyz = args.shift if args.first.is_a?(NMatrix) || args.first.is_a?(Array)
  default_dtype = xyz.respond_to?(:dtype) ? xyz.dtype : NMatrix.guess_dtype(xyz)
  opts = {dtype: default_dtype}
  opts = opts.merge(args.pop) if args.size > 0 && args.last.is_a?(Hash)
  xyz ||= args

  n = if args.size > 0
    NMatrix.eye(4, opts)
  else
    NMatrix.eye(4, opts)
  end
  n[0..2,3] = xyz
  n
end

.upcastObject

Singleton methods



170
# File 'ext/nmatrix/ruby_nmatrix.c', line 170

static VALUE nm_upcast(VALUE self, VALUE t1, VALUE t2);

.x_rotation(angle_in_radians, opts = {}) ⇒ Object

call-seq:

x_rotation(angle_in_radians) -> NMatrix
x_rotation(angle_in_radians, dtype: dtype) -> NMatrix
y_rotation(angle_in_radians) -> NMatrix
y_rotation(angle_in_radians, dtype: dtype) -> NMatrix
z_rotation(angle_in_radians) -> NMatrix
z_rotation(angle_in_radians, dtype: dtype) -> NMatrix

Generate a 4x4 homogeneous transformation matrix representing a rotation about the x, y, or z axis respectively.

  • Arguments :

    • angle_in_radians -> The angle of rotation in radians.

    • dtype -> (optional) Default is :float64

  • Returns :

    • A homogeneous transformation matrix consisting of a single rotation.

Examples:

NMatrix.x_rotation(Math::PI.quo(6)) # =>
                                          1.0      0.0       0.0       0.0
                                          0.0      0.866025 -0.499999  0.0
                                          0.0      0.499999  0.866025  0.0
                                          0.0      0.0       0.0       1.0

NMatrix.x_rotation(Math::PI.quo(6), dtype: :float32) # =>
                                          1.0      0.0       0.0       0.0
                                          0.0      0.866025 -0.5       0.0
                                          0.0      0.5       0.866025  0.0
                                          0.0      0.0       0.0       1.0


66
67
68
69
70
71
72
73
# File 'lib/nmatrix/homogeneous.rb', line 66

def x_rotation angle_in_radians, opts={}
  c = Math.cos(angle_in_radians)
  s = Math.sin(angle_in_radians)
  NMatrix.new(4, [1.0, 0.0, 0.0, 0.0,
                  0.0, c,   -s,  0.0,
                  0.0, s,    c,  0.0,
                  0.0, 0.0, 0.0, 1.0], {dtype: :float64}.merge(opts))
end

.y_rotation(angle_in_radians, opts = {}) ⇒ Object



75
76
77
78
79
80
81
82
# File 'lib/nmatrix/homogeneous.rb', line 75

def y_rotation angle_in_radians, opts={}
  c = Math.cos(angle_in_radians)
  s = Math.sin(angle_in_radians)
  NMatrix.new(4, [ c,  0.0,  s,  0.0,
                  0.0, 1.0, 0.0, 0.0,
                  -s,  0.0,  c,  0.0,
                  0.0, 0.0, 0.0, 1.0], {dtype: :float64}.merge(opts))
end

.z_rotation(angle_in_radians, opts = {}) ⇒ Object



84
85
86
87
88
89
90
91
# File 'lib/nmatrix/homogeneous.rb', line 84

def z_rotation angle_in_radians, opts={}
  c = Math.cos(angle_in_radians)
  s = Math.sin(angle_in_radians)
  NMatrix.new(4, [ c,  -s,  0.0, 0.0,
                   s,   c,  0.0, 0.0,
                  0.0, 0.0, 1.0, 0.0,
                  0.0, 0.0, 0.0, 1.0], {dtype: :float64}.merge(opts))
end

.zeros(shape, opts = {}) ⇒ Object Also known as: zeroes

call-seq:

zeros(shape) -> NMatrix
zeros(shape, dtype: dtype) -> NMatrix
zeros(shape, dtype: dtype, stype: stype) -> NMatrix

Creates a new matrix of zeros with the dimensions supplied as parameters.

  • Arguments :

    • shape -> Array (or integer for square matrix) specifying the dimensions.

    • dtype -> (optional) Default is :float64

    • stype -> (optional) Default is :dense.

  • Returns :

    • NMatrix filled with zeros.

Examples:

NMatrix.zeros(2) # =>  0.0   0.0
                       0.0   0.0

NMatrix.zeros([2, 3], dtype: :int32) # =>  0  0  0
                                           0  0  0

NMatrix.zeros([1, 5], dtype: :int32) # =>  0  0  0  0  0


146
147
148
# File 'lib/nmatrix/shortcuts.rb', line 146

def zeros(shape, opts = {})
  NMatrix.new(shape, 0, {:dtype => :float64}.merge(opts))
end

.zeros_like(nm) ⇒ NMatrix

call-seq:

zeros_like(nm) -> NMatrix

Creates a new matrix of zeros with the same stype, dtype, and shape as the provided matrix.

Parameters:

  • nm (NMatrix)

    the nmatrix whose stype, dtype, and shape will be used

Returns:

  • (NMatrix)

    a new nmatrix filled with zeros.



197
198
199
# File 'lib/nmatrix/shortcuts.rb', line 197

def zeros_like(nm)
  NMatrix.zeros(nm.shape, dtype: nm.dtype, stype: nm.stype, capacity: nm.capacity, default: 0)
end

Instance Method Details

#!~Object

#%Object

#*Object

#**Object

#+Object

#-Object

#-@Object

#/Object

#<Object

#<=Object

#==Object



152
# File 'ext/nmatrix/ruby_nmatrix.c', line 152

static VALUE nm_eqeq(VALUE left, VALUE right);

#=~Object

#>Object

#>=Object

#[]Object



60
# File 'ext/nmatrix/ruby_nmatrix.c', line 60

static VALUE nm_mref(int argc, VALUE* argv, VALUE self);

#[]=Object



58
# File 'ext/nmatrix/ruby_nmatrix.c', line 58

static VALUE nm_mset(int argc, VALUE* argv, VALUE self);

#__yale_ary__to_s(sym) ⇒ Object

:nodoc:



333
334
335
336
337
# File 'lib/nmatrix/nmatrix.rb', line 333

def __yale_ary__to_s(sym) #:nodoc:
  ary = self.send("__yale_#{sym.to_s}__".to_sym)

  '[' + ary.collect { |a| a ? a : 'nil'}.join(',') + ']'
end

#absObject

call-seq:

abs -> NMatrix

Maps all values in a matrix to their absolute values.



840
841
842
843
844
845
846
847
848
849
# File 'lib/nmatrix/math.rb', line 840

def abs
  if stype == :dense
    self.__dense_map__ { |v| v.abs }
  elsif stype == :list
    # FIXME: Need __list_map_stored__, but this will do for now.
    self.__list_map_merged_stored__(nil, nil) { |v,dummy| v.abs }
  else
    self.__yale_map_stored__ { |v| v.abs }
  end.cast(self.stype, abs_dtype)
end

#abs_dtypeObject

call-seq:

abs_dtype -> Symbol

Returns the dtype of the result of a call to #abs. In most cases, this is the same as dtype; it should only differ for :complex64 (where it’s :float32) and :complex128 (:float64).



824
825
826
827
828
829
830
831
832
# File 'lib/nmatrix/math.rb', line 824

def abs_dtype
  if self.dtype == :complex64
    :float32
  elsif self.dtype == :complex128
    :float64
  else
    self.dtype
  end
end

#acosObject

#acoshObject

#angle_vectorObject

call-seq:

angle_vector -> [angle, about_vector]

Find the angle vector for a quaternion. Assumes the quaternion has unit length.

Source: www.euclideanspace.com/maths/geometry/rotations/conversions/quaternionToAngle/

  • Returns :

    • An angle (in radians) describing the rotation about the about_vector.

    • A length-3 NMatrix representing the corresponding quaternion.

Examples:

q.angle_vector # => [1, 0, 0, 0]

Raises:



228
229
230
231
232
233
234
235
236
237
238
239
240
# File 'lib/nmatrix/homogeneous.rb', line 228

def angle_vector
  raise(ShapeError, "Expected length-4 vector or matrix (quaternion)") if self.shape[0] != 4
  raise("Expected unit quaternion") if self[0] > 1

  xyz = NMatrix.new([3], dtype: self.dtype)

  angle = 2 * Math.acos(self[0])
  s = Math.sqrt(1.0 - self[0]*self[0])

  xyz[0..2] = self[1..3]
  xyz /= s if s >= 0.001 # avoid divide by zero
  return [angle, xyz]
end

#asinObject

#asinhObject

#asum(incx = 1, n = nil) ⇒ Object Also known as: absolute_sum

call-seq:

absolute_sum -> Numeric

Arguments

- +incx+ -> the skip size (defaults to 1, no skip)
- +n+ -> the number of elements to include

Return the sum of the contents of the vector. This is the BLAS asum routine.



861
862
863
864
865
866
867
868
# File 'lib/nmatrix/math.rb', line 861

def asum incx=1, n=nil
  if self.shape == [1]
    return self[0].abs unless self.complex_dtype?
    return self[0].real.abs + self[0].imag.abs
  end
  return method_missing(:asum, incx, n) unless vector?
  NMatrix::BLAS::asum(self, incx, self.size / incx)
end

#atanObject

#atan2Object

#atanhObject

#binned_sorted_indicesObject

call-seq:

binned_sorted_indices -> Array

Returns an array of arrays of indices ordered by value sorted. Functions basically like sorted_indices, but groups indices together for those values that are the same.



915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
# File 'lib/nmatrix/nmatrix.rb', line 915

def binned_sorted_indices
  return method_missing(:sorted_indices) unless vector?
  ary = self.to_flat_array
  ary2 = []
  last_bin = ary.each_index.sort_by { |i| [ary[i]] }.inject([]) do |result, element|
    if result.empty? || ary[result[-1]] == ary[element]
      result << element
    else
      ary2 << result
      [element]
    end
  end
  ary2 << last_bin unless last_bin.empty?
  ary2
end

#capacityObject



50
# File 'ext/nmatrix/ruby_nmatrix.c', line 50

static VALUE nm_capacity(VALUE self);

#cast(*params) ⇒ Object

call-seq:

cast(stype, dtype, default) -> NMatrix
cast(stype, dtype) -> NMatrix
cast(stype) -> NMatrix
cast(options) -> NMatrix

This is a user-friendly helper for calling #cast_full. The easiest way to call this function is using an options hash, e.g.,

n.cast(:stype => :yale, :dtype => :int64, :default => false)

For list and yale, :default sets the “default value” or “init” of the matrix. List allows a bit more freedom since non-zeros are permitted. For yale, unpredictable behavior may result if the value is not false, nil, or some version of 0. Dense discards :default.

dtype and stype are inferred from the matrix upon which #cast is called – so you only really need to provide one. You can actually call this function with no arguments, in which case it functions like #clone.

If your dtype is :object and you are converting from :dense to a sparse type, it is recommended that you provide a :default, as 0 may behave differently from its Float or Complex equivalent. If no option is given, Fixnum 0 will be used.



226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
# File 'lib/nmatrix/nmatrix.rb', line 226

def cast(*params)
  if (params.size > 0 && params[0].is_a?(Hash))
    opts = {
        :stype => self.stype,
        :dtype => self.dtype,
        :default => self.stype == :dense ? 0 : self.default_value
    }.merge(params[0])

    self.cast_full(opts[:stype], opts[:dtype], opts[:default])
  else
    params << self.stype if params.size == 0
    params << self.dtype if params.size == 1
    #HACK: the default value can cause an exception if dtype is not complex
    #and default_value is. (The ruby C code apparently won't convert these.)
    #Perhaps this should be fixed in the C code (in rubyval_to_cval).
    default_value = maybe_get_noncomplex_default_value(params[1])
    params << (self.stype == :dense ? 0 : default_value) if params.size == 2
    self.cast_full(*params)
  end

end

#cast_fullObject

#cbrtObject

#ceilObject

#clone_structure(capacity = nil) ⇒ Object

call-seq:

clone_structure -> NMatrix

This function is like clone, but it only copies the structure and the default value. None of the other values are copied. It takes an optional capacity argument. This is mostly only useful for dense, where you may not want to initialize; for other types, you should probably use zeros_like.



991
992
993
994
995
# File 'lib/nmatrix/nmatrix.rb', line 991

def clone_structure(capacity = nil)
  opts = {stype: self.stype, default: self.default_value, dtype: self.dtype}
  opts = {capacity: capacity}.merge(opts) if self.yale?
  NMatrix.new(self.shape, opts)
end

#colsObject

call-seq:

cols -> Integer

This shortcut use #shape to return the number of columns (the second dimension) of the matrix.



267
268
269
# File 'lib/nmatrix/nmatrix.rb', line 267

def cols
  shape[1]
end

#column(column_number, get_by = :copy) ⇒ Object Also known as: col

call-seq:

column(column_number) -> NMatrix
column(column_number, get_by) -> NMatrix

Returns the column specified. Uses slicing by copy as default.

  • Arguments :

    • column_number -> Integer.

    • get_by -> Type of slicing to use, :copy or :reference.

  • Returns :

    • A NMatrix representing the requested column as a column vector.

Examples:

m = NMatrix.new(2, [1, 4, 9, 14], :int32) # =>  1   4
                                                9  14

m.column(1) # =>   4
                  14


516
517
518
# File 'lib/nmatrix/nmatrix.rb', line 516

def column(column_number, get_by = :copy)
  rank(1, column_number, get_by)
end

#complex_conjugate(new_stype = self.stype) ⇒ Object

call-seq:

complex_conjugate -> NMatrix
complex_conjugate(new_stype) -> NMatrix

Get the complex conjugate of this matrix. See also complex_conjugate! for an in-place operation (provided the dtype is already :complex64 or :complex128).

Doesn’t work on list matrices, but you can optionally pass in the stype you want to cast to if you’re dealing with a list matrix.

  • Arguments :

    • new_stype -> stype for the new matrix.

  • Returns :

    • If the original NMatrix isn’t complex, the result is a :complex128 NMatrix. Otherwise, it’s the original dtype.



552
553
554
# File 'lib/nmatrix/math.rb', line 552

def complex_conjugate(new_stype = self.stype)
  self.cast(new_stype, NMatrix::upcast(dtype, :complex64)).complex_conjugate!
end

#complex_conjugate!Object



161
# File 'ext/nmatrix/ruby_nmatrix.c', line 161

static VALUE nm_complex_conjugate_bang(VALUE self);

#complex_dtype?Boolean

call-seq:

complex_dtype?() -> Boolean

Checks if dtype is a complex type

Returns:

  • (Boolean)


364
365
366
# File 'lib/nmatrix/nmatrix.rb', line 364

def complex_dtype?
  [:complex64, :complex128].include?(self.dtype)
end

#concat(*matrices) ⇒ Object

call-seq:

matrix1.concat(*m2) -> NMatrix
matrix1.concat(*m2, rank) -> NMatrix
matrix1.hconcat(*m2) -> NMatrix
matrix1.vconcat(*m2) -> NMatrix
matrix1.dconcat(*m3) -> NMatrix

Joins two matrices together into a new larger matrix. Attempts to determine which direction to concatenate on by looking for the first common element of the matrix shape in reverse. In other words, concatenating two columns together without supplying rank will glue them into an n x 2 matrix.

You can also use hconcat, vconcat, and dconcat for the first three ranks. concat performs an hconcat when no rank argument is provided.

The two matrices must have the same dim.

  • Arguments :

    • matrices -> one or more matrices

    • rank -> Fixnum (for rank); alternatively, may use :row, :column, or

    :layer for 0, 1, 2, respectively



664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
# File 'lib/nmatrix/nmatrix.rb', line 664

def concat(*matrices)
  rank = nil
  rank = matrices.pop unless matrices.last.is_a?(NMatrix)

  # Find the first matching dimension and concatenate along that (unless rank is specified)
  if rank.nil?
    rank = self.dim-1
    self.shape.reverse_each.with_index do |s,i|
      matrices.each do |m|
        if m.shape[i] != s
          rank -= 1
          break
        end
      end
    end
  elsif rank.is_a?(Symbol) # Convert to numeric
    rank = {:row => 0, :column => 1, :col => 1, :lay => 2, :layer => 2}[rank]
  end

  # Need to figure out the new shape.
  new_shape = self.shape.dup
  new_shape[rank] = matrices.inject(self.shape[rank]) { |total,m| total + m.shape[rank] }

  # Now figure out the options for constructing the concatenated matrix.
  opts = {stype: self.stype, default: self.default_value, dtype: self.dtype}
  if self.yale?
    # We can generally predict the new capacity for Yale. Subtract out the number of rows
    # for each matrix being concatenated, and then add in the number of rows for the new
    # shape. That takes care of the diagonal. The rest of the capacity is represented by
    # the non-diagonal non-default values.
    new_cap = matrices.inject(self.capacity - self.shape[0]) do |total,m|
      total + m.capacity - m.shape[0]
    end - self.shape[0] + new_shape[0]
    opts = {capacity: new_cap}.merge(opts)
  end

  # Do the actual construction.
  n = NMatrix.new(new_shape, opts)

  # Figure out where to start and stop the concatenation. We'll use NMatrices instead of
  # Arrays because then we can do elementwise addition.
  ranges = self.shape.map.with_index { |s,i| 0...self.shape[i] }

  matrices.unshift(self)
  matrices.each do |m|
    n[*ranges] = m

    # move over by the requisite amount
    ranges[rank]  = (ranges[rank].first + m.shape[rank])...(ranges[rank].last + m.shape[rank])
  end

  n
end

#conjugate_transposeObject

call-seq:

conjugate_transpose -> NMatrix

Calculate the conjugate transpose of a matrix. If your dtype is already complex, this should only require one copy (for the transpose).

  • Returns :

    • The conjugate transpose of the matrix as a copy.



680
681
682
# File 'lib/nmatrix/math.rb', line 680

def conjugate_transpose
  self.transpose.complex_conjugate!
end

#corrObject

Calculate the correlation matrix.

Raises:

  • (NotImplementedError)


579
580
581
582
583
# File 'lib/nmatrix/math.rb', line 579

def corr
  raise NotImplementedError, "Does not work for complex dtypes" if complex_dtype?
  standard_deviation = std
  cov / (standard_deviation.transpose.dot(standard_deviation))
end

#cosObject

#coshObject

#cov(opts = {}) ⇒ Object

Calculate the variance co-variance matrix

Options

  • :for_sample_data - Default true. If set to false will consider the denominator for population data (i.e. N, as opposed to N-1 for sample data).

References

Raises:

  • (TypeError)


566
567
568
569
570
571
572
573
574
575
576
# File 'lib/nmatrix/math.rb', line 566

def cov(opts={})
  raise TypeError, "Only works for non-integer dtypes" if integer_dtype?
   opts = {
    for_sample_data: true
  }.merge(opts)
  
  denominator      = opts[:for_sample_data] ? rows - 1 : rows
  ones             = NMatrix.ones [rows,1] 
  deviation_scores = self - ones.dot(ones.transpose).dot(self) / rows
  deviation_scores.transpose.dot(deviation_scores) / denominator
end

#data_pointerObject

///////////////



68
# File 'ext/nmatrix/ruby_nmatrix.c', line 68

static VALUE nm_data_pointer(VALUE self);

#dconcat(*matrices) ⇒ Object

Depth concatenation with matrices.



729
730
731
# File 'lib/nmatrix/nmatrix.rb', line 729

def dconcat(*matrices)
  concat(*matrices, :layer)
end

#default_valueObject



43
# File 'ext/nmatrix/ruby_nmatrix.c', line 43

static VALUE nm_default_value(VALUE self);

#dense?Boolean

call-seq:

m.dense? -> true or false

Determine if m is a dense matrix.

Returns:

  • (Boolean)


41
# File 'lib/nmatrix/shortcuts.rb', line 41

def dense?; return stype == :dense; end

#detObject

call-seq:

det -> determinant

Calculate the determinant by way of LU decomposition. This is accomplished using clapack_getrf, and then by taking the product of the diagonal elements. There is a risk of underflow/overflow.

There are probably also more efficient ways to calculate the determinant. This method requires making a copy of the matrix, since clapack_getrf modifies its input.

For smaller matrices, you may be able to use #det_exact.

This function is guaranteed to return the same type of data in the matrix upon which it is called.

Integer matrices are converted to floating point matrices for the purposes of performing the calculation, as xGETRF can’t work on integer matrices.

  • Returns :

    • The determinant of the matrix. It’s the same type as the matrix’s dtype.

  • Raises :

    • ShapeError -> Must be used on square matrices.

Raises:



510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
# File 'lib/nmatrix/math.rb', line 510

def det
  raise(ShapeError, "determinant can be calculated only for square matrices") unless self.dim == 2 && self.shape[0] == self.shape[1]

  # Cast to a dtype for which getrf is implemented
  new_dtype = self.integer_dtype? ? :float64 : self.dtype
  copy = self.cast(:dense, new_dtype)

  # Need to know the number of permutations. We'll add up the diagonals of
  # the factorized matrix.
  pivot = copy.getrf!

  num_perm = 0 #number of permutations
  pivot.each_with_index do |swap, i|
    #pivot indexes rows starting from 1, instead of 0, so need to subtract 1 here
    num_perm += 1 if swap-1 != i
  end
  prod = num_perm % 2 == 1 ? -1 : 1 # odd permutations => negative
  [shape[0],shape[1]].min.times do |i|
    prod *= copy[i,i]
  end

  # Convert back to an integer if necessary
  new_dtype != self.dtype ? prod.round : prod #prevent rounding errors
end

#det_exactObject



157
# File 'ext/nmatrix/ruby_nmatrix.c', line 157

static VALUE nm_det_exact(VALUE self);

#diagonal(main_diagonal = true) ⇒ Object

Return the main diagonal or antidiagonal a matrix. Only works with 2D matrices.

Arguments

  • main_diagonal - Defaults to true. If passed ‘false’, then will return the antidiagonal of the matrix.

References



281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
# File 'lib/nmatrix/nmatrix.rb', line 281

def diagonal main_diagonal=true
  diag_size = [cols, rows].min
  diag = NMatrix.new [diag_size], dtype: dtype

  if main_diagonal
    0.upto(diag_size-1) do |i|
      diag[i] = self[i,i]
    end
  else
    row = 0
    (diag_size-1).downto(0) do |col|
      diag[row] = self[row,col]
      row += 1
    end
  end

  diag
end

#dimensionsObject Also known as: dim



46
# File 'ext/nmatrix/ruby_nmatrix.c', line 46

static VALUE nm_dim(VALUE self);

#dotObject

///////////////////////



156
# File 'ext/nmatrix/ruby_nmatrix.c', line 156

static VALUE nm_multiply(VALUE left_v, VALUE right_v);

#dtypeObject



41
# File 'ext/nmatrix/ruby_nmatrix.c', line 41

static VALUE nm_dtype(VALUE self);

#each(&bl) ⇒ Object

call-seq:

each -> Enumerator

Enumerate through the matrix. @see Enumerable#each

For dense, this actually calls a specialized each iterator (in C). For yale and list, it relies upon #each_with_indices (which is about as fast as reasonably possible for C code).



40
41
42
43
44
45
46
47
48
49
50
51
52
# File 'lib/nmatrix/enumerate.rb', line 40

def each &bl
  if self.stype == :dense
    self.__dense_each__(&bl)
  elsif block_given?
    self.each_with_indices(&bl)
  else # Handle case where no block is given
    Enumerator.new do |yielder|
      self.each_with_indices do |params|
        yielder.yield params
      end
    end
  end
end

#each_column(get_by = :reference) ⇒ Object

call-seq:

each_column { |column| block } -> NMatrix

Iterate through each column, referencing it as an NMatrix slice.



144
145
146
147
148
149
150
# File 'lib/nmatrix/enumerate.rb', line 144

def each_column(get_by=:reference)
  return enum_for(:each_column, get_by) unless block_given?
  (0...self.shape[1]).each do |j|
    yield self.column(j, get_by)
  end
  self
end

#each_layer(get_by = :reference) ⇒ Object

call-seq:

each_layer -> { |column| block } -> ...

Iterate through each layer, referencing it as an NMatrix slice.

Note: If you have a 3-dimensional matrix, the first dimension contains rows, the second contains columns, and the third contains layers.



160
161
162
163
164
165
166
# File 'lib/nmatrix/enumerate.rb', line 160

def each_layer(get_by=:reference)
  return enum_for(:each_layer, get_by) unless block_given?
  (0...self.shape[2]).each do |k|
    yield self.layer(k, get_by)
  end
  self
end

#each_ordered_stored_with_indicesObject



53
# File 'ext/nmatrix/ruby_nmatrix.c', line 53

static VALUE nm_each_ordered_stored_with_indices(VALUE nmatrix);

#each_rank(dimen = 0, get_by = :reference) ⇒ Object Also known as: each_along_dim

call-seq:

each_rank() -> NMatrix
each_rank() { |rank| block } -> NMatrix
each_rank(dimen) -> Enumerator
each_rank(dimen) { |rank| block } -> NMatrix

Generic for @each_row, @each_col

Iterate through each rank by reference.

Parameters:

  • dimen (Fixnum) (defaults to: 0)

    the rank being iterated over.



117
118
119
120
121
122
123
# File 'lib/nmatrix/enumerate.rb', line 117

def each_rank(dimen=0, get_by=:reference)
  return enum_for(:each_rank, dimen, get_by) unless block_given?
  (0...self.shape[dimen]).each do |idx|
    yield self.rank(dimen, idx, get_by)
  end
  self
end

#each_row(get_by = :reference) ⇒ Object

call-seq:

each_row { |row| block } -> NMatrix

Iterate through each row, referencing it as an NMatrix slice.



131
132
133
134
135
136
137
# File 'lib/nmatrix/enumerate.rb', line 131

def each_row(get_by=:reference)
  return enum_for(:each_row, get_by) unless block_given?
  (0...self.shape[0]).each do |i|
    yield self.row(i, get_by)
  end
  self
end

#each_stored_with_index(&block) ⇒ Object

call-seq:

each_stored_with_index -> Enumerator

Allow iteration across a vector NMatrix’s stored values. See also @each_stored_with_indices

Raises:

  • (NotImplementedError)


175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
# File 'lib/nmatrix/enumerate.rb', line 175

def each_stored_with_index(&block)
  raise(NotImplementedError, "only works for dim 2 vectors") unless self.dim <= 2
  return enum_for(:each_stored_with_index) unless block_given?

  self.each_stored_with_indices do |v, i, j|
    if shape[0] == 1
      yield(v,j)
    elsif shape[1] == 1
      yield(v,i)
    else
      method_missing(:each_stored_with_index, &block)
    end
  end
  self
end

#each_stored_with_indicesObject



52
# File 'ext/nmatrix/ruby_nmatrix.c', line 52

static VALUE nm_each_stored_with_indices(VALUE nmatrix);

#each_with_indicesObject

Iterators public methods



51
# File 'ext/nmatrix/ruby_nmatrix.c', line 51

static VALUE nm_each_with_indices(VALUE nmatrix);

#effective_dimensionsObject Also known as: effective_dim



45
# File 'ext/nmatrix/ruby_nmatrix.c', line 45

static VALUE nm_effective_dim(VALUE self);

#erfObject

#erfcObject

#expObject

#factorize_choleskyObject

call-seq:

factorize_cholesky -> [upper NMatrix, lower NMatrix]

Calculates the Cholesky factorization of a matrix and returns the upper and lower matrices such that A=LU and L=U*, where * is either the transpose or conjugate transpose.

Unlike potrf!, this makes method requires that the original is matrix is symmetric or Hermitian. However, it is still your responsibility to make sure it is positive-definite.



204
205
206
207
208
209
# File 'lib/nmatrix/math.rb', line 204

def factorize_cholesky
  raise "Matrix must be symmetric/Hermitian for Cholesky factorization" unless self.hermitian?
  l = self.clone.potrf_lower!.tril!
  u = l.conjugate_transpose
  [u,l]
end

#factorize_lu(with_permutation_matrix = nil) ⇒ Object

call-seq:

factorize_lu -> ...

LU factorization of a matrix. Optionally return the permutation matrix.

Note that computing the permutation matrix will introduce a slight memory
and time overhead.

Arguments

with_permutation_matrix - If set to true will return the permutation

matrix alongwith the LU factorization as a second return value.

Raises:

  • (NotImplementedError)


224
225
226
227
228
229
230
231
232
233
# File 'lib/nmatrix/math.rb', line 224

def factorize_lu with_permutation_matrix=nil
  raise(NotImplementedError, "only implemented for dense storage") unless self.stype == :dense
  raise(NotImplementedError, "matrix is not 2-dimensional") unless self.dimensions == 2

  t     = self.clone
  pivot = t.getrf!
  return t unless with_permutation_matrix

  [t, FactorizeLUMethods.permutation_matrix_from(pivot)]
end

#flat_mapObject

call-seq:

flat_map -> Enumerator
flat_map { |elem| block } -> Array

Maps using Enumerator (returns an Array or an Enumerator)



60
# File 'lib/nmatrix/enumerate.rb', line 60

alias_method :flat_map, :map

#float_dtype?Boolean

call-seq:

float_dtype?() -> Boolean

Checks if dtype is a floating point type

Returns:

  • (Boolean)


354
355
356
# File 'lib/nmatrix/nmatrix.rb', line 354

def float_dtype?
  [:float32, :float64].include?(dtype)
end

#floorObject

#gammaObject

#gesdd(workspace_size = nil) ⇒ Object

call-seq:

gesdd -> [u, sigma, v_transpose]
gesdd -> [u, sigma, v_conjugate_transpose] # complex

Compute the singular value decomposition of a matrix using LAPACK’s GESDD function. This uses a divide-and-conquer strategy. See also #gesvd.

Optionally accepts a workspace_size parameter, which will be honored only if it is larger than what LAPACK requires.



406
407
408
# File 'lib/nmatrix/math.rb', line 406

def gesdd(workspace_size=nil)
  self.clone.gesdd!(workspace_size)
end

#gesdd!(workspace_size = nil) ⇒ Object

call-seq:

gesdd! -> [u, sigma, v_transpose]
gesdd! -> [u, sigma, v_conjugate_transpose] # complex

Compute the singular value decomposition of a matrix using LAPACK’s GESDD function. This uses a divide-and-conquer strategy. This is destructive, modifying the source NMatrix. See also #gesvd.

Optionally accepts a workspace_size parameter, which will be honored only if it is larger than what LAPACK requires.



391
392
393
# File 'lib/nmatrix/math.rb', line 391

def gesdd!(workspace_size=nil)
  NMatrix::LAPACK::gesdd(self, workspace_size)
end

#gesvd(workspace_size = 1) ⇒ Object

call-seq:

gesvd -> [u, sigma, v_transpose]
gesvd -> [u, sigma, v_conjugate_transpose] # complex

Compute the singular value decomposition of a matrix using LAPACK’s GESVD function.

Optionally accepts a workspace_size parameter, which will be honored only if it is larger than what LAPACK requires.



374
375
376
# File 'lib/nmatrix/math.rb', line 374

def gesvd(workspace_size=1)
  self.clone.gesvd!(workspace_size)
end

#gesvd!(workspace_size = 1) ⇒ Object

call-seq:

gesvd! -> [u, sigma, v_transpose]
gesvd! -> [u, sigma, v_conjugate_transpose] # complex

Compute the singular value decomposition of a matrix using LAPACK’s GESVD function. This is destructive, modifying the source NMatrix. See also #gesdd.

Optionally accepts a workspace_size parameter, which will be honored only if it is larger than what LAPACK requires.



360
361
362
# File 'lib/nmatrix/math.rb', line 360

def gesvd!(workspace_size=1)
  NMatrix::LAPACK::gesvd(self, workspace_size)
end

#getrf!Object

call-seq:

getrf! -> Array

LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. The LU factorization is A = PLU, where P is a row permutation matrix, L is a lower triangular matrix with unit diagonals, and U is an upper triangular matrix (note that this convention is different from the clapack_getrf behavior, but matches the standard LAPACK getrf). A is overwritten with the elements of L and U (the unit diagonal elements of L are not saved). P is not returned directly and must be constructed from the pivot array ipiv. The row indices in ipiv are indexed starting from 1. Only works for dense matrices.

  • Returns :

    • The IPIV vector. The L and U matrices are stored in A.

  • Raises :

    • StorageTypeError -> ATLAS functions only work on dense matrices.

Raises:



135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
# File 'lib/nmatrix/math.rb', line 135

def getrf!
  raise(StorageTypeError, "ATLAS functions only work on dense matrices") unless self.dense?

  #For row-major matrices, clapack_getrf uses a different convention than
  #described above (U has unit diagonal elements instead of L and columns
  #are interchanged rather than rows). For column-major matrices, clapack
  #uses the stanard conventions. So we just transpose the matrix before
  #and after calling clapack_getrf.
  #Unfortunately, this is not a very good way, uses a lot of memory.
  temp = self.transpose
  ipiv = NMatrix::LAPACK::clapack_getrf(:col, self.shape[0], self.shape[1], temp, self.shape[0])
  temp = temp.transpose
  self[0...self.shape[0], 0...self.shape[1]] = temp

  #for some reason, in clapack_getrf, the indices in ipiv start from 0
  #instead of 1 as in LAPACK.
  ipiv.each_index { |i| ipiv[i]+=1 }

  return ipiv
end

#hconcat(*matrices) ⇒ Object

Horizontal concatenation with matrices.



719
720
721
# File 'lib/nmatrix/nmatrix.rb', line 719

def hconcat(*matrices)
  concat(*matrices, :column)
end

#hermitian?Boolean

Returns:

  • (Boolean)


150
# File 'ext/nmatrix/ruby_nmatrix.c', line 150

static VALUE nm_hermitian(VALUE self);

#hessenbergObject

Reduce self to upper hessenberg form using householder transforms.

References



241
242
243
# File 'lib/nmatrix/math.rb', line 241

def hessenberg
  clone.hessenberg!
end

#hessenberg!Object

Destructive version of #hessenberg

Raises:



246
247
248
249
250
251
252
253
254
255
256
257
# File 'lib/nmatrix/math.rb', line 246

def hessenberg!
  raise ShapeError, "Trying to reduce non 2D matrix to hessenberg form" if 
    shape.size != 2
  raise ShapeError, "Trying to reduce non-square matrix to hessenberg form" if 
    shape[0] != shape[1]
  raise StorageTypeError, "Matrix must be dense" if stype != :dense
  raise TypeError, "Works with float matrices only" unless 
    [:float64,:float32].include?(dtype)

  __hessenberg__(self)
  self
end

#hypotObject

#index(value) ⇒ Object

Returns the index of the first occurence of the specified value. Returns an array containing the position of the value, nil in case the value is not found.



968
969
970
971
972
973
974
975
976
977
978
979
980
# File 'lib/nmatrix/nmatrix.rb', line 968

def index(value)
  index = nil

  self.each_with_indices do |yields|
    if yields.first == value
      yields.shift
      index = yields
      break
    end
  end 

  index
end

#initialize_copyObject



35
# File 'ext/nmatrix/ruby_nmatrix.c', line 35

static VALUE nm_init_copy(VALUE copy, VALUE original);

#inject(sym) ⇒ Object

call-seq:

inject -> symbol

This overrides the inject function to use map_stored for yale matrices



960
961
962
963
# File 'lib/nmatrix/nmatrix.rb', line 960

def inject(sym)
  return super(sym) unless self.yale?
  return self.map_stored.inject(sym)
end

#inject_rank(dimen = 0, initial = nil, dtype = nil) ⇒ NMatrix Also known as: reduce_along_dim, inject_along_dim

call-seq:

inject_rank() -> Enumerator
inject_rank(dimen) -> Enumerator
inject_rank(dimen, initial) -> Enumerator
inject_rank(dimen, initial, dtype) -> Enumerator
inject_rank() { |elem| block } -> NMatrix
inject_rank(dimen) { |elem| block } -> NMatrix
inject_rank(dimen, initial) { |elem| block } -> NMatrix
inject_rank(dimen, initial, dtype) { |elem| block } -> NMatrix

Reduces an NMatrix using a supplied block over a specified dimension. The block should behave the same way as for Enumerable#reduce.

Parameters:

  • dimen (Integer) (defaults to: 0)

    the dimension being reduced

  • initial (Numeric) (defaults to: nil)

    the initial value for the reduction (i.e. the usual parameter to Enumerable#reduce). Supply nil or do not supply this argument to have it follow the usual Enumerable#reduce behavior of using the first element as the initial value.

  • dtype (Symbol) (defaults to: nil)

    if non-nil/false, forces the accumulated result to have this dtype

Returns:

  • (NMatrix)

    an NMatrix with the same number of dimensions as the input, but with the input dimension now having size 1. Each element is the result of the reduction at that position along the specified dimension.

Raises:

  • (RangeError)


217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
# File 'lib/nmatrix/enumerate.rb', line 217

def inject_rank(dimen=0, initial=nil, dtype=nil)

  raise(RangeError, "requested dimension (#{dimen}) does not exist (shape: #{shape})") if dimen > self.dim

  return enum_for(:inject_rank, dimen, initial, dtype) unless block_given?

  new_shape = shape
  new_shape[dimen] = 1

  first_as_acc = false

  if initial then
    acc = NMatrix.new(new_shape, initial, :dtype => dtype || self.dtype, stype: self.stype)
  else
    each_rank(dimen) do |sub_mat|
      acc = (sub_mat.is_a?(NMatrix) and !dtype.nil? and dtype != self.dtype) ? sub_mat.cast(self.stype, dtype) : sub_mat
      break
    end
    first_as_acc = true
  end

  each_rank(dimen) do |sub_mat|
    if first_as_acc
      first_as_acc = false
      next
    end
    acc = yield(acc, sub_mat)
  end

  acc
end

#inspectObject

:nodoc:



327
328
329
330
331
# File 'lib/nmatrix/nmatrix.rb', line 327

def inspect #:nodoc:
  original_inspect = super()
  original_inspect = original_inspect[0...original_inspect.size-1]
  original_inspect + " " + inspect_helper.join(" ") + ">"
end

#integer_dtype?Boolean

call-seq:

integer_dtype?() -> Boolean

Checks if dtype is an integer type

Returns:

  • (Boolean)


345
346
347
# File 'lib/nmatrix/nmatrix.rb', line 345

def integer_dtype?
  [:byte, :int8, :int16, :int32, :int64].include?(self.dtype)
end

#invertObject Also known as: inverse

call-seq:

invert -> NMatrix

Make a copy of the matrix, then invert using Gauss-Jordan elimination. Works without LAPACK.

  • Returns :

    • A dense NMatrix. Will be the same type as the input NMatrix,

    except if the input is an integral dtype, in which case it will be a :float64 NMatrix.

  • Raises :

    • StorageTypeError -> only implemented on dense matrices.

    • ShapeError -> matrix must be square.



102
103
104
105
106
107
108
109
110
111
112
# File 'lib/nmatrix/math.rb', line 102

def invert
  #write this in terms of invert! so plugins will only have to overwrite
  #invert! and not invert
  if self.integer_dtype?
    cloned = self.cast(dtype: :float64)
    cloned.invert!
  else
    cloned = self.clone
    cloned.invert!
  end
end

#invert!Object

call-seq:

invert! -> NMatrix

Use LAPACK to calculate the inverse of the matrix (in-place) if available. Only works on dense matrices. Alternatively uses in-place Gauss-Jordan elimination.

  • Raises :

    • StorageTypeError -> only implemented on dense matrices.

    • ShapeError -> matrix must be square.

    • DataTypeError -> cannot invert an integer matrix in-place.

Raises:



77
78
79
80
81
82
83
84
# File 'lib/nmatrix/math.rb', line 77

def invert!
  raise(StorageTypeError, "invert only works on dense matrices currently") unless self.dense?
  raise(ShapeError, "Cannot invert non-square matrix") unless self.dim == 2 && self.shape[0] == self.shape[1]
  raise(DataTypeError, "Cannot invert an integer matrix in-place") if self.integer_dtype?

  #No internal implementation of getri, so use this other function
  __inverse__(self, true)
end

#is_ref?Boolean

Returns:

  • (Boolean)


61
# File 'ext/nmatrix/ruby_nmatrix.c', line 61

static VALUE nm_is_ref(VALUE self);

#kron_prod(mat) ⇒ Object

Compute the Kronecker product of self and other NMatrix

Arguments

* +mat+ - A 2D NMatrix object

Usage

a = NMatrix.new([2,2],[1,2,
                       3,4])
b = NMatrix.new([2,3],[1,1,1,
                       1,1,1], dtype: :float64)
a.kron_prod(b) # => [ [1.0, 1.0, 1.0, 2.0, 2.0, 2.0]
                      [1.0, 1.0, 1.0, 2.0, 2.0, 2.0]
                      [3.0, 3.0, 3.0, 4.0, 4.0, 4.0]
                      [3.0, 3.0, 3.0, 4.0, 4.0, 4.0] ]


643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
# File 'lib/nmatrix/math.rb', line 643

def kron_prod(mat)
  unless self.dimensions==2 and mat.dimensions==2
    raise ShapeError, "Implemented for 2D NMatrix objects only."
  end

  # compute the shape [n,m] of the product matrix
  n, m = self.shape[0]*mat.shape[0], self.shape[1]*mat.shape[1]
  # compute the entries of the product matrix
  kron_prod_array = []
  if self.yale?
    # +:yale+ requires to get the row by copy in order to apply +#transpose+ to it
    self.each_row(getby=:copy) do |selfr|
      mat.each_row do |matr|
        kron_prod_array += (selfr.transpose.dot matr).to_flat_a
      end
    end
  else
    self.each_row do |selfr|
      mat.each_row do |matr|
        kron_prod_array += (selfr.transpose.dot matr).to_flat_a
      end
    end
  end

  NMatrix.new([n,m], kron_prod_array) 
end

#laswp(ary, opts = {}) ⇒ Object Also known as: permute_columns

call-seq:

laswp(ary) -> NMatrix

Permute the columns of a dense matrix using LASWP according to the order given in an array ary.

If :convention is :lapack, then ary represents a sequence of pair-wise permutations which are performed successively. That is, the i’th entry of ary is the index of the column to swap the i’th column with, having already applied all earlier swaps. This is the default.

If :convention is :intuitive, then ary represents the order of columns after the permutation. That is, the i’th entry of ary is the index of the column that will be in position i after the reordering (Matlab-like behaviour).

Not yet implemented for yale or list.

Arguments

  • ary - An Array specifying the order of the columns. See above for details.

Options

  • :covention - Possible values are :lapack and :intuitive. Default is :lapack. See above for details.



481
482
483
# File 'lib/nmatrix/math.rb', line 481

def laswp(ary, opts={})
  self.clone.laswp!(ary, opts)
end

#laswp!(ary, opts = {}) ⇒ Object Also known as: permute_columns!

call-seq:

laswp!(ary) -> NMatrix

In-place permute the columns of a dense matrix using LASWP according to the order given as an array ary.

If :convention is :lapack, then ary represents a sequence of pair-wise permutations which are performed successively. That is, the i’th entry of ary is the index of the column to swap the i’th column with, having already applied all earlier swaps.

If :convention is :intuitive, then ary represents the order of columns after the permutation. That is, the i’th entry of ary is the index of the column that will be in position i after the reordering (Matlab-like behaviour). This is the default.

Not yet implemented for yale or list.

Arguments

  • ary - An Array specifying the order of the columns. See above for details.

Options

  • :covention - Possible values are :lapack and :intuitive. Default is :intuitive. See above for details.

Raises:



434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
# File 'lib/nmatrix/math.rb', line 434

def laswp!(ary, opts={})
  raise(StorageTypeError, "ATLAS functions only work on dense matrices") unless self.dense?
  opts = { convention: :intuitive }.merge(opts)
  
  if opts[:convention] == :intuitive
    if ary.length != ary.uniq.length
      raise(ArgumentError, "No duplicated entries in the order array are allowed under convention :intuitive")
    end
    n = self.shape[1]
    p = []
    order = (0...n).to_a
    0.upto(n-2) do |i|
      p[i] = order.index(ary[i])
      order[i], order[p[i]] = order[p[i]], order[i]
    end
    p[n-1] = n-1
  else
    p = ary
  end

  NMatrix::LAPACK::laswp(self, p)
end

#layer(layer_number, get_by = :copy) ⇒ Object

call-seq:

layer(layer_number) -> NMatrix
row(layer_number, get_by) -> NMatrix
  • Arguments :

    • layer_number -> Integer.

    • get_by -> Type of slicing to use, :copy or :reference.

  • Returns :

    • A NMatrix representing the requested layer as a layer vector.



855
856
857
# File 'lib/nmatrix/nmatrix.rb', line 855

def layer(layer_number, get_by = :copy)
  rank(2, layer_number, get_by)
end

#ldexpObject

#list?Boolean

call-seq:

m.list? -> true or false

Determine if m is a list-of-lists matrix.

Returns:

  • (Boolean)


53
# File 'lib/nmatrix/shortcuts.rb', line 53

def list?;  return stype == :list; end

#logObject

#log10Object

#log2Object

#lower_triangle(k = 0) ⇒ Object Also known as: tril

call-seq:

lower_triangle -> NMatrix
lower_triangle(k) -> NMatrix
tril -> NMatrix
tril(k) -> NMatrix

Returns the lower triangular portion of a matrix. This is analogous to the tril method in MATLAB.

  • Arguments :

    • k -> Integer. How many extra diagonals to include in the lower triangular portion.

Raises:

  • (NotImplementedError)


802
803
804
805
806
807
808
809
810
811
812
813
814
815
# File 'lib/nmatrix/nmatrix.rb', line 802

def lower_triangle(k = 0)
  raise(NotImplementedError, "only implemented for 2D matrices") if self.shape.size > 2

  t = self.clone_structure
  (0...self.shape[0]).each do |i|
    if i + k >= shape[0]
      t[i, :*] = self[i, :*]
    else
      t[i, (i+k+1)...self.shape[1]] = 0
      t[i, 0..(i+k)] = self[i, 0..(i+k)]
    end
  end
  t
end

#lower_triangle!(k = 0) ⇒ Object Also known as: tril!

call-seq:

lower_triangle! -> NMatrix
lower_triangle!(k) -> NMatrix
tril! -> NMatrix
tril!(k) -> NMatrix

Deletes the upper triangular portion of the matrix (in-place) so only the lower portion remains.

  • Arguments :

    • k -> Integer. How many extra diagonals to include in the deletion.

Raises:

  • (NotImplementedError)


831
832
833
834
835
836
837
838
839
840
# File 'lib/nmatrix/nmatrix.rb', line 831

def lower_triangle!(k = 0)
  raise(NotImplementedError, "only implemented for 2D matrices") if self.shape.size > 2

  (0...self.shape[0]).each do |i|
    if i + k < shape[0]
      self[i, (i+k+1)...self.shape[1]] = 0
    end
  end
  self
end

#map(&bl) ⇒ Object

call-seq:

map -> Enumerator
map { |elem| block } -> NMatrix

Returns an NMatrix if a block is given. For an Array, use #flat_map

Note that #map will always return an :object matrix, because it has no way of knowing how to handle operations on the different dtypes.



72
73
74
75
76
77
# File 'lib/nmatrix/enumerate.rb', line 72

def map(&bl)
  return enum_for(:map) unless block_given?
  cp = self.cast(dtype: :object)
  cp.map!(&bl)
  cp
end

#map!Object

call-seq:

map! -> Enumerator
map! { |elem| block } -> NMatrix

Maps in place.

See Also:



87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
# File 'lib/nmatrix/enumerate.rb', line 87

def map!
  return enum_for(:map!) unless block_given?
  iterated = false
  self.each_stored_with_indices do |e, *i|
    iterated = true
    self[*i] = (yield e)
  end
  #HACK: if there's a single element in a non-dense matrix, it won't iterate and
  #won't change the default value; this ensures that it does get changed.
  unless iterated then
    self.each_with_indices do |e, *i|
      self[*i] = (yield e)
    end
  end
end

#map_storedObject



54
# File 'ext/nmatrix/ruby_nmatrix.c', line 54

static VALUE nm_map_stored(VALUE nmatrix);

#max(dimen = 0) ⇒ Object

call-seq:

max() -> NMatrix
max(dimen) -> NMatrix

Calculates the maximum along the specified dimension.

See Also:



768
769
770
771
772
773
774
775
776
# File 'lib/nmatrix/math.rb', line 768

def max(dimen=0)
  inject_rank(dimen) do |max, sub_mat|
    if max.is_a? NMatrix then
      max * (max >= sub_mat).cast(self.stype, self.dtype) + ((max)*0.0 + (max < sub_mat).cast(self.stype, self.dtype)) * sub_mat
    else
      max >= sub_mat ? max : sub_mat
    end
  end
end

#mean(dimen = 0) ⇒ Object

call-seq:

mean() -> NMatrix
mean(dimen) -> NMatrix

Calculates the mean along the specified dimension.

This will force integer types to float64 dtype.

See Also:



715
716
717
718
719
720
721
722
723
# File 'lib/nmatrix/math.rb', line 715

def mean(dimen=0)
  reduce_dtype = nil
  if integer_dtype? then
    reduce_dtype = :float64
  end
  inject_rank(dimen, 0.0, reduce_dtype) do |mean, sub_mat|
    mean + sub_mat
  end / shape[dimen]
end

#min(dimen = 0) ⇒ Object

call-seq:

min() -> NMatrix
min(dimen) -> NMatrix

Calculates the minimum along the specified dimension.

See Also:



749
750
751
752
753
754
755
756
757
# File 'lib/nmatrix/math.rb', line 749

def min(dimen=0)
  inject_rank(dimen) do |min, sub_mat|
    if min.is_a? NMatrix then
      min * (min <= sub_mat).cast(self.stype, self.dtype) + ((min)*0.0 + (min > sub_mat).cast(self.stype, self.dtype)) * sub_mat
    else
      min <= sub_mat ? min : sub_mat
    end
  end
end

#nrm2(incx = 1, n = nil) ⇒ Object Also known as: norm2

call-seq:

norm2 -> Numeric

Arguments

- +incx+ -> the skip size (defaults to 1, no skip)
- +n+ -> the number of elements to include

Return the 2-norm of the vector. This is the BLAS nrm2 routine.



880
881
882
883
# File 'lib/nmatrix/math.rb', line 880

def nrm2 incx=1, n=nil
  return method_missing(:nrm2, incx, n) unless vector?
  NMatrix::BLAS::nrm2(self, incx, self.size / incx)
end

#nvector?Boolean

call-seq:

nvector? -> true or false

Shortcut function for determining whether the effective dimension is less than the dimension. Useful when we take slices of n-dimensional matrices where n > 2.

Returns:

  • (Boolean)


429
430
431
# File 'lib/nmatrix/nmatrix.rb', line 429

def nvector?
  self.effective_dim < self.dim
end

#object_dtype?Boolean

call-seq:

object_dtype?() -> Boolean

Checks if dtype is a ruby object

Returns:

  • (Boolean)


374
375
376
# File 'lib/nmatrix/nmatrix.rb', line 374

def object_dtype?
  dtype == :object
end

#offsetObject



47
# File 'ext/nmatrix/ruby_nmatrix.c', line 47

static VALUE nm_offset(VALUE self);

#potrf!(which) ⇒ Object

call-seq:

potrf!(upper_or_lower) -> NMatrix

Cholesky factorization of a symmetric positive-definite matrix – or, if complex, a Hermitian positive-definite matrix A. The result will be written in either the upper or lower triangular portion of the matrix, depending on whether the argument is :upper or :lower. Also the function only reads in the upper or lower part of the matrix, so it doesn’t actually have to be symmetric/Hermitian. However, if the matrix (i.e. the symmetric matrix implied by the lower/upper half) is not positive-definite, the function will return nonsense.

This functions requires either the nmatrix-atlas or nmatrix-lapacke gem installed.

  • Returns : the triangular portion specified by the parameter

  • Raises :

    • StorageTypeError -> ATLAS functions only work on dense matrices.

    • ShapeError -> Must be square.

    • NotImplementedError -> If called without nmatrix-atlas or nmatrix-lapacke gem

Raises:

  • (NotImplementedError)


179
180
181
182
# File 'lib/nmatrix/math.rb', line 179

def potrf!(which)
  # The real implementation is in the plugin files.
  raise(NotImplementedError, "potrf! requires either the nmatrix-atlas or nmatrix-lapacke gem")
end

#potrf_lower!Object



188
189
190
# File 'lib/nmatrix/math.rb', line 188

def potrf_lower!
  potrf! :lower
end

#potrf_upper!Object



184
185
186
# File 'lib/nmatrix/math.rb', line 184

def potrf_upper!
  potrf! :upper
end

#pow(n) ⇒ Object

Raise a square matrix to a power. Be careful of numeric overflows! In case n is 0, an identity matrix of the same dimension is returned. In case of negative n, the matrix is inverted and the absolute value of n taken for computing the power.

Arguments

  • n - Integer to which self is to be raised.

References

  • R.G Dromey - How to Solve it by Computer. Link -

    http://www.amazon.com/Solve-Computer-Prentice-Hall-International-Science/dp/0134340019/ref=sr_1_1?ie=UTF8&qid=1422605572&sr=8-1&keywords=how+to+solve+it+by+computer
    

Raises:



598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
# File 'lib/nmatrix/math.rb', line 598

def pow n
  raise ShapeError, "Only works with 2D square matrices." if 
    shape[0] != shape[1] or shape.size != 2
  raise TypeError, "Only works with integer powers" unless n.is_a?(Integer)

  sequence = (integer_dtype? ? self.cast(dtype: :int64) : self).clone
  product  = NMatrix.eye shape[0], dtype: sequence.dtype, stype: sequence.stype 

  if n == 0
    return NMatrix.eye(shape, dtype: dtype, stype: stype)
  elsif n == 1
    return sequence
  elsif n < 0
    n = n.abs
    sequence.invert!
    product = NMatrix.eye shape[0], dtype: sequence.dtype, stype: sequence.stype
  end

  # Decompose n to reduce the number of multiplications.
  while n > 0
    product = product.dot(sequence) if n % 2 == 1
    n = n / 2
    sequence = sequence.dot(sequence)
  end

  product
end

#pretty_print(q) ⇒ Object

TODO: Make this actually pretty.



164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
# File 'lib/nmatrix/nmatrix.rb', line 164

def pretty_print(q) #:nodoc:
  if self.shape.size > 1 and self.shape[1] > 100
    self.inspect.pretty_print(q)
  elsif self.dim > 3 || self.dim == 1
    self.to_a.pretty_print(q)
  else
    # iterate through the whole matrix and find the longest number
    longest = Array.new(self.shape[1], 0)
    self.each_column.with_index do |col, j|
      col.each do |elem|
        elem_len   = elem.inspect.size
        longest[j] = elem_len if longest[j] < elem_len
      end
    end

    if self.dim == 3
      q.group(0, "\n{ layers:", "}") do
        self.each_layer.with_index do |layer,k|
          q.group(0, "\n  [\n", "  ]\n") do
            layer.each_row.with_index do |row,i|
              q.group(0, "    [", "]\n") do
                q.seplist(self[i,0...self.shape[1],k].to_flat_array, lambda { q.text ", "}, :each_with_index) { |v,j| q.text v.inspect.rjust(longest[j]) }
              end
            end
          end
        end
      end
    else # dim 2
      q.group(0, "\n[\n", "]") do
        self.each_row.with_index do |row,i|
          q.group(1, "  [", "]") do
            q.seplist(self.dim > 2 ? row.to_a[0] : row.to_a, lambda { q.text ", " }, :each_with_index) { |v,j| q.text v.inspect.rjust(longest[j]) }
          end
          q.breakable
        end
      end
    end
  end
end

#quaternionObject

call-seq:

quaternion -> NMatrix

Find the quaternion for a 3D rotation matrix.

Code borrowed from: courses.cms.caltech.edu/cs171/quatut.pdf

  • Returns :

    • A length-4 NMatrix representing the corresponding quaternion.

Examples:

n.quaternion # => [1, 0, 0, 0]

Raises:



159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
# File 'lib/nmatrix/homogeneous.rb', line 159

def quaternion
  raise(ShapeError, "Expected square matrix") if self.shape[0] != self.shape[1]
  raise(ShapeError, "Expected 3x3 rotation (or 4x4 homogeneous) matrix") if self.shape[0] > 4 || self.shape[0] < 3

  q = NMatrix.new([4], dtype: self.dtype == :float32 ? :float32: :float64)
  rotation_trace = self[0,0] + self[1,1] + self[2,2]
  if rotation_trace >= 0
    self_w = self.shape[0] == 4 ? self[3,3] : 1.0
    root_of_homogeneous_trace = Math.sqrt(rotation_trace + self_w)
    q[0] = root_of_homogeneous_trace * 0.5
    s = 0.5 / root_of_homogeneous_trace
    q[1] = (self[2,1] - self[1,2]) * s
    q[2] = (self[0,2] - self[2,0]) * s
    q[3] = (self[1,0] - self[0,1]) * s
  else
    h = 0
    h = 1 if self[1,1] > self[0,0]
    h = 2 if self[2,2] > self[h,h]

    case_macro = Proc.new do |i,j,k,ii,jj,kk|
      qq = NMatrix.new([4], dtype: :float64)
      self_w = self.shape[0] == 4 ? self[3,3] : 1.0
      s = Math.sqrt( (self[ii,ii] - (self[jj,jj] + self[kk,kk])) + self_w)
      qq[i] = s*0.5
      s = 0.5 / s
      qq[j] = (self[ii,jj] + self[jj,ii]) * s
      qq[k] = (self[kk,ii] + self[ii,kk]) * s
      qq[0] = (self[kk,jj] - self[jj,kk]) * s
      qq
    end

    case h
    when 0
      q = case_macro.call(1,2,3, 0,1,2)
    when 1
      q = case_macro.call(2,3,1, 1,2,0)
    when 2
      q = case_macro.call(3,1,2, 2,0,1)
    end

    self_w = self.shape[0] == 4 ? self[3,3] : 1.0
    if self_w != 1
      s = 1.0 / Math.sqrt(self_w)
      q[0] *= s
      q[1] *= s
      q[2] *= s
      q[3] *= s
    end
  end

  q
end

#rank(shape_idx, rank_idx, meth = :copy) ⇒ Object

call-seq:

rank(dimension, row_or_column_number) -> NMatrix
rank(dimension, row_or_column_number, :reference) -> NMatrix reference slice

Returns the rank (e.g., row, column, or layer) specified, using slicing by copy as default.

See @row (dimension = 0), @column (dimension = 1)



481
482
483
484
485
486
487
488
489
490
491
492
493
# File 'lib/nmatrix/nmatrix.rb', line 481

def rank(shape_idx, rank_idx, meth = :copy)

  if shape_idx > (self.dim-1)
    raise(RangeError, "#rank call was out of bounds")
  end

  params = Array.new(self.dim)
  params.each.with_index do |v,d|
    params[d] = d == shape_idx ? rank_idx : 0...self.shape[d]
  end

  meth == :reference ? self[*params] : self.slice(*params)
end

#repeat(count, axis) ⇒ Object

call-seq:

repeat(count, axis) -> NMatrix
  • Arguments :

    • count -> how many times NMatrix should be repeated

    • axis -> index of axis along which NMatrix should be repeated

  • Returns :

    • NMatrix created by repeating the existing one along an axis

  • Examples :

    m = NMatrix.new([2, 2], [1, 2, 3, 4])
    m.repeat(2, 0).to_a #<= [[1, 2], [3, 4], [1, 2], [3, 4]]
    m.repeat(2, 1).to_a #<= [[1, 2, 1, 2], [3, 4, 3, 4]]
    

Raises:

  • (ArgumentError)


1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
# File 'lib/nmatrix/nmatrix.rb', line 1010

def repeat(count, axis)
  raise(ArgumentError, 'Matrix should be repeated at least 2 times.') if count < 2
  new_shape = shape
  new_shape[axis] *= count
  new_matrix = NMatrix.new(new_shape)
  slice = new_shape.map { |axis_size| 0...axis_size }
  start = 0
  count.times do
    slice[axis] = start...(start += shape[axis])
    new_matrix[*slice] = self
  end
  new_matrix
end

#reshape(new_shape, *shapes) ⇒ Object

call-seq:

reshape(new_shape) -> NMatrix

Clone a matrix, changing the shape in the process. Note that this function does not do a resize; the product of the new and old shapes’ components must be equal.

  • Arguments :

    • new_shape -> Array of positive Fixnums.

  • Returns :

    • A copy with a different shape.



550
551
552
553
554
555
556
557
558
559
560
561
# File 'lib/nmatrix/nmatrix.rb', line 550

def reshape new_shape,*shapes
  if new_shape.is_a?Fixnum
    newer_shape =  [new_shape]+shapes
  else  # new_shape is an Array
    newer_shape = new_shape
  end
  t = reshape_clone_structure(newer_shape)
  left_params  = [:*]*newer_shape.size
  right_params = [:*]*self.shape.size
  t[*left_params] = self[*right_params]
  t
end

#reshape!(new_shape, *shapes) ⇒ Object

call-seq:

reshape!(new_shape) -> NMatrix
reshape! new_shape  -> NMatrix

Reshapes the matrix (in-place) to the desired shape. Note that this function does not do a resize; the product of the new and old shapes’ components must be equal.

  • Arguments :

    • new_shape -> Array of positive Fixnums.



575
576
577
578
579
580
581
582
583
584
585
586
# File 'lib/nmatrix/nmatrix.rb', line 575

def reshape! new_shape,*shapes
  if self.is_ref?
    raise(ArgumentError, "This operation cannot be performed on reference slices")
  else
    if new_shape.is_a?Fixnum
      shape =  [new_shape]+shapes
    else  # new_shape is an Array
      shape = new_shape
    end
    self.reshape_bang(shape)
  end
end

#respond_to?(method, include_all = false) ⇒ Boolean

:nodoc:

Returns:

  • (Boolean)


943
944
945
946
947
948
949
950
951
# File 'lib/nmatrix/nmatrix.rb', line 943

def respond_to?(method, include_all = false) #:nodoc:
  if [:shuffle, :shuffle!, :each_with_index, :sorted_indices, :binned_sorted_indices, :nrm2, :asum].include?(method.intern) # vector-only methods
    return vector?
  elsif [:each_layer, :layer].include?(method.intern) # 3-or-more dimensions only
    return dim > 2
  else
    super
  end
end

#roundObject



143
# File 'ext/nmatrix/ruby_nmatrix.c', line 143

static VALUE nm_unary_round(int argc, VALUE* argv, VALUE self);

#row(row_number, get_by = :copy) ⇒ Object

call-seq:

row(row_number) -> NMatrix
row(row_number, get_by) -> NMatrix
  • Arguments :

    • row_number -> Integer.

    • get_by -> Type of slicing to use, :copy or :reference.

  • Returns :

    • An NMatrix representing the requested row as a row vector.



533
534
535
# File 'lib/nmatrix/nmatrix.rb', line 533

def row(row_number, get_by = :copy)
  rank(0, row_number, get_by)
end

#rowsObject

call-seq:

rows -> Integer

This shortcut use #shape to return the number of rows (the first dimension) of the matrix.



256
257
258
# File 'lib/nmatrix/nmatrix.rb', line 256

def rows
  shape[0]
end

#shapeObject

handles list and dense, which are n-dimensional



48
# File 'ext/nmatrix/ruby_nmatrix.c', line 48

static VALUE nm_shape(VALUE self);

#shuffle(*args) ⇒ Object

call-seq:

shuffle -> ...
shuffle(rng) -> ...

Re-arranges the contents of an NVector.

TODO: Write more efficient version for Yale, list. TODO: Generalize for more dimensions.



888
889
890
891
892
# File 'lib/nmatrix/nmatrix.rb', line 888

def shuffle(*args)
  method_missing(:shuffle!, *args) if self.effective_dim > 1
  t = self.clone
  t.shuffle!(*args)
end

#shuffle!(*args) ⇒ Object

call-seq:

shuffle! -> ...
shuffle!(random: rng) -> ...

Re-arranges the contents of an NVector.

TODO: Write more efficient version for Yale, list. TODO: Generalize for more dimensions.



870
871
872
873
874
875
876
# File 'lib/nmatrix/nmatrix.rb', line 870

def shuffle!(*args)
  method_missing(:shuffle!, *args) if self.effective_dim > 1
  ary = self.to_flat_a
  ary.shuffle!(*args)
  ary.each.with_index { |v,idx| self[idx] = v }
  self
end

#sinObject

#sinhObject

#sizeObject

call-seq:

size -> Fixnum

Returns the total size of the NMatrix based on its shape.



413
414
415
# File 'lib/nmatrix/nmatrix.rb', line 413

def size
  NMatrix.size(self.shape)
end

#sliceObject



59
# File 'ext/nmatrix/ruby_nmatrix.c', line 59

static VALUE nm_mget(int argc, VALUE* argv, VALUE self);

#solve(b, opts = {}) ⇒ Object

Solve the matrix equation AX = B, where A is self, B is the first argument, and X is returned. A must be a nxn square matrix, while B must be nxm. Only works with dense matrices and non-integer, non-object data types.

Arguments

  • b - the right hand side

Options

  • form - Signifies the form of the matrix A in the linear system AX=B. If not set then it defaults to :general, which uses an LU solver. Other possible values are :lower_tri, :upper_tri and :pos_def (alternatively, non-abbreviated symbols :lower_triangular, :upper_triangular, and :positive_definite can be used. If :lower_tri or :upper_tri is set, then a specialized linear solver for linear systems AX=B with a lower or upper triangular matrix A is used. If :pos_def is chosen, then the linear system is solved via the Cholesky factorization. Note that when :lower_tri or :upper_tri is used, then the algorithm just assumes that all entries in the lower/upper triangle of the matrix are zeros without checking (which can be useful in certain applications).

Usage

a = NMatrix.new [2,2], [3,1,1,2], dtype: dtype
b = NMatrix.new [2,1], [9,8], dtype: dtype
a.solve(b)

# solve an upper triangular linear system more efficiently:
require 'benchmark'
require 'nmatrix/lapacke'
rand_mat = NMatrix.random([10000, 10000], dtype: :float64)
a = rand_mat.triu
b = NMatrix.random([10000, 10], dtype: :float64)
Benchmark.bm(10) do |bm|
  bm.report('general') { a.solve(b) }
  bm.report('upper_tri') { a.solve(b, form: :upper_tri) }
end
#                   user     system      total        real
#  general     73.170000   0.670000  73.840000 ( 73.810086)
#  upper_tri    0.180000   0.000000   0.180000 (  0.182491)

Raises:



302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
# File 'lib/nmatrix/math.rb', line 302

def solve(b, opts = {})
  raise(ShapeError, "Must be called on square matrix") unless self.dim == 2 && self.shape[0] == self.shape[1]
  raise(ShapeError, "number of rows of b must equal number of cols of self") if 
    self.shape[1] != b.shape[0]
  raise(ArgumentError, "only works with dense matrices") if self.stype != :dense
  raise(ArgumentError, "only works for non-integer, non-object dtypes") if 
    integer_dtype? or object_dtype? or b.integer_dtype? or b.object_dtype?

  opts = { form: :general }.merge(opts)
  x    = b.clone
  n    = self.shape[0]
  nrhs = b.shape[1]

  case opts[:form] 
  when :general
    clone = self.clone
    ipiv = NMatrix::LAPACK.clapack_getrf(:row, n, n, clone, n)
    # When we call clapack_getrs with :row, actually only the first matrix
    # (i.e. clone) is interpreted as row-major, while the other matrix (x)
    # is interpreted as column-major. See here: http://math-atlas.sourceforge.net/faq.html#RowSolve
    # So we must transpose x before and after
    # calling it.
    x = x.transpose
    NMatrix::LAPACK.clapack_getrs(:row, :no_transpose, n, nrhs, clone, n, ipiv, x, n)
    x.transpose
  when :upper_tri, :upper_triangular
    raise(ArgumentError, "upper triangular solver does not work with complex dtypes") if
      complex_dtype? or b.complex_dtype?
    # this is the correct function call; see https://github.com/SciRuby/nmatrix/issues/374
    NMatrix::BLAS::cblas_trsm(:row, :left, :upper, false, :nounit, n, nrhs, 1.0, self, n, x, nrhs)
    x
  when :lower_tri, :lower_triangular
    raise(ArgumentError, "lower triangular solver does not work with complex dtypes") if
      complex_dtype? or b.complex_dtype?
    # this is a workaround; see https://github.com/SciRuby/nmatrix/issues/422
    x = x.transpose
    NMatrix::BLAS::cblas_trsm(:row, :right, :lower, :transpose, :nounit, nrhs, n, 1.0, self, n, x, n)
    x.transpose
  when :pos_def, :positive_definite
    u, l = self.factorize_cholesky
    z = l.solve(b, form: :lower_tri)
    u.solve(z, form: :upper_tri)
  else
    raise(ArgumentError, "#{opts[:form]} is not a valid form option")
  end
end

#sorted_indicesObject

call-seq:

sorted_indices -> Array

Returns an array of the indices ordered by value sorted.



901
902
903
904
905
# File 'lib/nmatrix/nmatrix.rb', line 901

def sorted_indices
  return method_missing(:sorted_indices) unless vector?
  ary = self.to_flat_array
  ary.each_index.sort_by { |i| ary[i] }  # from: http://stackoverflow.com/a/17841159/170300
end

#sqrtObject

#std(dimen = 0) ⇒ Object

call-seq:

std() -> NMatrix
std(dimen) -> NMatrix

Calculates the sample standard deviation along the specified dimension.

This will force integer types to float64 dtype.

See Also:



813
814
815
# File 'lib/nmatrix/math.rb', line 813

def std(dimen=0)
  variance(dimen).sqrt
end

#stypeObject



42
# File 'ext/nmatrix/ruby_nmatrix.c', line 42

static VALUE nm_stype(VALUE self);

#sum(dimen = 0) ⇒ Object

call-seq:

sum() -> NMatrix
sum(dimen) -> NMatrix

Calculates the sum along the specified dimension.

See Also:



733
734
735
736
737
# File 'lib/nmatrix/math.rb', line 733

def sum(dimen=0)
  inject_rank(dimen, 0.0) do |sum, sub_mat|
    sum + sub_mat
  end
end

#supershapeObject



49
# File 'ext/nmatrix/ruby_nmatrix.c', line 49

static VALUE nm_supershape(VALUE self);

#symmetric?Boolean

Returns:

  • (Boolean)


149
# File 'ext/nmatrix/ruby_nmatrix.c', line 149

static VALUE nm_symmetric(VALUE self);

#tanObject

#tanhObject

#to_a(dimen = nil) ⇒ Object

call-seq:

to_a -> Array

Converts an NMatrix to an array of arrays, or an NMatrix of effective dimension 1 to an array.

Does not yet work for dimensions > 2



451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
# File 'lib/nmatrix/nmatrix.rb', line 451

def to_a(dimen=nil)
  if self.dim == 2

    return self.to_flat_a if self.shape[0] == 1

    ary = []
    begin
      self.each_row do |row|
        ary << row.to_flat_a
      end
    #rescue NotImplementedError # Oops. Try copying instead
    #  self.each_row(:copy) do |row|
    #    ary << row.to_a.flatten
    #  end
    end
    ary
  else
    to_a_rec(0)
  end
end

#to_fObject

call-seq:

to_f -> Float

Converts an nmatrix with a single element (but any number of dimensions)

to a float.

Raises an IndexError if the matrix does not have just a single element.

Raises:

  • (IndexError)


388
389
390
391
# File 'lib/nmatrix/nmatrix.rb', line 388

def to_f
  raise IndexError, 'to_f only valid for matrices with a single element' unless shape.all? { |e| e == 1 }
  self[*Array.new(shape.size, 0)]
end

#to_flat_arrayObject Also known as: to_flat_a

call-seq:

to_flat_array -> Array
to_flat_a -> Array

Converts an NMatrix to a one-dimensional Ruby Array.



400
401
402
403
404
# File 'lib/nmatrix/nmatrix.rb', line 400

def to_flat_array
  ary = Array.new(self.size)
  self.each.with_index { |v,i| ary[i] = v }
  ary
end

#to_hashObject Also known as: to_h

call-seq:

to_hash -> Hash

Create a Ruby Hash from an NMatrix.



306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
# File 'lib/nmatrix/nmatrix.rb', line 306

def to_hash
  if stype == :yale
    h = {}
    each_stored_with_indices do |val,i,j|
      next if val == 0 # Don't bother storing the diagonal zero values -- only non-zeros.
      if h.has_key?(i)
        h[i][j] = val
      else
        h[i] = {j => val}
      end
    end
    h
  else # dense and list should use a C internal function.
    # FIXME: Write a C internal to_h function.
    m = stype == :dense ? self.cast(:list, self.dtype) : self
    m.__list_to_hash__
  end
end

#to_sObject

:nodoc:



418
419
420
# File 'lib/nmatrix/nmatrix.rb', line 418

def to_s #:nodoc:
  self.to_flat_array.to_s
end

#traceObject

call-seq:

trace -> Numeric

Calculates the trace of an nxn matrix.

  • Raises :

    • ShapeError -> Expected square matrix

  • Returns :

    • The trace of the matrix (a numeric value)

Raises:



696
697
698
699
700
701
702
# File 'lib/nmatrix/math.rb', line 696

def trace
  raise(ShapeError, "Expected square matrix") unless self.shape[0] == self.shape[1] && self.dim == 2

  (0...self.shape[0]).inject(0) do |total,i|
    total + self[i,i]
  end
end

#transpose(permute = nil) ⇒ Object

call-seq:

transpose -> NMatrix
transpose(permutation) -> NMatrix

Clone a matrix, transposing it in the process. If the matrix is two-dimensional, the permutation is taken to be [1,0] automatically (switch dimension 0 with dimension 1). If the matrix is n-dimensional, you must provide a permutation of 0...n.

  • Arguments :

    • permutation -> Optional Array giving a permutation.

  • Returns :

    • A copy of the matrix, but transposed.



602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
# File 'lib/nmatrix/nmatrix.rb', line 602

def transpose(permute = nil)
  if permute.nil?
    if self.dim == 1
      return self.clone
    elsif self.dim == 2
      new_shape = [self.shape[1], self.shape[0]]
    else
      raise(ArgumentError, "need permutation array of size #{self.dim}")
    end
  elsif !permute.is_a?(Array) || permute.sort.uniq != (0...self.dim).to_a
    raise(ArgumentError, "invalid permutation array")
  else
    # Figure out the new shape based on the permutation given as an argument.
    new_shape = permute.map { |p| self.shape[p] }
  end

  if self.dim > 2 # FIXME: For dense, several of these are basically equivalent to reshape.

    # Make the new data structure.
    t = self.reshape_clone_structure(new_shape)

    self.each_stored_with_indices do |v,*indices|
      p_indices = permute.map { |p| indices[p] }
      t[*p_indices] = v
    end
    t
  elsif self.list? # TODO: Need a C list transposition algorithm.
    # Make the new data structure.
    t = self.reshape_clone_structure(new_shape)

    self.each_column.with_index do |col,j|
      t[j,:*] = col.to_flat_array
    end
    t
  else
    # Call C versions of Yale and List transpose, which do their own copies
    self.clone_transpose
  end
end

#upper_triangle(k = 0) ⇒ Object Also known as: triu

call-seq:

upper_triangle -> NMatrix
upper_triangle(k) -> NMatrix
triu -> NMatrix
triu(k) -> NMatrix

Returns the upper triangular portion of a matrix. This is analogous to the triu method in MATLAB.

  • Arguments :

    • k -> Positive integer. How many extra diagonals to include in the upper triangular portion.

Raises:

  • (NotImplementedError)


747
748
749
750
751
752
753
754
755
756
757
758
759
760
# File 'lib/nmatrix/nmatrix.rb', line 747

def upper_triangle(k = 0)
  raise(NotImplementedError, "only implemented for 2D matrices") if self.shape.size > 2

  t = self.clone_structure
  (0...self.shape[0]).each do |i|
    if i - k < 0
      t[i, :*] = self[i, :*]
    else
      t[i, 0...(i-k)]             = 0
      t[i, (i-k)...self.shape[1]] = self[i, (i-k)...self.shape[1]]
    end
  end
  t
end

#upper_triangle!(k = 0) ⇒ Object Also known as: triu!

call-seq:

upper_triangle! -> NMatrix
upper_triangle!(k) -> NMatrix
triu! -> NMatrix
triu!(k) -> NMatrix

Deletes the lower triangular portion of the matrix (in-place) so only the upper portion remains.

  • Arguments :

    • k -> Integer. How many extra diagonals to include in the deletion.

Raises:

  • (NotImplementedError)


776
777
778
779
780
781
782
783
784
785
# File 'lib/nmatrix/nmatrix.rb', line 776

def upper_triangle!(k = 0)
  raise(NotImplementedError, "only implemented for 2D matrices") if self.shape.size > 2

  (0...self.shape[0]).each do |i|
    if i - k >= 0
      self[i, 0...(i-k)] = 0
    end
  end
  self
end

#variance(dimen = 0) ⇒ Object

call-seq:

variance() -> NMatrix
variance(dimen) -> NMatrix

Calculates the sample variance along the specified dimension.

This will force integer types to float64 dtype.

See Also:



790
791
792
793
794
795
796
797
798
799
# File 'lib/nmatrix/math.rb', line 790

def variance(dimen=0)
  reduce_dtype = nil
  if integer_dtype? then
    reduce_dtype = :float64
  end
  m = mean(dimen)
  inject_rank(dimen, 0.0, reduce_dtype) do |var, sub_mat|
    var + (m - sub_mat)*(m - sub_mat)/(shape[dimen]-1)
  end
end

#vconcat(*matrices) ⇒ Object

Vertical concatenation with matrices.



724
725
726
# File 'lib/nmatrix/nmatrix.rb', line 724

def vconcat(*matrices)
  concat(*matrices, :row)
end

#vector?Boolean

call-seq:

vector? -> true or false

Shortcut function for determining whether the effective dimension is 1. See also #nvector?

Returns:

  • (Boolean)


439
440
441
# File 'lib/nmatrix/nmatrix.rb', line 439

def vector?
  self.effective_dim == 1
end

#writeObject



38
# File 'ext/nmatrix/ruby_nmatrix.c', line 38

static VALUE nm_write(int argc, VALUE* argv, VALUE self);

#yale?Boolean

call-seq:

m.yale? -> true or false

Determine if m is a Yale matrix.

Returns:

  • (Boolean)


47
# File 'lib/nmatrix/shortcuts.rb', line 47

def yale?;  return stype == :yale; end