Extended maintenance of Ruby versions 1.8.7 and 1.9.2 ended on July 31, 2014. Read more

In Files

  • bigdecimal/lib/bigdecimal/ludcmp.rb

Class/Module Index [+]

Quicksearch

LUSolve

Solves a*x = b for x, using LU decomposition.

Public Instance Methods

ludecomp(a,n,zero=0,one=1) click to toggle source

Performs LU decomposition of the n by n matrix a.

 
               # File bigdecimal/lib/bigdecimal/ludcmp.rb, line 10
  def ludecomp(a,n,zero=0,one=1)
    prec = BigDecimal.limit(nil)
    ps     = []
    scales = []
    for i in 0...n do  # pick up largest(abs. val.) element in each row.
      ps <<= i
      nrmrow  = zero
      ixn = i*n
      for j in 0...n do
        biggst = a[ixn+j].abs
        nrmrow = biggst if biggst>nrmrow
      end
      if nrmrow>zero then
        scales <<= one.div(nrmrow,prec)
      else
        raise "Singular matrix"
      end
    end
    n1          = n - 1
    for k in 0...n1 do # Gaussian elimination with partial pivoting.
      biggst  = zero;
      for i in k...n do
        size = a[ps[i]*n+k].abs*scales[ps[i]]
        if size>biggst then
          biggst = size
          pividx  = i
        end
      end
      raise "Singular matrix" if biggst<=zero
      if pividx!=k then
        j = ps[k]
        ps[k] = ps[pividx]
        ps[pividx] = j
      end
      pivot   = a[ps[k]*n+k]
      for i in (k+1)...n do
        psin = ps[i]*n
        a[psin+k] = mult = a[psin+k].div(pivot,prec)
        if mult!=zero then
          pskn = ps[k]*n
          for j in (k+1)...n do
            a[psin+j] -= mult.mult(a[pskn+j],prec)
          end
        end
      end
    end
    raise "Singular matrix" if a[ps[n1]*n+n1] == zero
    ps
  end

  # Solves a*x = b for x, using LU decomposition.
  #
  # a is a matrix, b is a constant vector, x is the solution vector.
  #
  # ps is the pivot, a vector which indicates the permutation of rows performed
  # during LU decomposition.
  def lusolve(a,b,ps,zero=0.0)
    prec = BigDecimal.limit(nil)
    n = ps.size
    x = []
    for i in 0...n do
      dot = zero
      psin = ps[i]*n
      for j in 0...i do
        dot = a[psin+j].mult(x[j],prec) + dot
      end
      x <<= b[ps[i]] - dot
    end
    (n-1).downto(0) do |i|
      dot = zero
      psin = ps[i]*n
      for j in (i+1)...n do
        dot = a[psin+j].mult(x[j],prec) + dot
      end
      x[i]  = (x[i]-dot).div(a[psin+i],prec)
    end
    x
  end
end

            
lusolve(a,b,ps,zero=0.0) click to toggle source

Solves a*x = b for x, using LU decomposition.

a is a matrix, b is a constant vector, x is the solution vector.

ps is the pivot, a vector which indicates the permutation of rows performed during LU decomposition.

 
               # File bigdecimal/lib/bigdecimal/ludcmp.rb, line 66
  def lusolve(a,b,ps,zero=0.0)
    prec = BigDecimal.limit(nil)
    n = ps.size
    x = []
    for i in 0...n do
      dot = zero
      psin = ps[i]*n
      for j in 0...i do
        dot = a[psin+j].mult(x[j],prec) + dot
      end
      x <<= b[ps[i]] - dot
    end
    (n-1).downto(0) do |i|
      dot = zero
      psin = ps[i]*n
      for j in (i+1)...n do
        dot = a[psin+j].mult(x[j],prec) + dot
      end
      x[i]  = (x[i]-dot).div(a[psin+i],prec)
    end
    x
  end
end
            

Commenting is here to help enhance the documentation. For example, code samples, or clarification of the documentation.

If you have questions about Ruby or the documentation, please post to one of the Ruby mailing lists. You will get better, faster, help that way.

If you wish to post a correction of the docs, please do so, but also file bug report so that it can be corrected for the next release. Thank you.

If you want to help improve the Ruby documentation, please visit Documenting-ruby.org.

blog comments powered by Disqus