class N

NMatrix is a matrix class that supports both multidimensional arrays (`:dense` stype) and sparse storage (`:list` or `:yale` stypes) and 13 data types, including complex numbers, various integer and floating-point sizes and ruby objects.

Public Class Methods

NMatrix[Numeric, ..., Numeric, dtype: Symbol] → NMatrix click to toggle source
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 ::[] 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] ]
# 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_diag(*params)
Alias for: block_diagonal
block_diagonal(*params) click to toggle source

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]
      ]
# 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
Also aliased as: block_diag
diag(entries, opts={})
Alias for: diagonal
diagonals(array) → NMatrix click to toggle source
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
# 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
Also aliased as: diag, diagonals
diagonals(entries, opts={})
Alias for: diagonal
eye(shape) → NMatrix click to toggle source
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
# 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
Also aliased as: identity
guess_dtype(p1) click to toggle source

Guess the dtype given a Ruby VALUE and return it as a symbol.

Not to be confused with nm_dtype_guess, which returns an nm::dtype_t. (This calls that.)

static VALUE nm_guess_dtype(VALUE self, VALUE v) {
  return ID2SYM(rb_intern(DTYPE_NAMES[nm_dtype_guess(v)]));
}
identity(shape, opts={})
Alias for: eye
load_matlab_file(path) → Mat5Reader click to toggle source
  • Arguments :

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

  • Returns :

    • A Mat5Reader object.

# 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(path) → PointCloudReader::MetaReader click to toggle source
  • Arguments :

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

  • Returns :

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

# File lib/nmatrix/nmatrix.rb, line 99
def load_pcd_file(file_path)
  NMatrix::IO::PointCloudReader::MetaReader.new(file_path)
end
meshgrid(arrs) → Array of NMatrix click to toggle source
meshgrid(arrs, options) → Array of NMatrix

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.

  • 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 :

  • 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]]
    
# 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_dtype(p1) click to toggle source

Get the minimum allowable dtype for a Ruby VALUE and return it as a symbol.

static VALUE nm_min_dtype(VALUE self, VALUE v) {
  return ID2SYM(rb_intern(DTYPE_NAMES[nm_dtype_min(v)]));
}
new(shape) → NMatrix click to toggle source
new(shape, initial_value) → NMatrix
new(shape, initial_array) → NMatrix
new(shape, initial_value, options) → NMatrix
new(shape, initial_array, options) → NMatrix

Create a new NMatrix.

The only mandatory argument is shape, which may be a positive integer or an array of positive integers.

It is recommended that you supply an initialization value or array of values. Without one, Yale and List matrices will be initialized to 0; and dense matrices will be undefined.

Additional options may be provided using keyword arguments. The keywords are +:dtype, :stype, :capacity, and :default. Only Yale uses a capacity argument, which is used to reserve the initial size of its storage vectors. List and Yale both accept a default value (which itself defaults to 0). This default is taken from the initial value if such a value is given; it is more likely to be required when an initial array is provided.

The storage type, or stype, is used to specify whether we want a :dense, :list, or :yale matrix; dense is the default.

The data type, or dtype, can be one of: :byte, :int8, :int16, :int32, :int64, :float32, :float64, :complex64, :complex128, or :object. The constructor will attempt to guess it from the initial value/array/default provided, if any. Otherwise, the default is :object, which stores any type of Ruby object.

In addition to the above, there is a legacy constructor from the alpha version. To use that version, you must be providing exactly four arguments. It is now deprecated.

There is one additional constructor for advanced users, which takes seven arguments and is only for creating Yale matrices with known IA, JA, and A arrays. This is used primarily internally for IO, e.g., reading Matlab matrices, which are stored in old Yale (not our Yale) format. But be careful; there are no overflow warnings. All of these constructors are defined for power-users. Everyone else should probably resort to the shortcut functions defined in shortcuts.rb.

static VALUE nm_init(int argc, VALUE* argv, VALUE nm) {
  NM_CONSERVATIVE(nm_register_value(&nm));
  NM_CONSERVATIVE(nm_register_values(argv, argc));

  if (argc <= 3) { // Call the new constructor unless all four arguments are given (or the 7-arg version is given)
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&nm));
        return nm_init_new_version(argc, argv, nm);
  }

  /* First, determine stype (dense by default) */
  nm::stype_t stype;
  size_t  offset = 0;

  if (!SYMBOL_P(argv[0]) && TYPE(argv[0]) != T_STRING) {
    stype = nm::DENSE_STORE;

  } else {
    // 0: String or Symbol
    stype  = interpret_stype(argv[0]);
    offset = 1;
  }

  // If there are 7 arguments and Yale, refer to a different init function with fewer sanity checks.
  if (argc == 7) {
    if (stype == nm::YALE_STORE) {
      NM_CONSERVATIVE(nm_unregister_values(argv, argc));
      NM_CONSERVATIVE(nm_unregister_value(&nm));
      return nm_init_yale_from_old_yale(argv[1], argv[2], argv[3], argv[4], argv[5], argv[6], nm);

    } else {
      NM_CONSERVATIVE(nm_unregister_values(argv, argc));
      NM_CONSERVATIVE(nm_unregister_value(&nm));
      rb_raise(rb_eArgError, "Expected 2-4 arguments (or 7 for internal Yale creation)");
    }
  }

  // 1: Array or Fixnum
  size_t dim;
  size_t* shape = interpret_shape(argv[offset], &dim);

  // 2-3: dtype
  nm::dtype_t dtype = interpret_dtype(argc-1-offset, argv+offset+1, stype);

  size_t init_cap = 0, init_val_len = 0;
  void* init_val  = NULL;
  if (!SYMBOL_P(argv[1+offset]) || TYPE(argv[1+offset]) == T_ARRAY) {
        // Initial value provided (could also be initial capacity, if yale).

    if (stype == nm::YALE_STORE && NM_RUBYVAL_IS_NUMERIC(argv[1+offset])) {
      init_cap = FIX2UINT(argv[1+offset]);

    } else {
        // 4: initial value / dtype
      init_val = interpret_initial_value(argv[1+offset], dtype);

      if (TYPE(argv[1+offset]) == T_ARRAY)      init_val_len = RARRAY_LEN(argv[1+offset]);
      else                                  init_val_len = 1;
    }

  } else {
        // DType is RUBYOBJ.

    if (stype == nm::DENSE_STORE) {
        /*
         * No need to initialize dense with any kind of default value unless it's
         * an RUBYOBJ matrix.
         */
      if (dtype == nm::RUBYOBJ) {
        // Pretend [nil] was passed for RUBYOBJ.
        init_val = NM_ALLOC(VALUE);
        *(VALUE*)init_val = Qnil;

        init_val_len = 1;

      } else {
        init_val = NULL;
      }
    } else if (stype == nm::LIST_STORE) {
      init_val = NM_ALLOC_N(char, DTYPE_SIZES[dtype]);
      std::memset(init_val, 0, DTYPE_SIZES[dtype]);
    }
  }

  if (dtype == nm::RUBYOBJ) {
    nm_register_values(reinterpret_cast<VALUE*>(init_val), init_val_len);
  }

  // TODO: Update to allow an array as the initial value.
  NMATRIX* nmatrix;
  UnwrapNMatrix(nm, nmatrix);

  nmatrix->stype = stype;

  switch (stype) {
    case nm::DENSE_STORE:
      nmatrix->storage = (STORAGE*)nm_dense_storage_create(dtype, shape, dim, init_val, init_val_len);
      break;

    case nm::LIST_STORE:
      nmatrix->storage = (STORAGE*)nm_list_storage_create(dtype, shape, dim, init_val);
      break;

    case nm::YALE_STORE:
      nmatrix->storage = (STORAGE*)nm_yale_storage_create(dtype, shape, dim, init_cap);
      nm_yale_storage_init((YALE_STORAGE*)(nmatrix->storage), NULL);
      break;
  }

  if (dtype == nm::RUBYOBJ) {
    nm_unregister_values(reinterpret_cast<VALUE*>(init_val), init_val_len);
  }

  NM_CONSERVATIVE(nm_unregister_values(argv, argc));
  NM_CONSERVATIVE(nm_unregister_value(&nm));

  return nm;
}
ones(shape) → NMatrix click to toggle source
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 :

Examples:

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

NMatrix.ones([2, 3], dtype: :int32) # =>  1  1  1
                                          1  1  1
# 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 click to toggle source

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

@param [NMatrix] nm the nmatrix whose dtype and shape will be used @return [NMatrix] a new nmatrix filled with ones.

# 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
rand(shape, opts={})
Alias for: random
random(shape) → NMatrix click to toggle source

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] ]
# 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
Also aliased as: rand
read(p1, p2 = v2) click to toggle source

Binary file reader for NMatrix standard format. file should be a path, which we aren't going to check very carefully (in other words, this function should generally be called from a Ruby helper method).

Note that currently, this function will by default refuse to read files that are newer than your version of NMatrix. To force an override, set the second argument to anything other than nil.

Returns an NMatrix Ruby object.

static VALUE nm_read(int argc, VALUE* argv, VALUE self) {
  using std::ifstream;

  NM_CONSERVATIVE(nm_register_values(argv, argc));
  NM_CONSERVATIVE(nm_register_value(&self));

  VALUE file, force_;

  // Read the arguments
  rb_scan_args(argc, argv, "11", &file, &force_);
  bool force   = (force_ != Qnil && force_ != Qfalse);


  if (!RB_FILE_EXISTS(file)) { // FIXME: Errno::ENOENT
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&self));
    rb_raise(rb_get_errno_exc("ENOENT"), "%s", RSTRING_PTR(file));
  }

  // Open a file stream
  ifstream f(RSTRING_PTR(file), std::ios::in | std::ios::binary);

  uint16_t major, minor, release;
  get_version_info(major, minor, release); // compare to NMatrix version

  uint16_t fmajor, fminor, frelease, null16;

  // READ FIRST 64-BIT BLOCK
  f.read(reinterpret_cast<char*>(&fmajor),   sizeof(uint16_t));
  f.read(reinterpret_cast<char*>(&fminor),   sizeof(uint16_t));
  f.read(reinterpret_cast<char*>(&frelease), sizeof(uint16_t));
  f.read(reinterpret_cast<char*>(&null16),   sizeof(uint16_t));

  int ver  = major * 10000 + minor * 100 + release,
      fver = fmajor * 10000 + fminor * 100 + release;
  if (fver > ver && force == false) {
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&self));
    rb_raise(rb_eIOError, "File was created in newer version of NMatrix than current (%u.%u.%u)", fmajor, fminor, frelease);
  }
  if (null16 != 0) rb_warn("nm_read: Expected zero padding was not zero (0)\n");

  uint8_t dt, st, it, sm;
  uint16_t dim;

  // READ SECOND 64-BIT BLOCK
  f.read(reinterpret_cast<char*>(&dt), sizeof(uint8_t));
  f.read(reinterpret_cast<char*>(&st), sizeof(uint8_t));
  f.read(reinterpret_cast<char*>(&it), sizeof(uint8_t)); // FIXME: should tell how few bytes indices are stored as
  f.read(reinterpret_cast<char*>(&sm), sizeof(uint8_t));
  f.read(reinterpret_cast<char*>(&null16), sizeof(uint16_t));
  f.read(reinterpret_cast<char*>(&dim), sizeof(uint16_t));

  if (null16 != 0) rb_warn("nm_read: Expected zero padding was not zero (1)");
  nm::stype_t stype = static_cast<nm::stype_t>(st);
  nm::dtype_t dtype = static_cast<nm::dtype_t>(dt);
  nm::symm_t  symm  = static_cast<nm::symm_t>(sm);
  //nm::itype_t itype = static_cast<nm::itype_t>(it);

  // READ NEXT FEW 64-BIT BLOCKS
  size_t* shape = NM_ALLOC_N(size_t, dim);
  read_padded_shape(f, dim, shape);

  STORAGE* s;
  if (stype == nm::DENSE_STORE) {
    s = nm_dense_storage_create(dtype, shape, dim, NULL, 0);
    nm_register_storage(stype, s);

    read_padded_dense_elements(f, reinterpret_cast<DENSE_STORAGE*>(s), symm, dtype);

  } else if (stype == nm::YALE_STORE) {
    uint32_t ndnz, length;

    // READ YALE-SPECIFIC 64-BIT BLOCK
    f.read(reinterpret_cast<char*>(&ndnz),     sizeof(uint32_t));
    f.read(reinterpret_cast<char*>(&length),   sizeof(uint32_t));

    s = nm_yale_storage_create(dtype, shape, dim, length); // set length as init capacity

    nm_register_storage(stype, s);

    read_padded_yale_elements(f, reinterpret_cast<YALE_STORAGE*>(s), length, symm, dtype);
  } else {
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&self));
    rb_raise(nm_eStorageTypeError, "please convert to yale or dense before saving");
  }

  NMATRIX* nm = nm_create(stype, s);

  // Return the appropriate matrix object (Ruby VALUE)
  // FIXME: This should probably return CLASS_OF(self) instead of cNMatrix, but I don't know how that works for
  // FIXME: class methods.
  nm_register_nmatrix(nm);
  VALUE to_return = Data_Wrap_Struct(cNMatrix, nm_mark, nm_delete, nm);

  nm_unregister_nmatrix(nm);
  NM_CONSERVATIVE(nm_unregister_values(argv, argc));
  NM_CONSERVATIVE(nm_unregister_value(&self));
  nm_unregister_storage(stype, s);

  switch(stype) {
  case nm::DENSE_STORE:
  case nm::YALE_STORE:
    return to_return;
  default: // this case never occurs (due to earlier rb_raise)
    return Qnil;
  }

}
register_lapack_extension(name) click to toggle source
# File lib/nmatrix/lapack_ext_common.rb, line 30
def NMatrix.register_lapack_extension(name)
  if (defined? @@lapack_extension)
    raise "Attempting to load #{name} when #{@@lapack_extension} is already loaded. You can only load one LAPACK extension."
  end

  @@lapack_extension = name
end
seq(shape) → NMatrix click to toggle source
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
# 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) click to toggle source

Calculate the size of an NMatrix of a given shape.

# 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(x, y, z) → NMatrix click to toggle source
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
# 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
upcast(first_dtype, second_dtype) → Symbol click to toggle source

Given a binary operation between types t1 and t2, what type will be returned?

This is a singleton method on NMatrix, e.g., ::upcast(:int32, :int64)

static VALUE nm_upcast(VALUE self, VALUE t1, VALUE t2) {
  nm::dtype_t d1    = nm_dtype_from_rbsymbol(t1),
              d2    = nm_dtype_from_rbsymbol(t2);

  return ID2SYM(rb_intern( DTYPE_NAMES[ Upcast[d1][d2] ] ));
}
x_rotation(angle_in_radians) → NMatrix click to toggle source
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
# 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={}) click to toggle source
# 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={}) click to toggle source
# 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
zeroes(shape, opts = {})
Alias for: zeros
zeros(shape) → NMatrix click to toggle source
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 :

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
# File lib/nmatrix/shortcuts.rb, line 146
def zeros(shape, opts = {})
  NMatrix.new(shape, 0, {:dtype => :float64}.merge(opts))
end
Also aliased as: zeroes
zeros_like(nm) → NMatrix click to toggle source

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

@param [NMatrix] nm the nmatrix whose stype, dtype, and shape will be used @return [NMatrix] a new nmatrix filled with zeros.

# 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

Public Instance Methods

==(p1) click to toggle source

Equality operator. Returns a single true or false value indicating whether the matrices are equivalent.

For elementwise, use =~ instead.

This method will raise an exception if dimensions do not match.

When stypes differ, this function calls a protected Ruby method.

static VALUE nm_eqeq(VALUE left, VALUE right) {
  NM_CONSERVATIVE(nm_register_value(&left));
  NM_CONSERVATIVE(nm_register_value(&right));

  NMATRIX *l, *r;

  CheckNMatrixType(left);
  CheckNMatrixType(right);

  UnwrapNMatrix(left, l);
  UnwrapNMatrix(right, r);

  bool result = false;

  // Check that the shapes match before going any further.
  if (l->storage->dim != r->storage->dim) {
    NM_CONSERVATIVE(nm_unregister_value(&left));
    NM_CONSERVATIVE(nm_unregister_value(&right));
    rb_raise(nm_eShapeError, "cannot compare matrices with different dimension");
  }

  size_t dim = l->storage->dim;
  for (size_t i=0; i<dim; i++) {
    if (l->storage->shape[i] != r->storage->shape[i]) {
      NM_CONSERVATIVE(nm_unregister_value(&left));
      NM_CONSERVATIVE(nm_unregister_value(&right));
      rb_raise(nm_eShapeError, "cannot compare matrices with different shapes");
    }
  }

  if (l->stype != r->stype) { // DIFFERENT STYPES

    if (l->stype == nm::DENSE_STORE)
      result = rb_funcall(left, rb_intern("dense_eql_sparse?"), 1, right);
    else if (r->stype == nm::DENSE_STORE)
      result = rb_funcall(right, rb_intern("dense_eql_sparse?"), 1, left);
    else
      result = rb_funcall(left, rb_intern("sparse_eql_sparse?"), 1, right);

  } else {

    switch(l->stype) {       // SAME STYPES
    case nm::DENSE_STORE:
      result = nm_dense_storage_eqeq(l->storage, r->storage);
      break;
    case nm::LIST_STORE:
      result = nm_list_storage_eqeq(l->storage, r->storage);
      break;
    case nm::YALE_STORE:
      result = nm_yale_storage_eqeq(l->storage, r->storage);
      break;
    }
  }

  NM_CONSERVATIVE(nm_unregister_value(&left));
  NM_CONSERVATIVE(nm_unregister_value(&right));

  return result ? Qtrue : Qfalse;
}
matrix[indices] → ... click to toggle source

Access the contents of an NMatrix at given coordinates by reference.

n[3,3]  # => 5.0
n[0..1,0..1] #=> matrix [2,2]
static VALUE nm_mref(int argc, VALUE* argv, VALUE self) {
  static void* (*ttable[nm::NUM_STYPES])(const STORAGE*, SLICE*) = {
    nm_dense_storage_ref,
    nm_list_storage_ref,
    nm_yale_storage_ref
  };
  nm::stype_t stype = NM_STYPE(self);
  return nm_xslice(argc, argv, ttable[stype], nm_delete_ref, self);
}
[]=(*args) click to toggle source

Modify the contents of an NMatrix in the given cell

n[3,3] = 5.0

Also returns the new contents, so you can chain:

n[3,3] = n[2,3] = 5.0
static VALUE nm_mset(int argc, VALUE* argv, VALUE self) {

  size_t dim = NM_DIM(self); // last arg is the value

  VALUE to_return = Qnil;

  if ((size_t)(argc) > NM_DIM(self)+1) {
    rb_raise(rb_eArgError, "wrong number of arguments (%d for %lu)", argc, effective_dim(NM_STORAGE(self))+1);
  } else {
    NM_CONSERVATIVE(nm_register_value(&self));
    NM_CONSERVATIVE(nm_register_values(argv, argc));

    SLICE* slice = get_slice(dim, argc-1, argv, NM_STORAGE(self)->shape);

    static void (*ttable[nm::NUM_STYPES])(VALUE, SLICE*, VALUE) = {
      nm_dense_storage_set,
      nm_list_storage_set,
      nm_yale_storage_set
    };

    ttable[NM_STYPE(self)](self, slice, argv[argc-1]);

    free_slice(slice);

    to_return = argv[argc-1];

    NM_CONSERVATIVE(nm_unregister_value(&self));
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
  }

  return to_return;
}
abs → NMatrix click to toggle source

Maps all values in a matrix to their absolute values.

# File lib/nmatrix/math.rb, line 784
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_dtype → Symbol click to toggle source

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).

# File lib/nmatrix/math.rb, line 768
def abs_dtype
  if self.dtype == :complex64
    :float32
  elsif self.dtype == :complex128
    :float64
  else
    self.dtype
  end
end
absolute_sum(incx=1, n=nil)
Alias for: asum
angle_vector → [angle, about_vector] click to toggle source

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]
# 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
absolute_sum → Numeric click to toggle source

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.

# File lib/nmatrix/math.rb, line 805
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
Also aliased as: absolute_sum
binned_sorted_indices → Array click to toggle source

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.

# File lib/nmatrix/nmatrix.rb, line 913
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
capacity() click to toggle source

Find the capacity of an NMatrix. The capacity only differs from the size for Yale matrices, which occasionally allocate more space than they need. For list and dense, capacity gives the number of elements in the matrix.

If you call this on a slice, it may behave unpredictably. Most likely it'll just return the original matrix's capacity.

static VALUE nm_capacity(VALUE self) {
  NM_CONSERVATIVE(nm_register_value(&self));
  VALUE cap;

  switch(NM_STYPE(self)) {
  case nm::YALE_STORE:
    cap = UINT2NUM(reinterpret_cast<YALE_STORAGE*>(NM_STORAGE_YALE(self)->src)->capacity);
    break;

  case nm::DENSE_STORE:
    cap = UINT2NUM(nm_storage_count_max_elements( NM_STORAGE_DENSE(self) ));
    break;

  case nm::LIST_STORE:
    cap = UINT2NUM(nm_list_storage_count_elements( NM_STORAGE_LIST(self) ));
    break;

  default:
    NM_CONSERVATIVE(nm_unregister_value(&self));
    rb_raise(nm_eStorageTypeError, "unrecognized stype in nm_capacity()");
  }

  NM_CONSERVATIVE(nm_unregister_value(&self));
  return cap;
}
cast(stype, dtype, default) → NMatrix click to toggle source
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.

# 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_full(stype) → NMatrix click to toggle source
cast_full(stype, dtype, sparse_basis) → NMatrix

Copy constructor for changing dtypes and stypes.

VALUE nm_cast(VALUE self, VALUE new_stype_symbol, VALUE new_dtype_symbol, VALUE init) {
  NM_CONSERVATIVE(nm_register_value(&self));
  NM_CONSERVATIVE(nm_register_value(&init));

  nm::dtype_t new_dtype = nm_dtype_from_rbsymbol(new_dtype_symbol);
  nm::stype_t new_stype = nm_stype_from_rbsymbol(new_stype_symbol);

  CheckNMatrixType(self);
  NMATRIX *rhs;

  UnwrapNMatrix( self, rhs );

  void* init_ptr = NM_ALLOCA_N(char, DTYPE_SIZES[new_dtype]);
  rubyval_to_cval(init, new_dtype, init_ptr);

  NMATRIX* m = nm_cast_with_ctype_args(rhs, new_stype, new_dtype, init_ptr);
  nm_register_nmatrix(m);

  VALUE to_return = Data_Wrap_Struct(CLASS_OF(self), nm_mark, nm_delete, m);

  nm_unregister_nmatrix(m);
  NM_CONSERVATIVE(nm_unregister_value(&self));
  NM_CONSERVATIVE(nm_unregister_value(&init));
  return to_return;

}
clone_structure → NMatrix click to toggle source

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.

# File lib/nmatrix/nmatrix.rb, line 989
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
col(column_number, get_by = :copy)
Alias for: column
cols → Integer click to toggle source

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

# File lib/nmatrix/nmatrix.rb, line 267
def cols
  shape[1]
end
column(column_number) → NMatrix click to toggle source
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
# File lib/nmatrix/nmatrix.rb, line 516
def column(column_number, get_by = :copy)
  rank(1, column_number, get_by)
end
Also aliased as: col
complex_conjugate → NMatrix click to toggle source

Transform the matrix (non-in-place) to its complex conjugate. Only works on complex matrices.

static VALUE nm_complex_conjugate(VALUE self) {
  VALUE copy;
  return nm_complex_conjugate_bang(nm_init_copy(copy,self));
}
complex_conjugate_bang → NMatrix click to toggle source

Transform the matrix (in-place) to its complex conjugate. Only works on complex matrices.

Bang should imply that no copy is being made, even temporarily.

static VALUE nm_complex_conjugate_bang(VALUE self) {

  NMATRIX* m;
  void* elem;
  size_t size, p;

  UnwrapNMatrix(self, m);

  if (m->stype == nm::DENSE_STORE) {

    size = nm_storage_count_max_elements(NM_STORAGE(self));
    elem = NM_STORAGE_DENSE(self)->elements;

  } else if (m->stype == nm::YALE_STORE) {

    size = nm_yale_storage_get_size(NM_STORAGE_YALE(self));
    elem = NM_STORAGE_YALE(self)->a;

  } else {
    rb_raise(rb_eNotImpError, "please cast to yale or dense (complex) first");
  }

  // Walk through and negate the imaginary component
  if (NM_DTYPE(self) == nm::COMPLEX64) {

    for (p = 0; p < size; ++p) {
      reinterpret_cast<nm::Complex64*>(elem)[p].i = -reinterpret_cast<nm::Complex64*>(elem)[p].i;
    }

  } else if (NM_DTYPE(self) == nm::COMPLEX128) {

    for (p = 0; p < size; ++p) {
      reinterpret_cast<nm::Complex128*>(elem)[p].i = -reinterpret_cast<nm::Complex128*>(elem)[p].i;
    }

  }
  return self;
}
complex_dtype?() → Boolean click to toggle source

Checks if dtype is a complex type

# File lib/nmatrix/nmatrix.rb, line 364
def complex_dtype?
  [:complex64, :complex128].include?(self.dtype)
end
concat(*m2) → NMatrix click to toggle source
concat(*m2, rank) → NMatrix
hconcat(*m2) → NMatrix
vconcat(*m2) → NMatrix
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

# File lib/nmatrix/nmatrix.rb, line 662
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_transpose → NMatrix click to toggle source

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.

# File lib/nmatrix/math.rb, line 624
def conjugate_transpose
  self.transpose.complex_conjugate!
end
corr() click to toggle source

Calculate the correlation matrix.

# File lib/nmatrix/math.rb, line 523
def corr
  raise NotImplementedError, "Does not work for complex dtypes" if complex_dtype?
  standard_deviation = std
  cov / (standard_deviation.transpose.dot(standard_deviation))
end
cov(opts={}) click to toggle source

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

# File lib/nmatrix/math.rb, line 510
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_pointer() click to toggle source

Returns the pointer to the matrix storage's data. This is useful primarily when you are using FFI with NMatrix – say, for example, you want to pass a float* to some function, and your NMatrix is a :float32 :dense matrix. Then you can call this function and get that pointer directly instead of copying the data.

static VALUE nm_data_pointer(VALUE self) {
  //if (NM_DTYPE(self) == nm::LIST_STORE)
  //  rb_warn("pointer requested for list storage, which may be meaningless");

  // This is actually pretty easy, since all of the storage types have their elements positioned in the same place
  // relative to one another. So yes, believe it or not, this should work just as well for Yale or list storage as for
  // dense.
  return INT2FIX(NM_STORAGE_DENSE(self)->elements);
}
dconcat(*matrices) click to toggle source

Depth concatenation with matrices.

# File lib/nmatrix/nmatrix.rb, line 727
def dconcat(*matrices)
  concat(*matrices, :layer)
end
default_value → ... click to toggle source

Get the default value for the matrix. For dense, this is undefined and will return Qnil. For list, it is user-defined. For yale, it's going to be some variation on zero, but may be Qfalse or Qnil.

static VALUE nm_default_value(VALUE self) {
  switch(NM_STYPE(self)) {
  case nm::YALE_STORE:
    return nm_yale_default_value(self);
  case nm::LIST_STORE:
    return nm_list_default_value(self);
  case nm::DENSE_STORE:
  default:
    return Qnil;
  }
}
dense? → true or false click to toggle source

Determine if m is a dense matrix.

# File lib/nmatrix/shortcuts.rb, line 41
def dense?; return stype == :dense; end
det → determinant click to toggle source

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.

# File lib/nmatrix/math.rb, line 454
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_exact() click to toggle source

Calculate the exact determinant of a dense matrix.

Returns nil for dense matrices which are not square or number of dimensions other than 2.

Note: Currently only implemented for 2x2 and 3x3 matrices.

static VALUE nm_det_exact(VALUE self) {

  if (NM_STYPE(self) != nm::DENSE_STORE) {
    rb_raise(rb_eNotImpError, "can only calculate exact determinant for dense matrices");
    return Qnil;
  }
  if (NM_DIM(self) != 2 || NM_SHAPE0(self) != NM_SHAPE1(self)) {
    rb_raise(nm_eShapeError, "matrices must be square to have a determinant defined");
    return Qnil;
  }

  NM_CONSERVATIVE(nm_register_value(&self));

  // Calculate the determinant and then assign it to the return value
  void* result = NM_ALLOCA_N(char, DTYPE_SIZES[NM_DTYPE(self)]);
  nm::dtype_t dtype = NM_DTYPE(self);
  nm_math_det_exact(NM_SHAPE0(self), NM_STORAGE_DENSE(self)->elements, NM_SHAPE0(self), NM_DTYPE(self), result);

  if (dtype == nm::RUBYOBJ) {
    nm_register_values(reinterpret_cast<VALUE*>(result), 1);
  }
  VALUE to_return = rubyobj_from_cval(result, NM_DTYPE(self)).rval;
  if (dtype == nm::RUBYOBJ) {
    nm_unregister_values(reinterpret_cast<VALUE*>(result), 1);
  }
  NM_CONSERVATIVE(nm_unregister_value(&self));

  return to_return;
}
diagonal(main_diagonal=true) click to toggle source

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

# 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
dim()
Alias for: dimensions
dim → Integer click to toggle source

Get the number of dimensions of a matrix.

In other words, if you set your matrix to be 3x4, the dim is 2. If the matrix was initialized as 3x4x3, the dim is 3.

Use effective_dim to get the dimension of an NMatrix which acts as a vector (e.g., a column or row).

static VALUE nm_dim(VALUE self) {
  return INT2FIX(NM_STORAGE(self)->dim);
}
Also aliased as: dim
dot(p1) click to toggle source

Matrix multiply (dot product): against another matrix or a vector.

For elementwise, use * instead.

The two matrices must be of the same stype (for now). If dtype differs, an upcast will occur.

static VALUE nm_multiply(VALUE left_v, VALUE right_v) {
  NM_CONSERVATIVE(nm_register_value(&left_v));
  NM_CONSERVATIVE(nm_register_value(&right_v));

  NMATRIX *left, *right;

  UnwrapNMatrix( left_v, left );

  if (NM_RUBYVAL_IS_NUMERIC(right_v)) {
    NM_CONSERVATIVE(nm_unregister_value(&left_v));
    NM_CONSERVATIVE(nm_unregister_value(&right_v));
    return matrix_multiply_scalar(left, right_v);
  }

  else if (TYPE(right_v) == T_ARRAY) {
    NM_CONSERVATIVE(nm_unregister_value(&left_v));
    NM_CONSERVATIVE(nm_unregister_value(&right_v));
    rb_raise(rb_eNotImpError, "please convert array to nx1 or 1xn NMatrix first");
  }

  else { // both are matrices (probably)
    CheckNMatrixType(right_v);
    UnwrapNMatrix( right_v, right );

    // work like vector dot product for 1dim
    if (left->storage->dim == 1 && right->storage->dim == 1) {
      if (left->storage->shape[0] != right->storage->shape[0]) {
        NM_CONSERVATIVE(nm_unregister_value(&left_v));
        NM_CONSERVATIVE(nm_unregister_value(&right_v));
        rb_raise(rb_eArgError, "The left- and right-hand sides of the operation must have the same dimensionality.");
      } else {
        VALUE result = elementwise_op(nm::EW_MUL, left_v, right_v);
        VALUE to_return = rb_funcall(result, rb_intern("sum"),0);
        NM_CONSERVATIVE(nm_unregister_value(&left_v));
        NM_CONSERVATIVE(nm_unregister_value(&right_v));
        return to_return;
      }
    }

    if (left->storage->shape[1] != right->storage->shape[0]) {
      NM_CONSERVATIVE(nm_unregister_value(&left_v));
      NM_CONSERVATIVE(nm_unregister_value(&right_v));
      rb_raise(rb_eArgError, "incompatible dimensions");
    }

    if (left->stype != right->stype) {
      NM_CONSERVATIVE(nm_unregister_value(&left_v));
      NM_CONSERVATIVE(nm_unregister_value(&right_v));
      rb_raise(rb_eNotImpError, "matrices must have same stype");
    }

    NM_CONSERVATIVE(nm_unregister_value(&left_v));
    NM_CONSERVATIVE(nm_unregister_value(&right_v));
    return matrix_multiply(left, right);

  }

  NM_CONSERVATIVE(nm_unregister_value(&left_v));
  NM_CONSERVATIVE(nm_unregister_value(&right_v));

  return Qnil;
}
Also aliased as: internal_dot
dtype → Symbol click to toggle source

Get the data type (dtype) of a matrix, e.g., :byte, :int8, :int16, :int32, :int64, :float32, :float64, :complex64, :complex128, or :object (the last is a Ruby object).

static VALUE nm_dtype(VALUE self) {
  ID dtype = rb_intern(DTYPE_NAMES[NM_DTYPE(self)]);
  return ID2SYM(dtype);
}
each → Enumerator click to toggle source

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).

# 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_along_dim(dimen=0, get_by=:reference)
Alias for: each_rank
each_column { |column| block } → NMatrix click to toggle source

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

# 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 -> { |column| block } → ... click to toggle source

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.

# 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_indices → Enumerator click to toggle source

Very similar to each_stored_with_indices. The key difference is that it enforces matrix ordering rather than storage ordering, which only matters if your matrix is Yale.

static VALUE nm_each_ordered_stored_with_indices(VALUE nmatrix) {
  NM_CONSERVATIVE(nm_register_value(&nmatrix));
  VALUE to_return = Qnil;

  switch(NM_STYPE(nmatrix)) {
  case nm::YALE_STORE:
    to_return = nm_yale_each_ordered_stored_with_indices(nmatrix);
    break;
  case nm::DENSE_STORE:
    to_return = nm_dense_each_with_indices(nmatrix);
    break;
  case nm::LIST_STORE:
    to_return = nm_list_each_with_indices(nmatrix, true);
    break;
  default:
    NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
    rb_raise(nm_eDataTypeError, "Not a proper storage type");
  }

  NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
  return to_return;
}
each_rank() → NMatrix click to toggle source
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.

@param [Fixnum] dimen the rank being iterated over.

# 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
Also aliased as: each_along_dim
each_row { |row| block } → NMatrix click to toggle source

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

# 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 → Enumerator click to toggle source

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

# 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_indices → Enumerator click to toggle source

Iterate over the stored entries of any matrix. For dense and yale, this iterates over non-zero entries; for list, this iterates over non-default entries. Yields dim+1 values for each entry: i, j, …, and the entry itself.

static VALUE nm_each_stored_with_indices(VALUE nmatrix) {
  NM_CONSERVATIVE(nm_register_value(&nmatrix));
  VALUE to_return = Qnil;

  switch(NM_STYPE(nmatrix)) {
  case nm::YALE_STORE:
    to_return = nm_yale_each_stored_with_indices(nmatrix);
    break;
  case nm::DENSE_STORE:
    to_return = nm_dense_each_with_indices(nmatrix);
    break;
  case nm::LIST_STORE:
    to_return = nm_list_each_with_indices(nmatrix, true);
    break;
  default:
    NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
    rb_raise(nm_eDataTypeError, "Not a proper storage type");
  }

  NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
  return to_return;
}
each_with_indices → Enumerator click to toggle source

Iterate over all entries of any matrix in standard storage order (as with each), and include the indices.

static VALUE nm_each_with_indices(VALUE nmatrix) {
  NM_CONSERVATIVE(nm_register_value(&nmatrix));
  VALUE to_return = Qnil;

  switch(NM_STYPE(nmatrix)) {
  case nm::YALE_STORE:
    to_return = nm_yale_each_with_indices(nmatrix);
    break;
  case nm::DENSE_STORE:
    to_return = nm_dense_each_with_indices(nmatrix);
    break;
  case nm::LIST_STORE:
    to_return = nm_list_each_with_indices(nmatrix, false);
    break;
  default:
    NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
    rb_raise(nm_eDataTypeError, "Not a proper storage type");
  }

  NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
  return to_return;
}
effective_dim()
effective_dim → Fixnum click to toggle source

Returns the number of dimensions that don't have length 1. Guaranteed to be less than or equal to dim.

static VALUE nm_effective_dim(VALUE self) {
  return INT2FIX(effective_dim(NM_STORAGE(self)));
}
Also aliased as: effective_dim
factorize_cholesky → [upper NMatrix, lower NMatrix] click to toggle source

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.

# 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 → ... click to toggle source

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.
# 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_map(&bl)

call-seq:

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

Maps using Enumerator (returns an Array or an Enumerator)

Alias for: map
float_dtype?() → Boolean click to toggle source

Checks if dtype is a floating point type

# File lib/nmatrix/nmatrix.rb, line 354
def float_dtype?
  [:float32, :float64].include?(dtype)
end
gesdd → [u, sigma, v_transpose] click to toggle source
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.

# File lib/nmatrix/math.rb, line 350
def gesdd(workspace_size=nil)
  self.clone.gesdd!(workspace_size)
end
gesdd! → [u, sigma, v_transpose] click to toggle source
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.

# File lib/nmatrix/math.rb, line 335
def gesdd!(workspace_size=nil)
  NMatrix::LAPACK::gesdd(self, workspace_size)
end
gesvd → [u, sigma, v_transpose] click to toggle source
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.

# File lib/nmatrix/math.rb, line 318
def gesvd(workspace_size=1)
  self.clone.gesvd!(workspace_size)
end
gesvd! → [u, sigma, v_transpose] click to toggle source
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.

# File lib/nmatrix/math.rb, line 304
def gesvd!(workspace_size=1)
  NMatrix::LAPACK::gesvd(self, workspace_size)
end
getrf!() click to toggle source
# File lib/nmatrix/lapacke.rb, line 169
def getrf!
  raise(StorageTypeError, "LAPACK functions only work on dense matrices") unless self.dense?

  ipiv = NMatrix::LAPACK::lapacke_getrf(:row, self.shape[0], self.shape[1], self, self.shape[1])

  return ipiv
end
hconcat(*matrices) click to toggle source

Horizontal concatenation with matrices.

# File lib/nmatrix/nmatrix.rb, line 717
def hconcat(*matrices)
  concat(*matrices, :column)
end
hermitian? → Boolean click to toggle source

Is this matrix hermitian?

Definition: en.wikipedia.org/wiki/Hermitian_matrix

For non-complex matrices, this function should return the same result as symmetric?.

static VALUE nm_hermitian(VALUE self) {
  return is_symmetric(self, true);
}
hessenberg() click to toggle source

Reduce self to upper hessenberg form using householder transforms.

References

# File lib/nmatrix/math.rb, line 241
def hessenberg
  clone.hessenberg!
end
hessenberg!() click to toggle source

Destructive version of hessenberg

# 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
index(value) click to toggle source

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.

# File lib/nmatrix/nmatrix.rb, line 966
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_copy(p1) click to toggle source

Copy constructor for no change of dtype or stype (used for initialize_copy hook).

static VALUE nm_init_copy(VALUE copy, VALUE original) {
  NM_CONSERVATIVE(nm_register_value(&copy));
  NM_CONSERVATIVE(nm_register_value(&original));

  NMATRIX *lhs, *rhs;

  CheckNMatrixType(original);

  if (copy == original) {
    NM_CONSERVATIVE(nm_unregister_value(&copy));
    NM_CONSERVATIVE(nm_unregister_value(&original));
    return copy;
  }

  UnwrapNMatrix( original, rhs );
  UnwrapNMatrix( copy,     lhs );

  lhs->stype = rhs->stype;

  // Copy the storage
  CAST_TABLE(ttable);
  lhs->storage = ttable[lhs->stype][rhs->stype](rhs->storage, rhs->storage->dtype, NULL);

  NM_CONSERVATIVE(nm_unregister_value(&copy));
  NM_CONSERVATIVE(nm_unregister_value(&original));

  return copy;
}
inject → symbol click to toggle source

This overrides the inject function to use #map_stored for yale matrices

Calls superclass method
# File lib/nmatrix/nmatrix.rb, line 958
def inject(sym)
  return super(sym) unless self.yale?
  return self.map_stored.inject(sym)
end
inject_along_dim(dimen=0, initial=nil, dtype=nil)
Alias for: inject_rank
inject_rank() → Enumerator click to toggle source
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.

@param [Integer] dimen the dimension being reduced @param [Numeric] initial 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.

@param [Symbol] dtype if non-nil/false, forces the accumulated result to have this dtype @return [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.
# 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
integer_dtype?() → Boolean click to toggle source

Checks if dtype is an integer type

# File lib/nmatrix/nmatrix.rb, line 345
def integer_dtype?
  [:byte, :int8, :int16, :int32, :int64].include?(self.dtype)
end
internal_dot(p1)
Alias for: dot
inverse()
Alias for: invert
invert → NMatrix click to toggle source

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

  • Returns :

    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.

# 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
Also aliased as: inverse
invert!() click to toggle source
# File lib/nmatrix/atlas.rb, line 186
def invert!
  raise(StorageTypeError, "invert only works on dense matrices currently") unless self.dense?
  raise(ShapeError, "Cannot invert non-square matrix") unless shape[0] == shape[1]
  raise(DataTypeError, "Cannot invert an integer matrix in-place") if self.integer_dtype?

  # Even though we are using the ATLAS plugin, we still might be missing
  # CLAPACK (and thus clapack_getri) if we are on OS X.
  if NMatrix.has_clapack?
    # Get the pivot array; factor the matrix
    # We can't used getrf! here since it doesn't have the clapack behavior,
    # so it doesn't play nicely with clapack_getri
    n = self.shape[0]
    pivot = NMatrix::LAPACK::clapack_getrf(:row, n, n, self, n)
    # Now calculate the inverse using the pivot array
    NMatrix::LAPACK::clapack_getri(:row, n, self, n, pivot)
    self
  else
    __inverse__(self,true)
  end
end
is_ref?() click to toggle source

Check to determine whether matrix is a reference to another matrix.

static VALUE nm_is_ref(VALUE self) {
  if (NM_SRC(self) == NM_STORAGE(self)) return Qfalse;
  return Qtrue;
}
kron_prod(mat) click to toggle source

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] ]
# File lib/nmatrix/math.rb, line 587
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) → NMatrix click to toggle source

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.

# File lib/nmatrix/math.rb, line 425
def laswp(ary, opts={})
  self.clone.laswp!(ary, opts)
end
Also aliased as: permute_columns
laswp!(ary) → NMatrix click to toggle source

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.

# File lib/nmatrix/math.rb, line 378
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
Also aliased as: permute_columns!
layer(layer_number) → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 853
def layer(layer_number, get_by = :copy)
  rank(2, layer_number, get_by)
end
list? → true or false click to toggle source

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

# File lib/nmatrix/shortcuts.rb, line 53
def list?;  return stype == :list; end
log(*args) click to toggle source
static VALUE nm_unary_log(int argc, VALUE* argv, VALUE self) {
  NM_CONSERVATIVE(nm_register_values(argv, argc));
  const double default_log_base = exp(1.0);
  NMATRIX* left;
  UnwrapNMatrix(self, left);
  std::string sym;

  switch(left->stype) {
  case nm::DENSE_STORE:
    sym = "__dense_unary_log__";
    break;
  case nm::YALE_STORE:
    sym = "__yale_unary_log__";
    break;
  case nm::LIST_STORE:
    sym = "__list_unary_log__";
    break;
  }
  NM_CONSERVATIVE(nm_unregister_values(argv, argc));
  if (argc > 0) { //supplied a base
    return rb_funcall(self, rb_intern(sym.c_str()), 1, argv[0]);
  }
  return rb_funcall(self, rb_intern(sym.c_str()), 1, nm::RubyObject(default_log_base).rval);
}
lower_triangle → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 800
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
Also aliased as: tril
lower_triangle! → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 829
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
Also aliased as: tril!
map → Enumerator click to toggle source
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.

# 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
Also aliased as: flat_map
map! → Enumerator click to toggle source
map! { |elem| block } → NMatrix

Maps in place. @see map

# 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_stored → Enumerator click to toggle source

Iterate over the stored entries of any matrix. For dense and yale, this iterates over non-zero entries; for list, this iterates over non-default entries. Yields dim+1 values for each entry: i, j, …, and the entry itself.

static VALUE nm_map_stored(VALUE nmatrix) {
  NM_CONSERVATIVE(nm_register_value(&nmatrix));
  VALUE to_return = Qnil;

  switch(NM_STYPE(nmatrix)) {
  case nm::YALE_STORE:
    to_return = nm_yale_map_stored(nmatrix);
    break;
  case nm::DENSE_STORE:
    to_return = nm_dense_map(nmatrix);
    break;
  case nm::LIST_STORE:
    to_return = nm_list_map_stored(nmatrix, Qnil);
    break;
  default:
    NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
    rb_raise(nm_eDataTypeError, "Not a proper storage type");
  }

  NM_CONSERVATIVE(nm_unregister_value(&nmatrix));
  return to_return;
}
max() → NMatrix click to toggle source
max(dimen) → NMatrix

Calculates the maximum along the specified dimension.

@see inject_rank

# File lib/nmatrix/math.rb, line 712
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() → NMatrix click to toggle source
mean(dimen) → NMatrix

Calculates the mean along the specified dimension.

This will force integer types to float64 dtype.

@see inject_rank

# File lib/nmatrix/math.rb, line 659
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() → NMatrix click to toggle source
min(dimen) → NMatrix

Calculates the minimum along the specified dimension.

@see inject_rank

# File lib/nmatrix/math.rb, line 693
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
norm2(incx=1, n=nil)
Alias for: nrm2
norm2 → Numeric click to toggle source

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.

# File lib/nmatrix/math.rb, line 824
def nrm2 incx=1, n=nil
  return method_missing(:nrm2, incx, n) unless vector?
  NMatrix::BLAS::nrm2(self, incx, self.size / incx)
end
Also aliased as: norm2
nvector? → true or false click to toggle source

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.

# File lib/nmatrix/nmatrix.rb, line 429
def nvector?
  self.effective_dim < self.dim
end
object_dtype?() → Boolean click to toggle source

Checks if dtype is a ruby object

# File lib/nmatrix/nmatrix.rb, line 374
def object_dtype?
  dtype == :object
end
offset → Array click to toggle source

Get the offset (slice position) of a matrix. Typically all zeros, unless you have a reference slice.

static VALUE nm_offset(VALUE self) {
  NM_CONSERVATIVE(nm_register_value(&self));
  STORAGE* s   = NM_STORAGE(self);

  // Copy elements into a VALUE array and then use those to create a Ruby array with rb_ary_new4.
  VALUE* offset = NM_ALLOCA_N(VALUE, s->dim);
  nm_register_values(offset, s->dim);
  for (size_t index = 0; index < s->dim; ++index)
    offset[index] = INT2FIX(s->offset[index]);

  nm_unregister_values(offset, s->dim);
  NM_CONSERVATIVE(nm_unregister_value(&self));
  return rb_ary_new4(s->dim, offset);
}
permute_columns(ary, opts={})
Alias for: laswp
permute_columns!(ary, opts={})
Alias for: laswp!
potrf!(which) click to toggle source
# File lib/nmatrix/atlas.rb, line 207
def potrf!(which)
  raise(StorageTypeError, "ATLAS functions only work on dense matrices") unless self.dense?
  raise(ShapeError, "Cholesky decomposition only valid for square matrices") unless self.dim == 2 && self.shape[0] == self.shape[1]

  NMatrix::LAPACK::clapack_potrf(:row, which, self.shape[0], self, self.shape[1])
end
potrf_lower!() click to toggle source
# File lib/nmatrix/math.rb, line 188
def potrf_lower!
  potrf! :lower
end
potrf_upper!() click to toggle source
# File lib/nmatrix/math.rb, line 184
def potrf_upper!
  potrf! :upper
end
pow(n) click to toggle source

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
# File lib/nmatrix/math.rb, line 542
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
quaternion → NMatrix click to toggle source

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]
# 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(dimension, row_or_column_number) → NMatrix click to toggle source
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)

# 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
reduce_along_dim(dimen=0, initial=nil, dtype=nil)
Alias for: inject_rank
repeat(count, axis) → NMatrix click to toggle source
  • 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]]
    
# File lib/nmatrix/nmatrix.rb, line 1008
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) → NMatrix click to toggle source

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.

# 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) → NMatrix click to toggle source
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.

# 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
round(*args) click to toggle source
static VALUE nm_unary_round(int argc, VALUE* argv, VALUE self) {
  NM_CONSERVATIVE(nm_register_values(argv, argc));
  const int default_precision = 0;
  NMATRIX* left;
  UnwrapNMatrix(self, left);
  std::string sym;

  switch(left->stype) {
  case nm::DENSE_STORE:
    sym = "__dense_unary_round__";
    break;
  case nm::YALE_STORE:
    sym = "__yale_unary_round__";
    break;
  case nm::LIST_STORE:
    sym = "__list_unary_round__";
    break;
  }
  NM_CONSERVATIVE(nm_unregister_values(argv, argc));
  if (argc > 0) { //supplied precision
    return rb_funcall(self, rb_intern(sym.c_str()), 1, argv[0]);
  }
  return rb_funcall(self, rb_intern(sym.c_str()), 1, nm::RubyObject(default_precision).rval);
}
row(row_number) → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 533
def row(row_number, get_by = :copy)
  rank(0, row_number, get_by)
end
rows → Integer click to toggle source

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

# File lib/nmatrix/nmatrix.rb, line 256
def rows
  shape[0]
end
shape → Array click to toggle source

Get the shape (dimensions) of a matrix.

static VALUE nm_shape(VALUE self) {
  NM_CONSERVATIVE(nm_register_value(&self));
  STORAGE* s   = NM_STORAGE(self);

  // Copy elements into a VALUE array and then use those to create a Ruby array with rb_ary_new4.
  VALUE* shape = NM_ALLOCA_N(VALUE, s->dim);
  nm_register_values(shape, s->dim);
  for (size_t index = 0; index < s->dim; ++index)
    shape[index] = INT2FIX(s->shape[index]);

  nm_unregister_values(shape, s->dim);
  NM_CONSERVATIVE(nm_unregister_value(&self));
  return rb_ary_new4(s->dim, shape);
}
shuffle → ... click to toggle source
shuffle(rng) → ...

Re-arranges the contents of an NVector.

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

# File lib/nmatrix/nmatrix.rb, line 886
def shuffle(*args)
  method_missing(:shuffle!, *args) if self.effective_dim > 1
  t = self.clone
  t.shuffle!(*args)
end
shuffle! → ... click to toggle source
shuffle!(random: rng) → ...

Re-arranges the contents of an NVector.

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

# File lib/nmatrix/nmatrix.rb, line 868
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
size → Fixnum click to toggle source

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

# File lib/nmatrix/nmatrix.rb, line 413
def size
  NMatrix.size(self.shape)
end
slice → ... click to toggle source

Access the contents of an NMatrix at given coordinates, using copying.

n.slice(3,3)  # => 5.0
n.slice(0..1,0..1) #=> matrix [2,2]
static VALUE nm_mget(int argc, VALUE* argv, VALUE self) {
  static void* (*ttable[nm::NUM_STYPES])(const STORAGE*, SLICE*) = {
    nm_dense_storage_get,
    nm_list_storage_get,
    nm_yale_storage_get
  };
  nm::stype_t stype = NM_STYPE(self);
  return nm_xslice(argc, argv, ttable[stype], nm_delete, self);
}
solve(b) click to toggle source
# File lib/nmatrix/lapacke.rb, line 198
def solve b
  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?

  x     = b.clone
  clone = self.clone
  n = self.shape[0]
  ipiv = NMatrix::LAPACK.lapacke_getrf(:row, n, n, clone, n)
  NMatrix::LAPACK.lapacke_getrs(:row, :no_transpose, n, b.shape[1], clone, n, ipiv, x, b.shape[1])
  x
end
sorted_indices → Array click to toggle source

Returns an array of the indices ordered by value sorted.

# File lib/nmatrix/nmatrix.rb, line 899
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
std() → NMatrix click to toggle source
std(dimen) → NMatrix

Calculates the sample standard deviation along the specified dimension.

This will force integer types to float64 dtype.

@see inject_rank

# File lib/nmatrix/math.rb, line 757
def std(dimen=0)
  variance(dimen).sqrt
end
stype → Symbol click to toggle source

Get the storage type (stype) of a matrix, e.g., :yale, :dense, or :list.

static VALUE nm_stype(VALUE self) {
  NM_CONSERVATIVE(nm_register_value(&self));
  VALUE stype = ID2SYM(rb_intern(STYPE_NAMES[NM_STYPE(self)]));
  NM_CONSERVATIVE(nm_unregister_value(&self));
  return stype;
}
sum() → NMatrix click to toggle source
sum(dimen) → NMatrix

Calculates the sum along the specified dimension.

@see inject_rank

# File lib/nmatrix/math.rb, line 677
def sum(dimen=0)
  inject_rank(dimen, 0.0) do |sum, sub_mat|
    sum + sub_mat
  end
end
supershape → Array click to toggle source

Get the shape of a slice's parent.

static VALUE nm_supershape(VALUE self) {

  STORAGE* s   = NM_STORAGE(self);
  if (s->src == s) {
    return nm_shape(self); // easy case (not a slice)
  }
  else s = s->src;

  NM_CONSERVATIVE(nm_register_value(&self));

  VALUE* shape = NM_ALLOCA_N(VALUE, s->dim);
  nm_register_values(shape, s->dim);
  for (size_t index = 0; index < s->dim; ++index)
    shape[index] = INT2FIX(s->shape[index]);

  nm_unregister_values(shape, s->dim);
  NM_CONSERVATIVE(nm_unregister_value(&self));
  return rb_ary_new4(s->dim, shape);
}
symmetric? → Boolean click to toggle source

Is this matrix symmetric?

static VALUE nm_symmetric(VALUE self) {
  return is_symmetric(self, false);
}
to_a → Array click to toggle source

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

# 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_f → Float click to toggle source

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.

# 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_a()
Alias for: to_flat_array
to_flat_array → Array click to toggle source
to_flat_a → Array

Converts an NMatrix to a one-dimensional Ruby Array.

# 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
Also aliased as: to_flat_a
to_h()
Alias for: to_hash
to_hash → Hash click to toggle source

Create a Ruby Hash from an NMatrix.

# 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
Also aliased as: to_h
trace → Numeric click to toggle source

Calculates the trace of an nxn matrix.

  • Raises :

    • ShapeError -> Expected square matrix

  • Returns :

    • The trace of the matrix (a numeric value)

# File lib/nmatrix/math.rb, line 640
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 → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 602
def transpose(permute = nil)
  if self.dim == 1
    return self.clone
  elsif self.dim == 2
    new_shape = [self.shape[1], self.shape[0]]
  elsif permute.nil?
    raise(ArgumentError, "need permutation array of size #{self.dim}")
  elsif 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
tril(k = 0)
Alias for: lower_triangle
tril!(k = 0)
Alias for: lower_triangle!
triu(k = 0)
Alias for: upper_triangle
triu!(k = 0)
Alias for: upper_triangle!
upper_triangle → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 745
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
Also aliased as: triu
upper_triangle! → NMatrix click to toggle source
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.

# File lib/nmatrix/nmatrix.rb, line 774
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
Also aliased as: triu!
variance() → NMatrix click to toggle source
variance(dimen) → NMatrix

Calculates the sample variance along the specified dimension.

This will force integer types to float64 dtype.

@see inject_rank

# File lib/nmatrix/math.rb, line 734
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) click to toggle source

Vertical concatenation with matrices.

# File lib/nmatrix/nmatrix.rb, line 722
def vconcat(*matrices)
  concat(*matrices, :row)
end
vector? → true or false click to toggle source

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

# File lib/nmatrix/nmatrix.rb, line 439
def vector?
  self.effective_dim == 1
end
write(*args) click to toggle source

Binary file writer for NMatrix standard format. file should be a path, which we aren't going to check very carefully (in other words, this function should generally be called from a Ruby helper method). Function also takes a symmetry argument, which allows us to specify that we only want to save the upper triangular portion of the matrix (or if the matrix is a lower triangular matrix, only the lower triangular portion). nil means regular storage.

static VALUE nm_write(int argc, VALUE* argv, VALUE self) {
  using std::ofstream;

  if (argc < 1 || argc > 2) {
    rb_raise(rb_eArgError, "Expected one or two arguments");
  }

  NM_CONSERVATIVE(nm_register_values(argv, argc));
  NM_CONSERVATIVE(nm_register_value(&self));

  VALUE file = argv[0],
        symm = argc == 1 ? Qnil : argv[1];

  NMATRIX* nmatrix;
  UnwrapNMatrix( self, nmatrix );

  nm::symm_t symm_ = interpret_symm(symm);

  if (nmatrix->storage->dtype == nm::RUBYOBJ) {
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&self));
    rb_raise(rb_eNotImpError, "Ruby Object writing is not implemented yet");
  }

  // Get the dtype, stype, itype, and symm and ensure they're the correct number of bytes.
  uint8_t st = static_cast<uint8_t>(nmatrix->stype),
          dt = static_cast<uint8_t>(nmatrix->storage->dtype),
          sm = static_cast<uint8_t>(symm_);
  uint16_t dim = nmatrix->storage->dim;

  //FIXME: Cast the matrix to the smallest possible index type. Write that in the place of IType.

  // Check arguments before starting to write.
  if (nmatrix->stype == nm::LIST_STORE) {
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&self));
    rb_raise(nm_eStorageTypeError, "cannot save list matrix; cast to yale or dense first");
  }
  if (symm_ != nm::NONSYMM) {
    NM_CONSERVATIVE(nm_unregister_values(argv, argc));
    NM_CONSERVATIVE(nm_unregister_value(&self));

    if (dim != 2) rb_raise(rb_eArgError, "symmetry/triangularity not defined for a non-2D matrix");
    if (nmatrix->storage->shape[0] != nmatrix->storage->shape[1])
      rb_raise(rb_eArgError, "symmetry/triangularity not defined for a non-square matrix");
    if (symm_ == nm::HERM &&
          dt != static_cast<uint8_t>(nm::COMPLEX64) && dt != static_cast<uint8_t>(nm::COMPLEX128) && dt != static_cast<uint8_t>(nm::RUBYOBJ))
      rb_raise(rb_eArgError, "cannot save a non-complex matrix as hermitian");
  }

  ofstream f(RSTRING_PTR(file), std::ios::out | std::ios::binary);

  // Get the NMatrix version information.
  uint16_t major, minor, release, null16 = 0;
  get_version_info(major, minor, release);

  // WRITE FIRST 64-BIT BLOCK
  f.write(reinterpret_cast<const char*>(&major),   sizeof(uint16_t));
  f.write(reinterpret_cast<const char*>(&minor),   sizeof(uint16_t));
  f.write(reinterpret_cast<const char*>(&release), sizeof(uint16_t));
  f.write(reinterpret_cast<const char*>(&null16),  sizeof(uint16_t));

  uint8_t ZERO = 0;
  // WRITE SECOND 64-BIT BLOCK
  f.write(reinterpret_cast<const char*>(&dt), sizeof(uint8_t));
  f.write(reinterpret_cast<const char*>(&st), sizeof(uint8_t));
  f.write(reinterpret_cast<const char*>(&ZERO),sizeof(uint8_t));
  f.write(reinterpret_cast<const char*>(&sm), sizeof(uint8_t));
  f.write(reinterpret_cast<const char*>(&null16), sizeof(uint16_t));
  f.write(reinterpret_cast<const char*>(&dim), sizeof(uint16_t));

  // Write shape (in 64-bit blocks)
  write_padded_shape(f, nmatrix->storage->dim, nmatrix->storage->shape);

  if (nmatrix->stype == nm::DENSE_STORE) {
    write_padded_dense_elements(f, reinterpret_cast<DENSE_STORAGE*>(nmatrix->storage), symm_, nmatrix->storage->dtype);
  } else if (nmatrix->stype == nm::YALE_STORE) {
    YALE_STORAGE* s = reinterpret_cast<YALE_STORAGE*>(nmatrix->storage);
    uint32_t ndnz   = s->ndnz,
             length = nm_yale_storage_get_size(s);
    f.write(reinterpret_cast<const char*>(&ndnz),   sizeof(uint32_t));
    f.write(reinterpret_cast<const char*>(&length), sizeof(uint32_t));

    write_padded_yale_elements(f, s, length, symm_, s->dtype);
  }

  f.close();

  NM_CONSERVATIVE(nm_unregister_values(argv, argc));
  NM_CONSERVATIVE(nm_unregister_value(&self));

  return Qtrue;
}
yale? → true or false click to toggle source

Determine if m is a Yale matrix.

# File lib/nmatrix/shortcuts.rb, line 47
def yale?;  return stype == :yale; end

Protected Instance Methods

__dense_unary_round__(precision) click to toggle source
# File lib/nmatrix/math.rb, line 921
def __dense_unary_round__(precision)
  if self.complex_dtype?
    self.__dense_map__ { |l| Complex(l.real.round(precision), l.imag.round(precision)) }
                                  .cast(stype, dtype)
  else
    self.__dense_map__ { |l| l.round(precision) }.cast(stype, dtype)
  end
end
__list_unary_round__(precision) click to toggle source

These are for rounding each value of a matrix. Takes an optional argument

# File lib/nmatrix/math.rb, line 903
def __list_unary_round__(precision)
  if self.complex_dtype?
    self.__list_map_stored__(nil) { |l| Complex(l.real.round(precision), l.imag.round(precision)) }
                                  .cast(stype, dtype)
  else
    self.__list_map_stored__(nil) { |l| l.round(precision) }.cast(stype, dtype)
  end
end
__yale_unary_round__(precision) click to toggle source
# File lib/nmatrix/math.rb, line 912
def __yale_unary_round__(precision)
  if self.complex_dtype?
    self.__yale_map_stored__ { |l| Complex(l.real.round(precision), l.imag.round(precision)) }
                                  .cast(stype, dtype)
  else
    self.__yale_map_stored__ { |l| l.round(precision) }.cast(stype, dtype)
  end
end
dtype_for_floor_or_ceil() click to toggle source

These are for calculating the floor or ceil of matrix

# File lib/nmatrix/math.rb, line 931
def dtype_for_floor_or_ceil
  if self.integer_dtype? or [:complex64, :complex128, :object].include?(self.dtype)
    return_dtype = dtype
  elsif [:float32, :float64].include?(self.dtype)
    return_dtype = :int64
  end

  return_dtype
end