diff options
Diffstat (limited to 'vendor/bundle/ruby/3.4.0/gems/bigdecimal-4.1.2/sample/nlsolve.rb')
| -rw-r--r-- | vendor/bundle/ruby/3.4.0/gems/bigdecimal-4.1.2/sample/nlsolve.rb | 57 |
1 files changed, 57 insertions, 0 deletions
diff --git a/vendor/bundle/ruby/3.4.0/gems/bigdecimal-4.1.2/sample/nlsolve.rb b/vendor/bundle/ruby/3.4.0/gems/bigdecimal-4.1.2/sample/nlsolve.rb new file mode 100644 index 0000000..7953653 --- /dev/null +++ b/vendor/bundle/ruby/3.4.0/gems/bigdecimal-4.1.2/sample/nlsolve.rb @@ -0,0 +1,57 @@ +# +# nlsolve.rb +# An example for solving nonlinear algebraic equation system. +# + +require "bigdecimal" +require_relative "linear" + +# Requires gem matrix +require "matrix" + +# :stopdoc: + +def func((x, y)) # defines functions solved + [ + x**2 + y**2 - 2, + (x - 1)**2 + (y + 1)**2 - 3 + ] +end + +def jacobian(x, f, delta, prec) + dim = x.size + dim.times.map do |i| + xplus = Array.new(dim) {|j| x[i] + (j == i ? delta : 0) } + xminus = Array.new(dim) {|j| x[i] - (j == i ? delta : 0) } + yplus = f.call(xplus) + yminus = f.call(xminus) + yplus.zip(yminus).map {|p, m| (p - m).div(2 * delta, prec) } + end.transpose +end + +def nlsolve(initial_x, prec:, max_iteration: 100, &f) + initial_x = initial_x.map {|v| BigDecimal(v) } + x = initial_x + delta = BigDecimal(0.01) + calc_prec = prec + 10 + max_iteration.times do |iteration| + # Newton step + jacobian = jacobian(x, f, delta, calc_prec) + matrix = Matrix[*jacobian.map {|row| row.map {|v| PrecisionSpecifiedValue.new(v, calc_prec) } }] + y = f.call(x) + vector = y.map {|v| PrecisionSpecifiedValue.new(v, calc_prec) } + dx = matrix.lup.solve(vector).map(&:value) + x_prev = x + x = x.zip(dx).map {|xi, di| xi.sub(di, prec) } + movement = x_prev.zip(x).map {|xn, xi| (xn - xi).abs }.max + delta = [movement, delta].min.mult(1, 10) + break if movement.zero? || movement.exponent < -prec + end + x +end + +initial_value = [1, 1] +ans = nlsolve(initial_value, prec: 100) {|x| func(x) } +diff = func(ans).map {|v| v.mult(1, 10) } +p(ans:) +p(diff:) |
