The Matrix
class represents a mathematical matrix. It provides
methods for creating matrices, operating on them arithmetically and
algebraically, and determining their mathematical properties (trace, rank,
inverse, determinant).
Method Catalogue
To create a matrix:
-
::rows(rows, copy = true)
-
::build(row_count, #column_count, &block)
-
::scalar(n, value)
-
Matrix.I(n)
To access Matrix elements/columns/rows/submatrices/properties:
-
[](i, j)
-
row_count (row_size)
-
column_count (column_size)
Properties of a matrix:
Matrix arithmetic:
Matrix functions:
Matrix decompositions:
Complex arithmetic:
-
conj
-
conjugate
-
imag
-
imaginary
-
real
-
rect
-
rectangular
Conversion to other data types:
String representations:
- #
- I
- #
- B
- C
- D
- E
- F
- H
- I
- L
- M
- N
- O
- P
- R
-
- rank,
- rank_e,
- real,
- real?,
- rect,
- rectangular,
- regular?,
- round,
- row,
- row_count,
- row_size,
- row_vector,
- row_vectors,
- rows
- S
- T
- U
- Z
SELECTORS | = | {all: true, diagonal: true, off_diagonal: true, lower: true, strict_lower: true, strict_upper: true, upper: true}.freeze |
[R] | column_count | Returns the number of columns. |
[R] | column_size | Returns the number of columns. |
[R] | rows | instance creations |
Creates a matrix where each argument is a row.
Matrix[ [25, 93], [-1, 66] ]
=> 25 93
-1 66
Creates a matrix of size row_count
x
column_count
. It fills the values by calling the given block,
passing the current row and column. Returns an enumerator if no block is
given.
m = Matrix.build(2, 4) {|row, col| col - row }
=> Matrix[[0, 1, 2, 3], [-1, 0, 1, 2]]
m = Matrix.build(3) { rand }
=> a 3x3 matrix with random elements
# File lib/matrix.rb, line 184 def Matrix.build(row_count, column_count = row_count) row_count = CoercionHelper.coerce_to_int(row_count) column_count = CoercionHelper.coerce_to_int(column_count) raise ArgumentError if row_count < 0 || column_count < 0 return to_enum :build, row_count, column_count unless block_given? rows = Array.new(row_count) do |i| Array.new(column_count) do |j| yield i, j end end new rows, column_count end
Creates a single-column matrix where the values of that column are as given
in column
.
Matrix.column_vector([4,5,6])
=> 4
5
6
Creates a matrix using columns
as an array of column vectors.
Matrix.columns([[25, 93], [-1, 66]])
=> 25 -1
93 66
Creates a matrix where the diagonal elements are composed of
values
.
Matrix.diagonal(9, 5, -3)
=> 9 0 0
0 5 0
0 0 -3
Creates a empty matrix of row_count
x
column_count
. At least one of row_count
or
column_count
must be 0.
m = Matrix.empty(2, 0)
m == Matrix[ [], [] ]
=> true
n = Matrix.empty(0, 3)
n == Matrix.columns([ [], [], [] ])
=> true
m * n
=> Matrix[[0, 0, 0], [0, 0, 0]]
Creates an n
by n
identity matrix.
Matrix.identity(2)
=> 1 0
0 1
# File lib/matrix.rb, line 297 def initialize(rows, column_count = rows[0].size) # No checking is done at this point. rows must be an Array of Arrays. # column_count must be the size of the first row, if there is one, # otherwise it *must* be specified and can be any integer >= 0 @rows = rows @column_count = column_count end
Creates a single-row matrix where the values of that row are as given in
row
.
Matrix.row_vector([4,5,6])
=> 4 5 6
Creates a matrix where rows
is an array of arrays, each of
which is a row of the matrix. If the optional argument copy
is false, use the given arrays as the internal structure of the matrix
without copying.
Matrix.rows([[25, 93], [-1, 66]])
=> 25 93
-1 66
# File lib/matrix.rb, line 151 def Matrix.rows(rows, copy = true) rows = convert_to_array(rows) rows.map! do |row| convert_to_array(row, copy) end size = (rows[0] || []).size rows.each do |row| raise ErrDimensionMismatch, "row size differs (#{row.size} should be #{size})" unless row.size == size end new rows, size end
Creates an n
by n
diagonal matrix where each
diagonal element is value
.
Matrix.scalar(2, 5)
=> 5 0
0 5
Creates a zero matrix.
Matrix.zero(2)
=> 0 0
0 0
Matrix multiplication.
Matrix[[2,4], [6,8]] * Matrix.identity(2)
=> 2 4
6 8
# File lib/matrix.rb, line 805 def *(m) # m is matrix or vector or number case(m) when Numeric rows = @rows.collect {|row| row.collect {|e| e * m } } return new_matrix rows, column_count when Vector m = self.class.column_vector(m) r = self * m return r.column(0) when Matrix Matrix.Raise ErrDimensionMismatch if column_count != m.row_count rows = Array.new(row_count) {|i| Array.new(m.column_count) {|j| (0 ... column_count).inject(0) do |vij, k| vij + self[i, k] * m[k, j] end } } return new_matrix rows, m.column_count else return apply_through_coercion(m, __method__) end end
Matrix exponentiation. Equivalent to multiplying the matrix by itself N times. Non integer exponents will be handled by diagonalizing the matrix.
Matrix[[7,6], [3,9]] ** 2
=> 67 96
48 99
# File lib/matrix.rb, line 972 def ** (other) case other when Integer x = self if other <= 0 x = self.inverse return self.class.identity(self.column_count) if other == 0 other = -other end z = nil loop do z = z ? z * x : x if other[0] == 1 return z if (other >>= 1).zero? x *= x end when Numeric v, d, v_inv = eigensystem v * self.class.diagonal(*d.each(:diagonal).map{|e| e ** other}) * v_inv else Matrix.Raise ErrOperationNotDefined, "**", self.class, other.class end end
Matrix addition.
Matrix.scalar(2,5) + Matrix[[1,0], [-4,7]]
=> 6 0
-4 12
# File lib/matrix.rb, line 838 def +(m) case m when Numeric Matrix.Raise ErrOperationNotDefined, "+", self.class, m.class when Vector m = self.class.column_vector(m) when Matrix else return apply_through_coercion(m, __method__) end Matrix.Raise ErrDimensionMismatch unless row_count == m.row_count and column_count == m.column_count rows = Array.new(row_count) {|i| Array.new(column_count) {|j| self[i, j] + m[i, j] } } new_matrix rows, column_count end
Matrix subtraction.
Matrix[[1,5], [4,2]] - Matrix[[9,3], [-4,1]]
=> -8 2
8 1
# File lib/matrix.rb, line 865 def -(m) case m when Numeric Matrix.Raise ErrOperationNotDefined, "-", self.class, m.class when Vector m = self.class.column_vector(m) when Matrix else return apply_through_coercion(m, __method__) end Matrix.Raise ErrDimensionMismatch unless row_count == m.row_count and column_count == m.column_count rows = Array.new(row_count) {|i| Array.new(column_count) {|j| self[i, j] - m[i, j] } } new_matrix rows, column_count end
Matrix division (multiplication by the inverse).
Matrix[[7,6], [3,9]] / Matrix[[2,9], [3,1]]
=> -7 1
-3 -6
Returns true
if and only if the two matrices contain equal
elements.
Returns element (i
,j
) of the matrix. That is:
row i
, column j
.
Returns a clone of the matrix, so that the contents of each do not reference identical objects. There should be no good reason to do this since Matrices are immutable.
The coerce method provides support for Ruby type coercion. This coercion mechanism is used by Ruby to handle mixed-type numeric operations: it is intended to find a compatible common type between the two operands of the operator. See also Numeric#coerce.
Returns a matrix that is the result of iteration of the given block over all elements of the matrix.
Matrix[ [1,2], [3,4] ].collect { |e| e**2 }
=> 1 4
9 16
Returns column vector number j
of the matrix as a Vector (starting at 0 like an array). When a block
is given, the elements of that vector are iterated.
# File lib/matrix.rb, line 358 def column(j) # :yield: e if block_given? return self if j >= column_count || j < -column_count row_count.times do |i| yield @rows[i][j] end self else return nil if j >= column_count || j < -column_count col = Array.new(row_count) {|i| @rows[i][j] } Vector.elements(col, false) end end
Returns the conjugate of the matrix.
Matrix[[Complex(1,2), Complex(0,1), 0], [1, 2, 3]]
=> 1+2i i 0
1 2 3
Matrix[[Complex(1,2), Complex(0,1), 0], [1, 2, 3]].conjugate
=> 1-2i -i 0
1 2 3
Returns the determinant of the matrix.
Beware that using Float values can yield erroneous results because of their lack of precision. Consider using exact types like Rational or BigDecimal instead.
Matrix[[7,6], [3,9]].determinant
=> 45
# File lib/matrix.rb, line 1009 def determinant Matrix.Raise ErrDimensionMismatch unless square? m = @rows case row_count # Up to 4x4, give result using Laplacian expansion by minors. # This will typically be faster, as well as giving good results # in case of Floats when 0 +1 when 1 + m[0][0] when 2 + m[0][0] * m[1][1] - m[0][1] * m[1][0] when 3 m0, m1, m2 = m + m0[0] * m1[1] * m2[2] - m0[0] * m1[2] * m2[1] - m0[1] * m1[0] * m2[2] + m0[1] * m1[2] * m2[0] + m0[2] * m1[0] * m2[1] - m0[2] * m1[1] * m2[0] when 4 m0, m1, m2, m3 = m + m0[0] * m1[1] * m2[2] * m3[3] - m0[0] * m1[1] * m2[3] * m3[2] - m0[0] * m1[2] * m2[1] * m3[3] + m0[0] * m1[2] * m2[3] * m3[1] + m0[0] * m1[3] * m2[1] * m3[2] - m0[0] * m1[3] * m2[2] * m3[1] - m0[1] * m1[0] * m2[2] * m3[3] + m0[1] * m1[0] * m2[3] * m3[2] + m0[1] * m1[2] * m2[0] * m3[3] - m0[1] * m1[2] * m2[3] * m3[0] - m0[1] * m1[3] * m2[0] * m3[2] + m0[1] * m1[3] * m2[2] * m3[0] + m0[2] * m1[0] * m2[1] * m3[3] - m0[2] * m1[0] * m2[3] * m3[1] - m0[2] * m1[1] * m2[0] * m3[3] + m0[2] * m1[1] * m2[3] * m3[0] + m0[2] * m1[3] * m2[0] * m3[1] - m0[2] * m1[3] * m2[1] * m3[0] - m0[3] * m1[0] * m2[1] * m3[2] + m0[3] * m1[0] * m2[2] * m3[1] + m0[3] * m1[1] * m2[0] * m3[2] - m0[3] * m1[1] * m2[2] * m3[0] - m0[3] * m1[2] * m2[0] * m3[1] + m0[3] * m1[2] * m2[1] * m3[0] else # For bigger matrices, use an efficient and general algorithm. # Currently, we use the Gauss-Bareiss algorithm determinant_bareiss end end
Returns true
is this is a diagonal matrix. Raises an error if
matrix is not square.
Yields all elements of the matrix, starting with those of the first row, or returns an Enumerator is no block given. Elements can be restricted by passing an argument:
-
:all (default): yields all elements
-
:diagonal: yields only elements on the diagonal
-
:off_diagonal: yields all elements except on the diagonal
-
:lower: yields only elements on or below the diagonal
-
:strict_lower: yields only elements below the diagonal
-
:strict_upper: yields only elements above the diagonal
-
:upper: yields only elements on or above the diagonal
Matrix[ [1,2], [3,4] ].each { |e| puts e }
# => prints the numbers 1 to 4
Matrix[ [1,2], [3,4] ].each(:strict_lower).to_a # => [3]
# File lib/matrix.rb, line 404 def each(which = :all) # :yield: e return to_enum :each, which unless block_given? last = column_count - 1 case which when :all block = Proc.new @rows.each do |row| row.each(&block) end when :diagonal @rows.each_with_index do |row, row_index| yield row.fetch(row_index){return self} end when :off_diagonal @rows.each_with_index do |row, row_index| column_count.times do |col_index| yield row[col_index] unless row_index == col_index end end when :lower @rows.each_with_index do |row, row_index| 0.upto([row_index, last].min) do |col_index| yield row[col_index] end end when :strict_lower @rows.each_with_index do |row, row_index| [row_index, column_count].min.times do |col_index| yield row[col_index] end end when :strict_upper @rows.each_with_index do |row, row_index| (row_index+1).upto(last) do |col_index| yield row[col_index] end end when :upper @rows.each_with_index do |row, row_index| row_index.upto(last) do |col_index| yield row[col_index] end end else raise ArgumentError, "expected #{which.inspect} to be one of :all, :diagonal, :off_diagonal, :lower, :strict_lower, :strict_upper or :upper" end self end
Same as each, but the row index and column index in addition to the element
Matrix[ [1,2], [3,4] ].each_with_index do |e, row, col|
puts "#{e} at #{row}, #{col}"
end
# => Prints:
# 1 at 0, 0
# 2 at 0, 1
# 3 at 1, 0
# 4 at 1, 1
# File lib/matrix.rb, line 465 def each_with_index(which = :all) # :yield: e, row, column return to_enum :each_with_index, which unless block_given? last = column_count - 1 case which when :all @rows.each_with_index do |row, row_index| row.each_with_index do |e, col_index| yield e, row_index, col_index end end when :diagonal @rows.each_with_index do |row, row_index| yield row.fetch(row_index){return self}, row_index, row_index end when :off_diagonal @rows.each_with_index do |row, row_index| column_count.times do |col_index| yield row[col_index], row_index, col_index unless row_index == col_index end end when :lower @rows.each_with_index do |row, row_index| 0.upto([row_index, last].min) do |col_index| yield row[col_index], row_index, col_index end end when :strict_lower @rows.each_with_index do |row, row_index| [row_index, column_count].min.times do |col_index| yield row[col_index], row_index, col_index end end when :strict_upper @rows.each_with_index do |row, row_index| (row_index+1).upto(last) do |col_index| yield row[col_index], row_index, col_index end end when :upper @rows.each_with_index do |row, row_index| row_index.upto(last) do |col_index| yield row[col_index], row_index, col_index end end else raise ArgumentError, "expected #{which.inspect} to be one of :all, :diagonal, :off_diagonal, :lower, :strict_lower, :strict_upper or :upper" end self end
Returns the Eigensystem of the matrix; see
EigenvalueDecomposition
.
m = Matrix[[1, 2], [3, 4]]
v, d, v_inv = m.eigensystem
d.diagonal? # => true
v.inv == v_inv # => true
(v * d * v_inv).round(5) == m # => true
Returns true
if this is an empty matrix, i.e. if the number of
rows or the number of columns is 0.
Returns a hash-code for the matrix.
Returns true
is this is an hermitian matrix. Raises an error
if matrix is not square.
Returns the imaginary part of the matrix.
Matrix[[Complex(1,2), Complex(0,1), 0], [1, 2, 3]]
=> 1+2i i 0
1 2 3
Matrix[[Complex(1,2), Complex(0,1), 0], [1, 2, 3]].imaginary
=> 2i i 0
0 0 0
The index method is specialized to return the index as [row, column] It
also accepts an optional selector
argument, see each for details.
Matrix[ [1,2], [3,4] ].index(&:even?) # => [0, 1]
Matrix[ [1,1], [1,1] ].index(1, :strict_lower) # => [1, 0]
# File lib/matrix.rb, line 528 def index(*args) raise ArgumentError, "wrong number of arguments(#{args.size} for 0-2)" if args.size > 2 which = (args.size == 2 || SELECTORS.include?(args.last)) ? args.pop : :all return to_enum :find_index, which, *args unless block_given? || args.size == 1 if args.size == 1 value = args.first each_with_index(which) do |e, row_index, col_index| return row_index, col_index if e == value end else each_with_index(which) do |e, row_index, col_index| return row_index, col_index if yield e end end nil end
Overrides Object#inspect
Returns the inverse of the matrix.
Matrix[[-1, -1], [0, -1]].inverse
=> -1 1
0 -1
Returns true
is this is a lower triangular matrix.
Returns the LUP decomposition of the matrix; see
LUPDecomposition
.
a = Matrix[[1, 2], [3, 4]]
l, u, p = a.lup
l.lower_triangular? # => true
u.upper_triangular? # => true
p.permutation? # => true
l * u == p * a # => true
a.lup.solve([2, 5]) # => Vector[(1/1), (1/2)]
Returns a section of the matrix. The parameters are either:
-
start_row, nrows, start_col, ncols; OR
-
row_range, col_range
Matrix.diagonal(9, 5, -3).minor(0..1, 0..2)
=> 9 0 0
0 5 0
Like Array#[], negative indices count backward from the end of the row or column (-1 is the last element). Returns nil if the starting row or column is greater than #row_count or #column_count respectively.
# File lib/matrix.rb, line 558 def minor(*param) case param.size when 2 row_range, col_range = param from_row = row_range.first from_row += row_count if from_row < 0 to_row = row_range.end to_row += row_count if to_row < 0 to_row += 1 unless row_range.exclude_end? size_row = to_row - from_row from_col = col_range.first from_col += column_count if from_col < 0 to_col = col_range.end to_col += column_count if to_col < 0 to_col += 1 unless col_range.exclude_end? size_col = to_col - from_col when 4 from_row, size_row, from_col, size_col = param return nil if size_row < 0 || size_col < 0 from_row += row_count if from_row < 0 from_col += column_count if from_col < 0 else raise ArgumentError, param.inspect end return nil if from_row > row_count || from_col > column_count || from_row < 0 || from_col < 0 rows = @rows[from_row, size_row].collect{|row| row[from_col, size_col] } new_matrix rows, [column_count - from_col, size_col].min end
Returns true
is this is a normal matrix. Raises an error if
matrix is not square.
# File lib/matrix.rb, line 634 def normal? Matrix.Raise ErrDimensionMismatch unless square? rows.each_with_index do |row_i, i| rows.each_with_index do |row_j, j| s = 0 rows.each_with_index do |row_k, k| s += row_i[k] * row_j[k].conj - row_k[i].conj * row_k[j] end return false unless s == 0 end end true end
Returns true
is this is an orthogonal matrix Raises an error
if matrix is not square.
Returns true
is this is a permutation matrix Raises an error
if matrix is not square.
# File lib/matrix.rb, line 670 def permutation? Matrix.Raise ErrDimensionMismatch unless square? cols = Array.new(column_count) rows.each_with_index do |row, i| found = false row.each_with_index do |e, j| if e == 1 return false if found || cols[j] found = cols[j] = true elsif e != 0 return false end end return false unless found end true end
Returns the rank of the matrix. Beware that using Float values can yield erroneous results because of their lack of precision. Consider using exact types like Rational or BigDecimal instead.
Matrix[[7,6], [3,9]].rank
=> 2
# File lib/matrix.rb, line 1106 def rank # We currently use Bareiss' multistep integer-preserving gaussian elimination # (see comments on determinant) a = to_a last_column = column_count - 1 last_row = row_count - 1 pivot_row = 0 previous_pivot = 1 0.upto(last_column) do |k| switch_row = (pivot_row .. last_row).find {|row| a[row][k] != 0 } if switch_row a[switch_row], a[pivot_row] = a[pivot_row], a[switch_row] unless pivot_row == switch_row pivot = a[pivot_row][k] (pivot_row+1).upto(last_row) do |i| ai = a[i] (k+1).upto(last_column) do |j| ai[j] = (pivot * ai[j] - ai[k] * a[pivot_row][j]) / previous_pivot end end pivot_row += 1 previous_pivot = pivot end end pivot_row end
Returns the real part of the matrix.
Matrix[[Complex(1,2), Complex(0,1), 0], [1, 2, 3]]
=> 1+2i i 0
1 2 3
Matrix[[Complex(1,2), Complex(0,1), 0], [1, 2, 3]].real
=> 1 0 0
1 2 3
Returns true
if all entries of the matrix are real.
Returns an array containing matrices corresponding to the real and imaginary parts of the matrix
m.rect == [m.real, m.imag] # ==> true for all matrices m
Returns true
if this is a regular (i.e. non-singular) matrix.
Returns a matrix with entries rounded to the given precision (see Float#round)
Returns row vector number i
of the matrix as a Vector (starting at 0 like an array). When a block
is given, the elements of that vector are iterated.
Returns true
is this is a singular matrix.
Returns true
is this is a square matrix.
Returns true
is this is a symmetric matrix. Raises an error if
matrix is not square.
Returns an array of arrays that describe the rows of the matrix.
Overrides Object#to_s
Returns the trace (sum of diagonal elements) of the matrix.
Matrix[[7,6], [3,9]].trace
=> 16
Returns the transpose of the matrix.
Matrix[[1,2], [3,4], [5,6]]
=> 1 2
3 4
5 6
Matrix[[1,2], [3,4], [5,6]].transpose
=> 1 3 5
2 4 6
Returns true
is this is a unitary matrix Raises an error if
matrix is not square.
Returns true
is this is an upper triangular matrix.