diff options
Diffstat (limited to 'ext/bigdecimal/lib')
| -rw-r--r-- | ext/bigdecimal/lib/bigdecimal.rb | 1 | ||||
| -rw-r--r-- | ext/bigdecimal/lib/bigdecimal/jacobian.rb | 90 | ||||
| -rw-r--r-- | ext/bigdecimal/lib/bigdecimal/ludcmp.rb | 89 | ||||
| -rw-r--r-- | ext/bigdecimal/lib/bigdecimal/math.rb | 232 | ||||
| -rw-r--r-- | ext/bigdecimal/lib/bigdecimal/newton.rb | 80 | ||||
| -rw-r--r-- | ext/bigdecimal/lib/bigdecimal/util.rb | 181 |
6 files changed, 0 insertions, 673 deletions
diff --git a/ext/bigdecimal/lib/bigdecimal.rb b/ext/bigdecimal/lib/bigdecimal.rb deleted file mode 100644 index 8fd2587c84..0000000000 --- a/ext/bigdecimal/lib/bigdecimal.rb +++ /dev/null @@ -1 +0,0 @@ -require 'bigdecimal.so' diff --git a/ext/bigdecimal/lib/bigdecimal/jacobian.rb b/ext/bigdecimal/lib/bigdecimal/jacobian.rb deleted file mode 100644 index 5e29304299..0000000000 --- a/ext/bigdecimal/lib/bigdecimal/jacobian.rb +++ /dev/null @@ -1,90 +0,0 @@ -# frozen_string_literal: false - -require 'bigdecimal' - -# require 'bigdecimal/jacobian' -# -# Provides methods to compute the Jacobian matrix of a set of equations at a -# point x. In the methods below: -# -# f is an Object which is used to compute the Jacobian matrix of the equations. -# It must provide the following methods: -# -# f.values(x):: returns the values of all functions at x -# -# f.zero:: returns 0.0 -# f.one:: returns 1.0 -# f.two:: returns 2.0 -# f.ten:: returns 10.0 -# -# f.eps:: returns the convergence criterion (epsilon value) used to determine whether two values are considered equal. If |a-b| < epsilon, the two values are considered equal. -# -# x is the point at which to compute the Jacobian. -# -# fx is f.values(x). -# -module Jacobian - module_function - - # Determines the equality of two numbers by comparing to zero, or using the epsilon value - def isEqual(a,b,zero=0.0,e=1.0e-8) - aa = a.abs - bb = b.abs - if aa == zero && bb == zero then - true - else - if ((a-b)/(aa+bb)).abs < e then - true - else - false - end - end - end - - - # Computes the derivative of f[i] at x[i]. - # fx is the value of f at x. - def dfdxi(f,fx,x,i) - nRetry = 0 - n = x.size - xSave = x[i] - ok = 0 - ratio = f.ten*f.ten*f.ten - dx = x[i].abs/ratio - dx = fx[i].abs/ratio if isEqual(dx,f.zero,f.zero,f.eps) - dx = f.one/f.ten if isEqual(dx,f.zero,f.zero,f.eps) - until ok>0 do - deriv = [] - nRetry += 1 - if nRetry > 100 - raise "Singular Jacobian matrix. No change at x[" + i.to_s + "]" - end - dx = dx*f.two - x[i] += dx - fxNew = f.values(x) - for j in 0...n do - if !isEqual(fxNew[j],fx[j],f.zero,f.eps) then - ok += 1 - deriv <<= (fxNew[j]-fx[j])/dx - else - deriv <<= f.zero - end - end - x[i] = xSave - end - deriv - end - - # Computes the Jacobian of f at x. fx is the value of f at x. - def jacobian(f,fx,x) - n = x.size - dfdx = Array.new(n*n) - for i in 0...n do - df = dfdxi(f,fx,x,i) - for j in 0...n do - dfdx[j*n+i] = df[j] - end - end - dfdx - end -end diff --git a/ext/bigdecimal/lib/bigdecimal/ludcmp.rb b/ext/bigdecimal/lib/bigdecimal/ludcmp.rb deleted file mode 100644 index dd265e482a..0000000000 --- a/ext/bigdecimal/lib/bigdecimal/ludcmp.rb +++ /dev/null @@ -1,89 +0,0 @@ -# frozen_string_literal: false -require 'bigdecimal' - -# -# Solves a*x = b for x, using LU decomposition. -# -module LUSolve - module_function - - # Performs LU decomposition of the n by n matrix a. - 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 diff --git a/ext/bigdecimal/lib/bigdecimal/math.rb b/ext/bigdecimal/lib/bigdecimal/math.rb deleted file mode 100644 index 0b9d0648bb..0000000000 --- a/ext/bigdecimal/lib/bigdecimal/math.rb +++ /dev/null @@ -1,232 +0,0 @@ -# frozen_string_literal: false -require 'bigdecimal' - -# -#-- -# Contents: -# sqrt(x, prec) -# sin (x, prec) -# cos (x, prec) -# atan(x, prec) Note: |x|<1, x=0.9999 may not converge. -# PI (prec) -# E (prec) == exp(1.0,prec) -# -# where: -# x ... BigDecimal number to be computed. -# |x| must be small enough to get convergence. -# prec ... Number of digits to be obtained. -#++ -# -# Provides mathematical functions. -# -# Example: -# -# require "bigdecimal/math" -# -# include BigMath -# -# a = BigDecimal((PI(100)/2).to_s) -# puts sin(a,100) # => 0.99999999999999999999......e0 -# -module BigMath - module_function - - # call-seq: - # sqrt(decimal, numeric) -> BigDecimal - # - # Computes the square root of +decimal+ to the specified number of digits of - # precision, +numeric+. - # - # BigMath.sqrt(BigDecimal('2'), 16).to_s - # #=> "0.1414213562373095048801688724e1" - # - def sqrt(x, prec) - x.sqrt(prec) - end - - # call-seq: - # sin(decimal, numeric) -> BigDecimal - # - # Computes the sine of +decimal+ to the specified number of digits of - # precision, +numeric+. - # - # If +decimal+ is Infinity or NaN, returns NaN. - # - # BigMath.sin(BigMath.PI(5)/4, 5).to_s - # #=> "0.70710678118654752440082036563292800375e0" - # - def sin(x, prec) - raise ArgumentError, "Zero or negative precision for sin" if prec <= 0 - return BigDecimal("NaN") if x.infinite? || x.nan? - n = prec + BigDecimal.double_fig - one = BigDecimal("1") - two = BigDecimal("2") - x = -x if neg = x < 0 - if x > (twopi = two * BigMath.PI(prec)) - if x > 30 - x %= twopi - else - x -= twopi while x > twopi - end - end - x1 = x - x2 = x.mult(x,n) - sign = 1 - y = x - d = y - i = one - z = one - while d.nonzero? && ((m = n - (y.exponent - d.exponent).abs) > 0) - m = BigDecimal.double_fig if m < BigDecimal.double_fig - sign = -sign - x1 = x2.mult(x1,n) - i += two - z *= (i-one) * i - d = sign * x1.div(z,m) - y += d - end - neg ? -y : y - end - - # call-seq: - # cos(decimal, numeric) -> BigDecimal - # - # Computes the cosine of +decimal+ to the specified number of digits of - # precision, +numeric+. - # - # If +decimal+ is Infinity or NaN, returns NaN. - # - # BigMath.cos(BigMath.PI(4), 16).to_s - # #=> "-0.999999999999999999999999999999856613163740061349e0" - # - def cos(x, prec) - raise ArgumentError, "Zero or negative precision for cos" if prec <= 0 - return BigDecimal("NaN") if x.infinite? || x.nan? - n = prec + BigDecimal.double_fig - one = BigDecimal("1") - two = BigDecimal("2") - x = -x if x < 0 - if x > (twopi = two * BigMath.PI(prec)) - if x > 30 - x %= twopi - else - x -= twopi while x > twopi - end - end - x1 = one - x2 = x.mult(x,n) - sign = 1 - y = one - d = y - i = BigDecimal("0") - z = one - while d.nonzero? && ((m = n - (y.exponent - d.exponent).abs) > 0) - m = BigDecimal.double_fig if m < BigDecimal.double_fig - sign = -sign - x1 = x2.mult(x1,n) - i += two - z *= (i-one) * i - d = sign * x1.div(z,m) - y += d - end - y - end - - # call-seq: - # atan(decimal, numeric) -> BigDecimal - # - # Computes the arctangent of +decimal+ to the specified number of digits of - # precision, +numeric+. - # - # If +decimal+ is NaN, returns NaN. - # - # BigMath.atan(BigDecimal('-1'), 16).to_s - # #=> "-0.785398163397448309615660845819878471907514682065e0" - # - def atan(x, prec) - raise ArgumentError, "Zero or negative precision for atan" if prec <= 0 - return BigDecimal("NaN") if x.nan? - pi = PI(prec) - x = -x if neg = x < 0 - return pi.div(neg ? -2 : 2, prec) if x.infinite? - return pi / (neg ? -4 : 4) if x.round(prec) == 1 - x = BigDecimal("1").div(x, prec) if inv = x > 1 - x = (-1 + sqrt(1 + x**2, prec))/x if dbl = x > 0.5 - n = prec + BigDecimal.double_fig - y = x - d = y - t = x - r = BigDecimal("3") - x2 = x.mult(x,n) - while d.nonzero? && ((m = n - (y.exponent - d.exponent).abs) > 0) - m = BigDecimal.double_fig if m < BigDecimal.double_fig - t = -t.mult(x2,n) - d = t.div(r,m) - y += d - r += 2 - end - y *= 2 if dbl - y = pi / 2 - y if inv - y = -y if neg - y - end - - # call-seq: - # PI(numeric) -> BigDecimal - # - # Computes the value of pi to the specified number of digits of precision, - # +numeric+. - # - # BigMath.PI(10).to_s - # #=> "0.3141592653589793238462643388813853786957412e1" - # - def PI(prec) - raise ArgumentError, "Zero or negative precision for PI" if prec <= 0 - n = prec + BigDecimal.double_fig - zero = BigDecimal("0") - one = BigDecimal("1") - two = BigDecimal("2") - - m25 = BigDecimal("-0.04") - m57121 = BigDecimal("-57121") - - pi = zero - - d = one - k = one - t = BigDecimal("-80") - while d.nonzero? && ((m = n - (pi.exponent - d.exponent).abs) > 0) - m = BigDecimal.double_fig if m < BigDecimal.double_fig - t = t*m25 - d = t.div(k,m) - k = k+two - pi = pi + d - end - - d = one - k = one - t = BigDecimal("956") - while d.nonzero? && ((m = n - (pi.exponent - d.exponent).abs) > 0) - m = BigDecimal.double_fig if m < BigDecimal.double_fig - t = t.div(m57121,n) - d = t.div(k,m) - pi = pi + d - k = k+two - end - pi - end - - # call-seq: - # E(numeric) -> BigDecimal - # - # Computes e (the base of natural logarithms) to the specified number of - # digits of precision, +numeric+. - # - # BigMath.E(10).to_s - # #=> "0.271828182845904523536028752390026306410273e1" - # - def E(prec) - raise ArgumentError, "Zero or negative precision for E" if prec <= 0 - BigMath.exp(1, prec) - end -end diff --git a/ext/bigdecimal/lib/bigdecimal/newton.rb b/ext/bigdecimal/lib/bigdecimal/newton.rb deleted file mode 100644 index 85bacb7f2e..0000000000 --- a/ext/bigdecimal/lib/bigdecimal/newton.rb +++ /dev/null @@ -1,80 +0,0 @@ -# frozen_string_literal: false -require "bigdecimal/ludcmp" -require "bigdecimal/jacobian" - -# -# newton.rb -# -# Solves the nonlinear algebraic equation system f = 0 by Newton's method. -# This program is not dependent on BigDecimal. -# -# To call: -# n = nlsolve(f,x) -# where n is the number of iterations required, -# x is the initial value vector -# f is an Object which is used to compute the values of the equations to be solved. -# It must provide the following methods: -# -# f.values(x):: returns the values of all functions at x -# -# f.zero:: returns 0.0 -# f.one:: returns 1.0 -# f.two:: returns 2.0 -# f.ten:: returns 10.0 -# -# f.eps:: returns the convergence criterion (epsilon value) used to determine whether two values are considered equal. If |a-b| < epsilon, the two values are considered equal. -# -# On exit, x is the solution vector. -# -module Newton - include LUSolve - include Jacobian - module_function - - def norm(fv,zero=0.0) # :nodoc: - s = zero - n = fv.size - for i in 0...n do - s += fv[i]*fv[i] - end - s - end - - # See also Newton - def nlsolve(f,x) - nRetry = 0 - n = x.size - - f0 = f.values(x) - zero = f.zero - one = f.one - two = f.two - p5 = one/two - d = norm(f0,zero) - minfact = f.ten*f.ten*f.ten - minfact = one/minfact - e = f.eps - while d >= e do - nRetry += 1 - # Not yet converged. => Compute Jacobian matrix - dfdx = jacobian(f,f0,x) - # Solve dfdx*dx = -f0 to estimate dx - dx = lusolve(dfdx,f0,ludecomp(dfdx,n,zero,one),zero) - fact = two - xs = x.dup - begin - fact *= p5 - if fact < minfact then - raise "Failed to reduce function values." - end - for i in 0...n do - x[i] = xs[i] - dx[i]*fact - end - f0 = f.values(x) - dn = norm(f0,zero) - end while(dn>=d) - d = dn - end - nRetry - end -end diff --git a/ext/bigdecimal/lib/bigdecimal/util.rb b/ext/bigdecimal/lib/bigdecimal/util.rb deleted file mode 100644 index cb645d2a71..0000000000 --- a/ext/bigdecimal/lib/bigdecimal/util.rb +++ /dev/null @@ -1,181 +0,0 @@ -# frozen_string_literal: false -# -#-- -# bigdecimal/util extends various native classes to provide the #to_d method, -# and provides BigDecimal#to_d and BigDecimal#to_digits. -#++ - -require 'bigdecimal' - -class Integer < Numeric - # call-seq: - # int.to_d -> bigdecimal - # - # Returns the value of +int+ as a BigDecimal. - # - # require 'bigdecimal' - # require 'bigdecimal/util' - # - # 42.to_d # => 0.42e2 - # - # See also BigDecimal::new. - # - def to_d - BigDecimal(self) - end -end - - -class Float < Numeric - # call-seq: - # float.to_d -> bigdecimal - # float.to_d(precision) -> bigdecimal - # - # Returns the value of +float+ as a BigDecimal. - # The +precision+ parameter is used to determine the number of - # significant digits for the result (the default is Float::DIG). - # - # require 'bigdecimal' - # require 'bigdecimal/util' - # - # 0.5.to_d # => 0.5e0 - # 1.234.to_d(2) # => 0.12e1 - # - # See also BigDecimal::new. - # - def to_d(precision=0) - BigDecimal(self, precision) - end -end - - -class String - # call-seq: - # str.to_d -> bigdecimal - # - # Returns the result of interpreting leading characters in +str+ - # as a BigDecimal. - # - # require 'bigdecimal' - # require 'bigdecimal/util' - # - # "0.5".to_d # => 0.5e0 - # "123.45e1".to_d # => 0.12345e4 - # "45.67 degrees".to_d # => 0.4567e2 - # - # See also BigDecimal::new. - # - def to_d - BigDecimal.interpret_loosely(self) - end -end - - -class BigDecimal < Numeric - # call-seq: - # a.to_digits -> string - # - # Converts a BigDecimal to a String of the form "nnnnnn.mmm". - # This method is deprecated; use BigDecimal#to_s("F") instead. - # - # require 'bigdecimal/util' - # - # d = BigDecimal("3.14") - # d.to_digits # => "3.14" - # - def to_digits - if self.nan? || self.infinite? || self.zero? - self.to_s - else - i = self.to_i.to_s - _,f,_,z = self.frac.split - i + "." + ("0"*(-z)) + f - end - end - - # call-seq: - # a.to_d -> bigdecimal - # - # Returns self. - # - # require 'bigdecimal/util' - # - # d = BigDecimal("3.14") - # d.to_d # => 0.314e1 - # - def to_d - self - end -end - - -class Rational < Numeric - # call-seq: - # rat.to_d(precision) -> bigdecimal - # - # Returns the value as a BigDecimal. - # - # The required +precision+ parameter is used to determine the number of - # significant digits for the result. - # - # require 'bigdecimal' - # require 'bigdecimal/util' - # - # Rational(22, 7).to_d(3) # => 0.314e1 - # - # See also BigDecimal::new. - # - def to_d(precision) - BigDecimal(self, precision) - end -end - - -class Complex < Numeric - # call-seq: - # cmp.to_d -> bigdecimal - # cmp.to_d(precision) -> bigdecimal - # - # Returns the value as a BigDecimal. - # - # The +precision+ parameter is required for a rational complex number. - # This parameter is used to determine the number of significant digits - # for the result. - # - # require 'bigdecimal' - # require 'bigdecimal/util' - # - # Complex(0.1234567, 0).to_d(4) # => 0.1235e0 - # Complex(Rational(22, 7), 0).to_d(3) # => 0.314e1 - # - # See also BigDecimal::new. - # - def to_d(*args) - BigDecimal(self) unless self.imag.zero? # to raise eerror - - if args.length == 0 - case self.real - when Rational - BigDecimal(self.real) # to raise error - end - end - self.real.to_d(*args) - end -end - - -class NilClass - # call-seq: - # nil.to_d -> bigdecimal - # - # Returns nil represented as a BigDecimal. - # - # require 'bigdecimal' - # require 'bigdecimal/util' - # - # nil.to_d # => 0.0 - # - def to_d - BigDecimal(0) - end -end |