# File bigdecimal/newton.rb, line 39
39:   def nlsolve(f,x)
40:     nRetry = 0
41:     n = x.size
42: 
43:     f0 = f.values(x)
44:     zero = f.zero
45:     one  = f.one
46:     two  = f.two
47:     p5 = one/two
48:     d  = norm(f0,zero)
49:     minfact = f.ten*f.ten*f.ten
50:     minfact = one/minfact
51:     e = f.eps
52:     while d >= e do
53:       nRetry += 1
54:       # Not yet converged. => Compute Jacobian matrix
55:       dfdx = jacobian(f,f0,x)
56:       # Solve dfdx*dx = -f0 to estimate dx
57:       dx = lusolve(dfdx,f0,ludecomp(dfdx,n,zero,one),zero)
58:       fact = two
59:       xs = x.dup
60:       begin
61:         fact *= p5
62:         if fact < minfact then
63:           raize "Failed to reduce function values."
64:         end
65:         for i in 0...n do
66:           x[i] = xs[i] - dx[i]*fact
67:         end
68:         f0 = f.values(x)
69:         dn = norm(f0,zero)
70:       end while(dn>=d)
71:       d = dn
72:     end
73:     nRetry
74:   end

    # File bigdecimal/newton.rb, line 30
30:   def norm(fv,zero=0.0)
31:     s = zero
32:     n = fv.size
33:     for i in 0...n do
34:       s += fv[i]*fv[i]
35:     end
36:     s
37:   end