class Matrix

Attributes

rows[R]
wrap[R]

Public Class Methods

build(n,m,val=0) { |i, j| ... } click to toggle source

Creates a matrix n x m If you provide a block, it will be used to set the values. If not, val will be used

# File lib/extendmatrix.rb, line 320
def build(n,m,val=0)
  f = (block_given?)? lambda {|i,j| yield(i, j)} : lambda {|i,j| val}
  Matrix.rows((0...n).collect {|i| (0...m).collect {|j| f.call(i,j)}})
end
diag(*args) click to toggle source

Creates a matrix with the given matrices as diagonal blocks

# File lib/extendmatrix.rb, line 328
def diag(*args)
  dsize = 0
  sizes = args.collect{|e| x = (e.is_a?(Matrix)) ? e.row_size : 1; dsize += x; x}
  m = Matrix.zero(dsize)
  count = 0

  sizes.size.times{|i|
    range = count..(count+sizes[i]-1)
    m[range, range] = args[i]
    count += sizes[i]
  }
  m
end
diag_in_delta?(m1, m0, eps = 1.0e-10) click to toggle source

Tests if all the diagonal elements of two matrix are equal in delta

# File lib/extendmatrix.rb, line 359
def diag_in_delta?(m1, m0, eps = 1.0e-10)
  n = m1.row_size
  return false if n != m0.row_size or m1.column_size != m0.column_size
  n.times{|i|
    return false if (m1[i,i]-m0[i,i]).abs > eps
  }
  true
end
equal_in_delta?(m0, m1, delta = 1.0e-10) click to toggle source

Tests if all the elements of two matrix are equal in delta

# File lib/extendmatrix.rb, line 345
def equal_in_delta?(m0, m1, delta = 1.0e-10)
  delta = delta.abs
  m0.row_size.times {|i|
    m0.column_size.times {|j|
      x=m0[i,j]; y=m1[i,j]
      return false if (x < y - delta or x > y + delta)
    }
  }
  true
end
robust_I(n) click to toggle source

Return an empty matrix on n=0 else, returns I(n)

# File lib/extendmatrix.rb, line 226
def self.robust_I(n)
  n>0 ? self.I(n) : self.[]
end

Public Instance Methods

/(v)

quo seems always desirable

Also aliased as: old_divition
Alias for: quo
[](i, j) click to toggle source

Return a value or a vector/matrix of values depending if the indexes are ranges or not m = Matrix.new_with_value(4, 3){|i, j| i * 3 + j} m: 0 1 2

3  4  5
6  7  8
9 10 11

m[1, 2] => 5 m => Vector[10, 11] m[0..1, 0..2] => Matrix[[0, 1, 2], [3, 4, 5]]

# File lib/extendmatrix.rb, line 241
def [](i, j)
  case i
  when Range
    case j
    when Range
      Matrix[*i.collect{|l| self.row(l)[j].to_a}]
    else
      column(j)[i]
    end
  else
    case j
    when Range
      row(i)[j]
    else
      ids(i, j)
    end
  end
end
Also aliased as: ids
[]=(i, j, v) click to toggle source

Set the values of a matrix m = Matrix.build(3, 3){|i, j| i * 3 + j} m: 0 1 2

3  4  5
6  7  8

m[1, 2] = 9 => Matrix[[0, 1, 2], [3, 4, 9], [6, 7, 8]] m = Vector[8, 8] => Matrix[[0, 1, 2], [3, 8, 8], [6, 7, 8]] m[0..1, 0..1] = Matrix[[0, 0, 0],[0, 0, 0]]

=> Matrix[[0, 0, 2], [0, 0, 8], [6, 7, 8]]
# File lib/extendmatrix.rb, line 270
def []=(i, j, v)
  case i
  when Range
    if i.entries.size == 1
      self[i.begin, j] = (v.is_a?(Matrix) ? v.row(0) : v)
    else
      case j
      when Range
        if j.entries.size == 1
          self[i, j.begin] = (v.is_a?(Matrix) ? v.column(0) : v)
        else
          i.each{|l| self.row= l, v.row(l - i.begin), j}
        end
      else
        self.column= j, v, i
      end
    end
  else
    case j
    when Range
      if j.entries.size == 1
        self[i, j.begin] = (v.is_a?(Vector) ? v[0] : v)
      else
        self.row= i, v, j
      end
    else
      @rows[i][j] = (v.is_a?(Vector) ? v[0] : v)

    end
  end
end
cols_len() click to toggle source

Returns a list with the maximum lengths

# File lib/extendmatrix.rb, line 429
def cols_len
  (0...column_size).collect {|j| max_len_column(j)}
end
dup() click to toggle source

Return a duplicate matrix, with all elements copied

# File lib/extendmatrix.rb, line 304
def dup
  (self.class).rows(self.rows.dup)
end
each() { |e| ... } click to toggle source

Iterate the elements of a matrix

# File lib/extendmatrix.rb, line 463
def each
  @rows.each {|x| x.each {|e| yield(e)}}
  nil
end
ids(i, j)
Alias for: []
initialize_copy(orig) click to toggle source
# File lib/extendmatrix.rb, line 308
def initialize_copy(orig)
  init_rows(orig.rows, true)
  self.wrap=(orig.wrap)
end
max_len_column(j) click to toggle source

Returns the maximum length of column elements

# File lib/extendmatrix.rb, line 422
def max_len_column(j)
  column_collect(j) {|x| x.to_s.length}.max
end
old_divition(v)
Alias for: /
quo(v) click to toggle source

Returns the matrix divided by a scalar

# File lib/extendmatrix.rb, line 372
def quo(v)
  map {|e| e.quo(v)}
end
Also aliased as: /
set(m) click to toggle source

Set de values of a matrix and the value of wrap property

# File lib/extendmatrix.rb, line 385
def set(m)
  0.upto(m.row_size - 1) do |i|
    0.upto(m.column_size - 1) do |j|
      self[i, j] = m[i, j]
    end
  end
  self.wrap = m.wrap
end
to_s(mode = :pretty, len_col = 3) click to toggle source

Returns a string for nice printing matrix

# File lib/extendmatrix.rb, line 437
def to_s(mode = :pretty, len_col = 3)
  return "Matrix[]" if empty?
  if mode == :pretty
    clen = cols_len
    to_a.collect {|r|
      i=0
      r.map {|x|
        l=clen[i]
        i+=1
        format("%#{l}s ", x.to_s)
      } << "\n"
    }.join("")
  else
    i = 0; s = ""; cs = column_size
    each do |e|
      i = (i + 1) % cs
      s += format("%#{len_col}s ", e.to_s)
      s += "\n" if i == 0
    end
    s
  end
end
Also aliased as: to_s_old
to_s_old(mode = :pretty, len_col = 3)
Alias for: to_s
wrap=(mode = :torus) click to toggle source

Set wrap feature of a matrix

# File lib/extendmatrix.rb, line 409
def wrap=(mode = :torus)
  case mode
  when :torus then eval(wraplate("i %= row_size; j %= column_size"))
  when :h_cylinder then eval(wraplate("i %= row_size"))
  when :v_cylinder then eval(wraplate("j %= column_size"))
  when :nil then eval(wraplate)
  end
  @wrap = mode
end
wraplate(ijwrap = "") click to toggle source
# File lib/extendmatrix.rb, line 394
def wraplate(ijwrap = "")
  "class << self
    def [](i, j)
    #{ijwrap}; @rows[i][j]
    end

    def []=(i, j, v)
    #{ijwrap}; @rows[i][j] = v
    end
  end"
end

Advanced methods

↑ top

Constants

EXTENSION_VERSION

Public Instance Methods

L() click to toggle source

Return the lower triangular matrix of LU factorization L = M_1^{-1} * … * M_{n-1}^{-1} = I + sum_{k=1}^{n-1} tau * e

# File lib/extendmatrix.rb, line 744
def L
  LU.factorization(self)[0]
end
U() click to toggle source

Return the upper triangular matrix of LU factorization M_{n-1} * … * M_1 * A = U

# File lib/extendmatrix.rb, line 737
def U
  LU.factorization(self)[1]
end
bidiagonal() click to toggle source

Returns the upper bidiagonal matrix obtained with Householder Bidiagonalization algorithm

# File lib/extendmatrix.rb, line 825
def bidiagonal
  ub, vb = Householder.bidiag(self)
  ub.t * self * vb
end
cJacobi(tol = 1.0e-10) click to toggle source

Classical Jacobi 8.4.3 Golub & van Loan

# File lib/extendmatrix.rb, line 1074
def cJacobi(tol = 1.0e-10)
  a = self.clone
  n = row_size
  v = Matrix.I(n)
  eps = tol * a.normF
  while Jacobi.off(a) > eps
    p, q = Jacobi.max(a)
    c, s = Jacobi.sym_schur2(a, p, q)
    #print "\np:#{p} q:#{q} c:#{c} s:#{s}\n"
    j = Jacobi.J(p, q, c, s, n)
    a = j.t * a * j
    v = v * j
  end
  return a, v
end
cJacobiA(tol = 1.0e-10) click to toggle source

Returns the aproximation matrix computed with Classical Jacobi algorithm. The aproximate eigenvalues values are in the diagonal of the matrix A.

# File lib/extendmatrix.rb, line 1094
def cJacobiA(tol = 1.0e-10)
  cJacobi(tol)[0]
end
cJacobiV(tol = 1.0e-10) click to toggle source

Returns the orthogonal matrix obtained with the Jacobi eigenvalue algorithm. The columns of V are the eigenvector.

# File lib/extendmatrix.rb, line 1111
def cJacobiV(tol = 1.0e-10)
  cJacobi(tol)[1]
end
column!(j) { |rows[j]| ... } click to toggle source

Returns column vector number “j” as a Vector. When the block is given, the elements of column “j” are mmodified

# File lib/extendmatrix.rb, line 591
def column!(j)
  if block_given?
    (0...row_size).collect { |i| @rows[i][j] = yield @rows[i][j] }
  else
    column(j)
  end
end
Also aliased as: column_collect!
column2matrix(c) click to toggle source

Returns the column/s of matrix as a Matrix

# File lib/extendmatrix.rb, line 665
def column2matrix(c)
  if c.is_a?(Range)
    a=c.map {|v| column(v).to_a}.find_all {|v| v.size>0}
  else
    a = column(c).to_a
  end
  if c.is_a?(Range)
    return Matrix.columns(a)
  else
    return Matrix[*a.collect{|x| [x]}]
  end
end
column=(args) click to toggle source

Set a certain column with the values of a Vector m = Matrix.build(3, 3){|i, j| i * 3 + j + 1} m.column= 1, Vector[1, 1, 1], 1..2 m => 1 2 3

4 1 6
7 1 9
# File lib/extendmatrix.rb, line 608
def column=(args)
  m = row_size
  c, v, r = MMatrix.id_vect_range(args, m)
  (m..r.begin - 1).each{|i| self[i, c] = 0}
  [v.size, r.entries.size].min.times{|i| self[i + r.begin, c] = v[i]}
  ((v.size + r.begin)..r.entries.last).each {|i| self[i, c] = 0}
end
column_collect(j, &block) click to toggle source

Returns an array with the elements collected from the column “j”. When a block is given, the elements of that vector are iterated.

# File lib/extendmatrix.rb, line 582
def column_collect(j, &block)
  f = MMatrix.default_block(block)
  (0...row_size).collect {|r| f.call(self[r, j])}
end
column_collect!(j)
Alias for: column!
column_sum() click to toggle source

Returns the marginal of columns

# File lib/extendmatrix.rb, line 685
def column_sum
  (0...column_size).collect {|i|
    column(i).sum
  }
end
eigenvaluesJacobi() click to toggle source

Returns a Vector with the eigenvalues aproximated values. The eigenvalues are computed with the Classic Jacobi Algorithm.

# File lib/extendmatrix.rb, line 1102
def eigenvaluesJacobi
  a = cJacobiA
  Vector[*(0...row_size).collect{|i| a[i, i]}]
end
empty?() click to toggle source

Tests if the matrix is empty or not

# File lib/extendmatrix.rb, line 641
def empty?
  @rows.empty? if @rows
end
givensQ() click to toggle source

Returns the orthogonal matrix Q of Givens QR factorization. Q = G_1 * … * G_t where G_j is the j’th Givens rotation and ‘t’ is the total number of rotations

# File lib/extendmatrix.rb, line 941
def givensQ
  Givens.QR(self)[1]
end
givensR() click to toggle source

Returns the upper triunghiular matrix R of a Givens QR factorization

# File lib/extendmatrix.rb, line 932
def givensR
  Givens.QR(self)[0]
end
gram_schmidt() click to toggle source

Modified Gram Schmidt QR factorization (MC, Golub, p. 232) A = Q_1 * R_1

# File lib/extendmatrix.rb, line 856
def gram_schmidt
  a = clone
  n = column_size
  m = row_size
  q = Matrix.build(m, n){0}
  r = Matrix.zero(n)
  for k in 0...n
    r[k,k] = a[0...m, k].norm
    q[0...m, k] = a[0...m, k] / r[k, k]
    for j in (k+1)...n
      r[k, j] = q[0...m, k].t * a[0...m, j]
      a[0...m, j] -= q[0...m, k] * r[k, j]
    end
  end
  return q, r
end
gram_schmidtQ() click to toggle source

Returns the Q_1 matrix of Modified Gram Schmidt algorithm Q_1 has orthonormal columns

# File lib/extendmatrix.rb, line 877
def gram_schmidtQ
  gram_schmidt[0]
end
gram_schmidtR() click to toggle source

Returns the R_1 upper triangular matrix of Modified Gram Schmidt algorithm

# File lib/extendmatrix.rb, line 884
def gram_schmidtR
  gram_schmidt[1]
end
hessenbergQ() click to toggle source

Returns the orthogonal matrix Q of Hessenberg QR factorization Q = G_1 G_(n-1) where G_j is the Givens rotation G_j = G(j, j+1, omega_j)

# File lib/extendmatrix.rb, line 967
def hessenbergQ
  Hessenberg.QR(self)[0]
end
hessenbergR() click to toggle source

Returns the upper triunghiular matrix R of a Hessenberg QR factorization

# File lib/extendmatrix.rb, line 974
def hessenbergR
  Hessenberg.QR(self)[1]
end
hessenberg_form_H() click to toggle source

Return an upper Hessenberg matrix obtained with Householder reduction to Hessenberg Form algorithm

# File lib/extendmatrix.rb, line 981
def hessenberg_form_H
  Householder.toHessenberg(self)[0]
end
houseQ() click to toggle source

Returns the orthogonal matrix Q of Householder QR factorization where Q = H_1 * H_2 * H_3 * … * H_n,

# File lib/extendmatrix.rb, line 834
def houseQ
  h = Householder.QR(self)
  q = h[0]
  (1...h.size).each{|i| q *= h[i]}
  q
end
houseR() click to toggle source

Returns the matrix R of Householder QR factorization R = H_n * H_n-1 * … * H_1 * A is an upper triangular matrix

# File lib/extendmatrix.rb, line 845
def houseR
  h = Householder.QR(self)
  r = self.clone
  h.size.times{|i| r = h[i] * r}
  r
end
norm(p = 2) click to toggle source
# File lib/extendmatrix.rb, line 629
def norm(p = 2)
  Vector::Norm.sqnorm(self, p) ** (1.quo(p))
end
normF()
Alias for: norm_frobenius
norm_frobenius() click to toggle source
# File lib/extendmatrix.rb, line 633
def norm_frobenius
  norm
end
Also aliased as: normF
realSchur(eps = 1.0e-10, steps = 100) click to toggle source

The real Schur decomposition. The eigenvalues are aproximated in diagonal elements of the real Schur decomposition matrix

# File lib/extendmatrix.rb, line 989
def realSchur(eps = 1.0e-10, steps = 100)
  h = self.hessenberg_form_H
  h1 = Matrix[]
  i = 0
  loop do
    h1 = h.hessenbergR * h.hessenbergQ
    break if Matrix.diag_in_delta?(h1, h, eps) or steps <= 0
    h = h1.clone
    steps -= 1
    i += 1
  end
  h1
end
row!(i) { |e| ... } click to toggle source

Returns row vector number “i” like Matrix.row as a Vector. When the block is given, the elements of row “i” are modified

# File lib/extendmatrix.rb, line 570
def row!(i)
  if block_given?
    @rows[i].collect! {|e| yield e }
  else
    Vector.elements(@rows[i], false)
  end
end
Also aliased as: row_collect!
row2matrix(r) click to toggle source

Returns row(s) of matrix as a Matrix

# File lib/extendmatrix.rb, line 649
def row2matrix(r)
  if r.is_a? Range
    a=r.map {|v| v<row_size ? row(v).to_a : nil }.find_all {|v| !v.nil?}
  else
    a = row(r).to_a
  end
  if r.is_a?(Range) and r.entries.size > 1
    return Matrix[*a]
  else
    return Matrix[a]
  end
end
row=(args) click to toggle source

Set a certain row with the values of a Vector m = Matrix.new(3, 3){|i, j| i * 3 + j + 1} m.row= 0, Vector[0, 0], 1..2 m => 1 0 0

4 5 6
7 8 9
# File lib/extendmatrix.rb, line 624
def row=(args)
  i, val, range = MMatrix.id_vect_range(args, column_size)
  row!(i)[range] = val
end
row_collect(i, &block) click to toggle source

Returns an array with the elements collected from the row “i”. When a block is given, the elements of that vector are iterated.

# File lib/extendmatrix.rb, line 562
def row_collect(i, &block)
  f = MMatrix.default_block(block)
  @rows[i].collect {|e| f.call(e)}
end
row_collect!(i)
Alias for: row!
row_sum() click to toggle source

Returns the marginal of rows

# File lib/extendmatrix.rb, line 679
def row_sum
  (0...row_size).collect {|i|
    row(i).sum
  }
end
total_sum() click to toggle source

Calculate sum of cells

# File lib/extendmatrix.rb, line 691
def total_sum
  row_sum.inject(&:+)
end

SPSS like methods

↑ top

Public Instance Methods

diagonal() click to toggle source

Diagonal of matrix. Returns an array with as many elements as the minimum of the number of rows and the number of columns in the argument

# File lib/extendmatrix.rb, line 526
def diagonal
  min = (row_size<column_size) ? row_size : column_size
  min.times.map {|i| self[i,i]}
end
e_mult(other) click to toggle source

Element wise multiplication

# File lib/extendmatrix.rb, line 479
def e_mult(other)
  elementwise_operation(:*,other)
end
e_quo(other) click to toggle source

Element wise multiplication

# File lib/extendmatrix.rb, line 484
def e_quo(other)
  elementwise_operation(:quo,other)
end
eigen() click to toggle source

Returns eigenvalues and eigenvectors of a matrix on a Hash like CALL EIGEN on SPSS.

  • :eigenvectors: contains the eigenvectors as columns of a new Matrix, ordered in descendent order

  • :eigenvalues: contains an array with the eigenvalues, ordered in descendent order

# File lib/extendmatrix.rb, line 505
def eigen
  ep=eigenpairs
  { :eigenvalues=>ep.map {|a| a[0]} ,
    :eigenvectors=>Matrix.columns(ep.map{|a| a[1]})
  }
end
eigenpairs() click to toggle source
# File lib/extendmatrix.rb, line 493
def eigenpairs
  eigval, eigvec= eigenvaluesJacobi, cJacobiV
  eigenpairs=eigval.size.times.map {|i|
    [eigval[i], eigvec.column(i)]
  }
  eigenpairs=eigenpairs.sort{|a,b| a[0]<=>b[0]}.reverse
end
elementwise_operation(op,other) click to toggle source

Element wise operation

# File lib/extendmatrix.rb, line 472
def elementwise_operation(op,other)
  Matrix.build(row_size,column_size) do |row, column|
    self[row,column].send(op,other[row,column])
  end
end
ln() click to toggle source

Returns a new matrix with the natural logarithm of initial values

# File lib/extendmatrix.rb, line 512
def ln
  map {|v| Math::log(v)}
end
mssq() click to toggle source

Matrix sum of squares

# File lib/extendmatrix.rb, line 489
def mssq
  @rows.inject(0){|ac,row| ac+(row.inject(0) {|acr,i| acr+(i**2)})}
end
sqrt() click to toggle source

Returns a new matrix, with the square root of initial values

# File lib/extendmatrix.rb, line 516
def sqrt
  map {|v| Math::sqrt(v)}
end
sscp() click to toggle source

Sum of squares and cross-products. Equal to m.t * m

# File lib/extendmatrix.rb, line 520
def sscp
  transpose*self
end