Module: Matrix::Jacobi

Defined in:
lib/extendmatrix.rb

Class Method Summary collapse

Class Method Details

.J(p, q, c, s, n) ⇒ Object

Returns the Jacobi rotation matrix



1022
1023
1024
1025
1026
1027
 
# File 'lib/extendmatrix.rb', line 1022

def self.J(p, q, c, s, n)
  j = Matrix.I(n)
  j[p,p] = c; j[p, q] = s
  j[q,p] = -s; j[q, q] = c
  j
end

.max(a) ⇒ Object

Returns the index pair (p, q) with 1<= p < q <= n and A[p, q] is the maximum in absolute value



983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
 
# File 'lib/extendmatrix.rb', line 983

def self.max(a)
        n = a.row_size
        max = 0
        p = 0
        q = 0
        n.times{|i|
((i+1)...n).each{|j|
  val = a[i, j].abs
  if val > max
    max = val
    p = i
    q = j
  end	}}
  return p, q
end

.off(a) ⇒ Object

Returns the nurm of the off-diagonal element



973
974
975
976
977
978
 
# File 'lib/extendmatrix.rb', line 973

def self.off(a)
  n = a.row_size
  sum = 0
  n.times{|i| n.times{|j| sum += a[i, j]**2 if j != i}}
  Math.sqrt(sum)
end

.sym_schur2(a, p, q) ⇒ Object

Compute the cosine-sine pair (c, s) for the element A[p, q]



1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
 
# File 'lib/extendmatrix.rb', line 1002

def self.sym_schur2(a, p, q)
  if a[p, q] != 0
    tau = Float(a[q, q] - a[p, p])/(2 * a[p, q])
    if tau >= 0
      t = 1./(tau + Math.sqrt(1 + tau ** 2))
    else
      t = -1./(-tau + Math.sqrt(1 + tau ** 2))
    end
    c = 1./Math.sqrt(1 + t ** 2)
    s = t * c
  else
    c = 1
    s = 0
  end
  return c, s
end