Rosetta Code

Multidimensional Newton-Raphson method: Difference between revisions

Page Discussion

Multidimensional Newton-Raphson method (view source)

Revision as of 14:10, 4 January 2024

48,222 bytes added ,  4 months ago
m
→‎{{header|Wren}}: Minor tidy and corrected a previous copy/paste error.
m (→‎{{header|Kotlin}}: Minor correction to pre-amble.)
m (→‎{{header|Wren}}: Minor tidy and corrected a previous copy/paste error.)
 
(14 intermediate revisions by 9 users not shown)
Line 3:
;Task:
Create a program that finds and outputs the root of a system of nonlinear equations
using Newton-Raphson metodmethod.
<br><br>
 
=={{header|C#}}==
=={{header|C sharp|C#}}==
For matrix inversion and matrix and vector definitions - see C# source from [[Gaussian elimination]]
<langsyntaxhighlight lang="csharp">
using System;
 
Line 42 ⟶ 43:
}
}
</syntaxhighlight>
</lang>
<langsyntaxhighlight lang="csharp">
using System;
 
Line 99 ⟶ 100:
}
}
</syntaxhighlight>
</lang>
{{out}}<pre>
2.54258545959024
1.06149981539336
</pre>
 
=={{header|Go}}==
{{trans|Kotlin}}
<br>
We follow the Kotlin example of coding our own matrix methods rather than using a third party library.
<syntaxhighlight lang="go">package main
 
import (
"fmt"
"math"
)
 
type vector = []float64
type matrix []vector
type fun = func(vector) float64
type funs = []fun
type jacobian = []funs
 
func (m1 matrix) mul(m2 matrix) matrix {
rows1, cols1 := len(m1), len(m1[0])
rows2, cols2 := len(m2), len(m2[0])
if cols1 != rows2 {
panic("Matrices cannot be multiplied.")
}
result := make(matrix, rows1)
for i := 0; i < rows1; i++ {
result[i] = make(vector, cols2)
for j := 0; j < cols2; j++ {
for k := 0; k < rows2; k++ {
result[i][j] += m1[i][k] * m2[k][j]
}
}
}
return result
}
 
func (m1 matrix) sub(m2 matrix) matrix {
rows, cols := len(m1), len(m1[0])
if rows != len(m2) || cols != len(m2[0]) {
panic("Matrices cannot be subtracted.")
}
result := make(matrix, rows)
for i := 0; i < rows; i++ {
result[i] = make(vector, cols)
for j := 0; j < cols; j++ {
result[i][j] = m1[i][j] - m2[i][j]
}
}
return result
}
 
func (m matrix) transpose() matrix {
rows, cols := len(m), len(m[0])
trans := make(matrix, cols)
for i := 0; i < cols; i++ {
trans[i] = make(vector, rows)
for j := 0; j < rows; j++ {
trans[i][j] = m[j][i]
}
}
return trans
}
 
func (m matrix) inverse() matrix {
le := len(m)
for _, v := range m {
if len(v) != le {
panic("Not a square matrix")
}
}
aug := make(matrix, le)
for i := 0; i < le; i++ {
aug[i] = make(vector, 2*le)
copy(aug[i], m[i])
// augment by identity matrix to right
aug[i][i+le] = 1
}
aug.toReducedRowEchelonForm()
inv := make(matrix, le)
// remove identity matrix to left
for i := 0; i < le; i++ {
inv[i] = make(vector, le)
copy(inv[i], aug[i][le:])
}
return inv
}
 
// note: this mutates the matrix in place
func (m matrix) toReducedRowEchelonForm() {
lead := 0
rowCount, colCount := len(m), len(m[0])
for r := 0; r < rowCount; r++ {
if colCount <= lead {
return
}
i := r
 
for m[i][lead] == 0 {
i++
if rowCount == i {
i = r
lead++
if colCount == lead {
return
}
}
}
 
m[i], m[r] = m[r], m[i]
if div := m[r][lead]; div != 0 {
for j := 0; j < colCount; j++ {
m[r][j] /= div
}
}
 
for k := 0; k < rowCount; k++ {
if k != r {
mult := m[k][lead]
for j := 0; j < colCount; j++ {
m[k][j] -= m[r][j] * mult
}
}
}
lead++
}
}
 
func solve(fs funs, jacob jacobian, guesses vector) vector {
size := len(fs)
var gu1 vector
gu2 := make(vector, len(guesses))
copy(gu2, guesses)
jac := make(matrix, size)
for i := 0; i < size; i++ {
jac[i] = make(vector, size)
}
tol := 1e-8
maxIter := 12
iter := 0
for {
gu1 = gu2
g := matrix{gu1}.transpose()
t := make(vector, size)
for i := 0; i < size; i++ {
t[i] = fs[i](gu1)
}
f := matrix{t}.transpose()
for i := 0; i < size; i++ {
for j := 0; j < size; j++ {
jac[i][j] = jacob[i][j](gu1)
}
}
g1 := g.sub(jac.inverse().mul(f))
gu2 = make(vector, size)
for i := 0; i < size; i++ {
gu2[i] = g1[i][0]
}
iter++
any := false
for i, v := range gu2 {
if math.Abs(v)-gu1[i] > tol {
any = true
break
}
}
if !any || iter >= maxIter {
break
}
}
return gu2
}
 
func main() {
/*
solve the two non-linear equations:
y = -x^2 + x + 0.5
y + 5xy = x^2
given initial guesses of x = y = 1.2
 
Example taken from:
http://www.fixoncloud.com/Home/LoginValidate/OneProblemComplete_Detailed.php?problemid=286
 
Expected results: x = 1.23332, y = 0.2122
*/
f1 := func(x vector) float64 { return -x[0]*x[0] + x[0] + 0.5 - x[1] }
f2 := func(x vector) float64 { return x[1] + 5*x[0]*x[1] - x[0]*x[0] }
fs := funs{f1, f2}
jacob := jacobian{
funs{
func(x vector) float64 { return -2*x[0] + 1 },
func(x vector) float64 { return -1 },
},
funs{
func(x vector) float64 { return 5*x[1] - 2*x[0] },
func(x vector) float64 { return 1 + 5*x[0] },
},
}
guesses := vector{1.2, 1.2}
sol := solve(fs, jacob, guesses)
fmt.Printf("Approximate solutions are x = %.7f, y = %.7f\n", sol[0], sol[1])
 
/*
solve the three non-linear equations:
9x^2 + 36y^2 + 4z^2 - 36 = 0
x^2 - 2y^2 - 20z = 0
x^2 - y^2 + z^2 = 0
given initial guesses of x = y = 1.0 and z = 0.0
 
Example taken from:
http://mathfaculty.fullerton.edu/mathews/n2003/FixPointNewtonMod.html (exercise 5)
 
Expected results: x = 0.893628, y = 0.894527, z = -0.0400893
*/
 
fmt.Println()
f3 := func(x vector) float64 { return 9*x[0]*x[0] + 36*x[1]*x[1] + 4*x[2]*x[2] - 36 }
f4 := func(x vector) float64 { return x[0]*x[0] - 2*x[1]*x[1] - 20*x[2] }
f5 := func(x vector) float64 { return x[0]*x[0] - x[1]*x[1] + x[2]*x[2] }
fs = funs{f3, f4, f5}
jacob = jacobian{
funs{
func(x vector) float64 { return 18 * x[0] },
func(x vector) float64 { return 72 * x[1] },
func(x vector) float64 { return 8 * x[2] },
},
funs{
func(x vector) float64 { return 2 * x[0] },
func(x vector) float64 { return -4 * x[1] },
func(x vector) float64 { return -20 },
},
funs{
func(x vector) float64 { return 2 * x[0] },
func(x vector) float64 { return -2 * x[1] },
func(x vector) float64 { return 2 * x[2] },
},
}
guesses = vector{1, 1, 0}
sol = solve(fs, jacob, guesses)
fmt.Printf("Approximate solutions are x = %.7f, y = %.7f, z = %.7f\n", sol[0], sol[1], sol[2])
}</syntaxhighlight>
 
{{out}}
<pre>
Approximate solutions are x = 1.2333178, y = 0.2122450
 
Approximate solutions are x = 0.8936282, y = 0.8945270, z = -0.0400893
</pre>
 
=={{header|Julia}}==
NLsolve is a Julia package for nonlinear systems of equations, with the Newton-Raphson method one of the choices for solvers.
<syntaxhighlight lang="julia"># from the NLSolve documentation: to solve
# (x, y) -> ((x+3)*(y^3-7)+18, sin(y*exp(x)-1))
using NLsolve
 
function f!(F, x)
F[1] = (x[1]+3)*(x[2]^3-7)+18
F[2] = sin(x[2]*exp(x[1])-1)
end
 
function j!(J, x)
J[1, 1] = x[2]^3-7
J[1, 2] = 3*x[2]^2*(x[1]+3)
u = exp(x[1])*cos(x[2]*exp(x[1])-1)
J[2, 1] = x[2]*u
J[2, 2] = u
end
 
println(nlsolve(f!, j!, [ 0.1; 1.2], method = :newton))
</syntaxhighlight>{{out}}
<pre>
Results of Nonlinear Solver Algorithm
* Algorithm: Newton with line-search
* Starting Point: [0.1, 1.2]
* Zero: [-3.7818e-16, 1.0]
* Inf-norm of residuals: 0.000000
* Iterations: 4
* Convergence: true
* |x - x'| < 0.0e+00: false
* |f(x)| < 1.0e-08: true
* Function Calls (f): 5
* Jacobian Calls (df/dx): 4
</pre>
 
Line 109 ⟶ 392:
 
As neither the JDK nor the Kotlin Standard Library have matrix functions built in, most of the functions used have been taken from other tasks.
<langsyntaxhighlight lang="scala">// Version 1.2.31
 
import kotlin.math.abs
Line 287 ⟶ 570:
val (xx2, yy2, zz2) = solve(funcs2, jacobian2, guesses2)
System.out.printf("Approximate solutions are x = %.7f, y = %.7f, z = %.7f\n", xx2, yy2, zz2)
}</langsyntaxhighlight>
 
{{out}}
Line 294 ⟶ 577:
 
Approximate solutions are x = 0.8936282, y = 0.8945270, z = -0.0400893
</pre>
 
=={{header|Nim}}==
{{trans|Kotlin}}
<syntaxhighlight lang="nim">import sequtils, strformat, sugar
 
type
Vector = seq[float]
Matrix = seq[Vector]
Func = (Vector) -> float
Funcs = seq[Func]
Jacobian = seq[Funcs]
 
 
func `*`(m1, m2: Matrix): Matrix =
let
rows1 = m1.len
cols1 = m1[0].len
rows2 = m2.len
cols2 = m2[0].len
doAssert cols1 == rows2
result = newSeqWith(rows1, newSeq[float](cols2))
for i in 0..<rows1:
for j in 0..<cols2:
for k in 0..<rows2:
result[i][j] += m1[i][k] * m2[k][j]
 
 
func `-`(m1, m2: Matrix): Matrix =
let
rows = m1.len
cols = m1[0].len
doAssert m2.len == rows and m2[0].len == cols
result = newSeqWith(rows, newSeq[float](cols))
for i in 0..<rows:
for j in 0..<cols:
result[i][j] = m1[i][j] - m2[i][j]
 
 
func transposed(m: Matrix): Matrix =
let
rows = m.len
cols = m[0].len
result = newSeqWith(cols, newSeq[float](rows))
for i in 0..<cols:
for j in 0..<rows:
result[i][j] = m[j][i]
 
 
func toReducedRowEchelonForm(m: var Matrix) =
var lead = 0
let rowCount = m.len
let colCount = m[0].len
for r in 0..<rowCount:
if colCount <= lead: return
var i = r
 
while m[i][lead] == 0:
inc i
if rowCount == i:
i = r
inc lead
if colCount == lead: return
 
swap m[i], m[r]
 
if m[r][lead] != 0:
let divisor = m[r][lead]
for j in 0..<colCount: m[r][j] /= divisor
 
for k in 0..<rowCount:
if k != r:
let mult = m[k][lead]
for j in 0..<colCount: m[k][j] -= m[r][j] * mult
 
inc lead
 
 
func inverse(m: Matrix): Matrix =
let size = m.len
doAssert m.allIt(it.len == size), "not a square matrix."
var aug = newSeqWith(size, newSeq[float](2 * size))
for i in 0..<size:
for j in 0..<size: aug[i][j] = m[i][j]
# Augment by identity matrix to right.
aug[i][i + size] = 1
aug.toReducedRowEchelonForm()
result = newSeqWith(size, newSeq[float](size))
# Remove identity matrix to left.
for i in 0..<size:
for j in 0..<size: result[i][j] = aug[i][j + size]
 
 
proc solve(funcs: Funcs; jacobian: Jacobian; guesses: Vector): Vector =
let size = funcs.len
result = guesses
var jac = newSeqWith(size, newSeq[float](size))
const Tol = 1e-8
let MaxIter = 12
var iter = 1
while true:
let gu = move(result)
let g = transposed(@[gu])
let f = transposed(@[funcs.mapIt(it(gu))])
for i in 0..<size:
for j in 0..<size:
jac[i][j] = jacobian[i][j](gu)
let g1 = g - inverse(jac) * f
result = g1.mapIt(it[0])
inc iter
if iter > MaxIter: break
var exit = true
for idx, val in result:
if abs(val - gu[idx]) > Tol:
exit = false
break
if exit: break
 
 
when isMainModule:
 
#[ Solve the two non-linear equations:
y = -x^2 + x + 0.5
y + 5xy = x^2
given initial guesses of x = y = 1.2
 
Example taken from:
http://www.fixoncloud.com/Home/LoginValidate/OneProblemComplete_Detailed.php?problemid=286
 
Expected results: x = 1.23332, y = 0.2122
]#
 
let
f1: Func = (x: Vector) => -x[0] * x[0] + x[0] + 0.5 - x[1]
f2: Func = (x: Vector) => x[1] + 5 * x[0] * x[1] - x[0] * x[0]
funcs1: Funcs = @[f1, f2]
jacobian1: Jacobian = @[@[Func((x: Vector) => - 2 * x[0] + 1),
Func((x: Vector) => -1.0)],
@[Func((x: Vector) => 5 * x[1] - 2 * x[0]),
Func((x: Vector) => 1 + 5 * x[0])]]
 
guesses1 = @[1.2, 1.2]
sol1 = solve(funcs1, jacobian1, guesses1)
echo &"Approximate solutions are x = {sol1[0]:.7f}, y = {sol1[1]:.7f}"
 
#[ Solve the three non-linear equations:
9x^2 + 36y^2 + 4z^2 - 36 = 0
x^2 - 2y^2 - 20z = 0
x^2 - y^2 + z^2 = 0
given initial guesses of x = y = 1.0 and z = 0.0
 
Example taken from:
http://mathfaculty.fullerton.edu/mathews/n2003/FixPointNewtonMod.html (exercise 5)
 
Expected results: x = 0.893628, y = 0.894527, z = -0.0400893
]#
 
echo()
let
f3: Func = (x: Vector) => 9 * x[0] * x[0] + 36 * x[1] * x[1] + 4 * x[2] * x[2] - 36
f4: Func = (x: Vector) => x[0] * x[0] - 2 * x[1] * x[1] - 20 * x[2]
f5: Func = (x: Vector) => x[0] * x[0] - x[1] * x[1] + x[2] * x[2]
funcs2: Funcs = @[f3, f4, f5]
jacobian2: Jacobian = @[@[Func((x: Vector) => 18 * x[0]),
Func((x: Vector) => 72 * x[1]),
Func((x: Vector) => 8 * x[2])],
@[Func((x: Vector) => 2 * x[0]),
Func((x: Vector) => -4 * x[1]),
Func((x: Vector) => -20.0)],
@[Func((x: Vector) => 2 * x[0]),
Func((x: Vector) => -2 * x[1]),
Func((x: Vector) => 2 * x[2])]]
guesses2 = @[1.0, 1.0, 0.0]
sol2 = solve(funcs2, jacobian2, guesses2)
echo &"Approximate solutions are x = {sol2[0]:.7f}, y = {sol2[1]:.7f}, z = {sol2[2]:.7f}"</syntaxhighlight>
 
{{out}}
<pre>Approximate solutions are x = 1.2333178, y = 0.2122450
 
Approximate solutions are x = 0.8936282, y = 0.8945270, z = -0.0400893</pre>
 
=={{header|Phix}}==
{{trans|Go}}
Uses code from [[Reduced_row_echelon_form#Phix]],
[[Gauss-Jordan_matrix_inversion#Phix]],
[[Matrix_transposition#Phix]], and
[[Matrix_multiplication#Phix]]<br>
See std distro for a complete runnable version.
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #000080;font-style:italic;">-- demo\rosetta\Multidimensional_Newton-Raphson_method.exw</span>
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">solve</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">fs</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">jacob</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">guesses</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">size</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">),</span>
<span style="color: #000000;">maxIter</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">12</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">iter</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">gu1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">g</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">t</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">f</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">g1</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">gu2</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">guesses</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">jac</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">size</span><span style="color: #0000FF;">),</span><span style="color: #000000;">size</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">tol</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1e-8</span>
<span style="color: #008080;">while</span> <span style="color: #004600;">true</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">gu1</span> <span style="color: #0000FF;">:=</span> <span style="color: #000000;">gu2</span>
<span style="color: #000000;">g</span> <span style="color: #0000FF;">:=</span> <span style="color: #000000;">matrix_transpose</span><span style="color: #0000FF;">({</span><span style="color: #000000;">gu1</span><span style="color: #0000FF;">})</span>
<span style="color: #000000;">t</span> <span style="color: #0000FF;">:=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">size</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">size</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">t</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">:=</span> <span style="color: #7060A8;">call_func</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">],{</span><span style="color: #000000;">gu1</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #000000;">f</span> <span style="color: #0000FF;">:=</span> <span style="color: #000000;">matrix_transpose</span><span style="color: #0000FF;">({</span><span style="color: #000000;">t</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">size</span> <span style="color: #008080;">do</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">size</span> <span style="color: #008080;">do</span>
<span style="color: #000000;">jac</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">:=</span> <span style="color: #7060A8;">call_func</span><span style="color: #0000FF;">(</span><span style="color: #000000;">jacob</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">],{</span><span style="color: #000000;">gu1</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #000000;">g1</span> <span style="color: #0000FF;">:=</span> <span style="color: #7060A8;">sq_sub</span><span style="color: #0000FF;">(</span><span style="color: #000000;">g</span><span style="color: #0000FF;">,</span><span style="color: #000000;">matrix_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">inverse</span><span style="color: #0000FF;">(</span><span style="color: #000000;">jac</span><span style="color: #0000FF;">),</span><span style="color: #000000;">f</span><span style="color: #0000FF;">))</span>
<span style="color: #000000;">gu2</span> <span style="color: #0000FF;">:=</span> <span style="color: #7060A8;">vslice</span><span style="color: #0000FF;">(</span><span style="color: #000000;">g1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">iter</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
<span style="color: #004080;">bool</span> <span style="color: #000000;">any</span> <span style="color: #0000FF;">:=</span> <span style="color: #7060A8;">find</span><span style="color: #0000FF;">(</span><span style="color: #004600;">true</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">sq_gt</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">sq_sub</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">sq_abs</span><span style="color: #0000FF;">(</span><span style="color: #000000;">gu2</span><span style="color: #0000FF;">),</span><span style="color: #000000;">gu1</span><span style="color: #0000FF;">),</span><span style="color: #000000;">tol</span><span style="color: #0000FF;">))!=</span><span style="color: #000000;">0</span>
<span style="color: #008080;">if</span> <span style="color: #008080;">not</span> <span style="color: #000000;">any</span> <span style="color: #008080;">or</span> <span style="color: #000000;">iter</span> <span style="color: #0000FF;">>=</span> <span style="color: #000000;">maxIter</span> <span style="color: #008080;">then</span> <span style="color: #008080;">exit</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">gu2</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">f1</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">x</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">+</span><span style="color: #000000;">x</span><span style="color: #0000FF;">+</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">-</span><span style="color: #000000;">y</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">f2</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">y</span><span style="color: #0000FF;">+</span><span style="color: #000000;">5</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">-</span><span style="color: #000000;">x</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">f3</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">9</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">+</span><span style="color: #000000;">36</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">+</span><span style="color: #000000;">4</span><span style="color: #0000FF;">*</span><span style="color: #000000;">z</span><span style="color: #0000FF;">*</span><span style="color: #000000;">z</span><span style="color: #0000FF;">-</span><span style="color: #000000;">36</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">f4</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">x</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">-</span><span style="color: #000000;">20</span><span style="color: #0000FF;">*</span><span style="color: #000000;">z</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">f5</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">x</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">-</span><span style="color: #000000;">y</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">+</span><span style="color: #000000;">z</span><span style="color: #0000FF;">*</span><span style="color: #000000;">z</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j1</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span><span style="color: #0000FF;">+</span><span style="color: #000000;">1</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j2</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000080;font-style:italic;">/*v*/</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">return</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j3</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">5</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span><span style="color: #0000FF;">-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j4</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">+</span><span style="color: #000000;">5</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j11</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">18</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j12</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{?,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">72</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j13</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{?,?,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">8</span><span style="color: #0000FF;">*</span><span style="color: #000000;">z</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j21</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j22</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{?,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">4</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j23</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000080;font-style:italic;">/*v*/</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">return</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">20</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j31</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">x</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j32</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{?,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">y</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">j33</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">v</span><span style="color: #0000FF;">)</span> <span style="color: #004080;">atom</span> <span style="color: #0000FF;">{?,?,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">v</span> <span style="color: #008080;">return</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">*</span><span style="color: #000000;">z</span> <span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">procedure</span> <span style="color: #000000;">main</span><span style="color: #0000FF;">()</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">fs</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">jacob</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">guesses</span>
<span style="color: #000080;font-style:italic;">/*
solve the two non-linear equations:
y = -x^2 + x + 0.5
y + 5xy = x^2
given initial guesses of x = y = 1.2
Example taken from:
http://www.fixoncloud.com/Home/LoginValidate/OneProblemComplete_Detailed.php?problemid=286
Expected results: x = 1.23332, y = 0.2122
*/</span>
<span style="color: #000000;">fs</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">f1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">f2</span><span style="color: #0000FF;">}</span>
<span style="color: #000000;">jacob</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{{</span><span style="color: #000000;">j1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j2</span><span style="color: #0000FF;">},</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">j3</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j4</span><span style="color: #0000FF;">}}</span>
<span style="color: #000000;">guesses</span> <span style="color: #0000FF;">:=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">1.2</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1.2</span><span style="color: #0000FF;">}</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Approximate solutions are x = %.7f, y = %.7f\n\n"</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">solve</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">jacob</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">guesses</span><span style="color: #0000FF;">))</span>
<span style="color: #000080;font-style:italic;">/*
solve the three non-linear equations:
9x^2 + 36y^2 + 4z^2 - 36 = 0
x^2 - 2y^2 - 20z = 0
x^2 - y^2 + z^2 = 0
given initial guesses of x = y = 1.0 and z = 0.0
Example taken from:
http://mathfaculty.fullerton.edu/mathews/n2003/FixPointNewtonMod.html (exercise 5)
Expected results: x = 0.893628, y = 0.894527, z = -0.0400893
*/</span>
<span style="color: #000000;">fs</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">f3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">f4</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">f5</span><span style="color: #0000FF;">}</span>
<span style="color: #000000;">jacob</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{{</span><span style="color: #000000;">j11</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j12</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j13</span><span style="color: #0000FF;">},</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">j21</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j22</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j23</span><span style="color: #0000FF;">},</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">j31</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j32</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j33</span><span style="color: #0000FF;">}}</span>
<span style="color: #000000;">guesses</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Approximate solutions are x = %.7f, y = %.7f, z = %.7f\n"</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">solve</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">jacob</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">guesses</span><span style="color: #0000FF;">))</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
<span style="color: #000000;">main</span><span style="color: #0000FF;">()</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
Approximate solutions are x = 1.2333178, y = 0.2122450
 
Approximate solutions are x = 0.8936282, y = 0.8945270, z = -0.04008929
</pre>
 
=={{header|Raku}}==
(formerly Perl 6)
<syntaxhighlight lang="raku" line># Reference:
# https://github.com/pierre-vigier/Perl6-Math-Matrix
# Mastering Algorithms with Perl
# By Jarkko Hietaniemi, John Macdonald, Jon Orwant
# Publisher: O'Reilly Media, ISBN-10: 1565923987
# https://resources.oreilly.com/examples/9781565923980/blob/master/ch16/solve
 
use v6;
 
sub solve_funcs ($funcs, @guesses, $iterations, $epsilon) {
my ($error_value, @values, @delta, @jacobian); my \ε = $epsilon;
for 1 .. $iterations {
for ^+$funcs { @values[$^i] = $funcs[$^i](|@guesses); }
$error_value = 0;
for ^+$funcs { $error_value += @values[$^i].abs }
return @guesses if $error_value ≤ ε;
for ^+$funcs { @delta[$^i] = -@values[$^i] }
@jacobian = jacobian $funcs, @guesses, ε;
@delta = solve_matrix @jacobian, @delta;
loop (my $j = 0, $error_value = 0; $j < +$funcs; $j++) {
$error_value += @delta[$j].abs ;
@guesses[$j] += @delta[$j];
}
return @guesses if $error_value ≤ ε;
}
return @guesses;
}
 
sub jacobian ($funcs is copy, @points is copy, $epsilon is copy) {
my ($Δ, @P, @M, @plusΔ, @minusΔ);
my Array @jacobian; my \ε = $epsilon;
for ^+@points -> $i {
@plusΔ = @minusΔ = @points;
$Δ = (ε * @points[$i].abs) || ε;
@plusΔ[$i] = @points[$i] + $Δ;
@minusΔ[$i] = @points[$i] - $Δ;
for ^+$funcs { @P[$^k] = $funcs[$^k](|@plusΔ); }
for ^+$funcs { @M[$^k] = $funcs[$^k](|@minusΔ); }
for ^+$funcs -> $j {
@jacobian[$j][$i] = (@P[$j] - @M[$j]) / (2 * $Δ);
}
}
return @jacobian;
}
 
sub solve_matrix (@matrix_array is copy, @delta is copy) {
# https://github.com/pierre-vigier/Perl6-Math-Matrix/issues/56
{ use Math::Matrix;
my $matrix = Math::Matrix.new(@matrix_array);
my $vector = Math::Matrix.new(@delta.map({.list}));
die "Matrix is not invertible" unless $matrix.is-invertible;
my @result = ( $matrix.inverted dot $vector ).transposed;
return @result.split(" ");
}
}
 
my $funcs = [
{ 9*$^x² + 36*$^y² + 4*$^z² - 36 }
{ $^x² - 2*$^y² - 20*$^z }
{ $^x² - $^y² + $^z² }
];
 
my @guesses = (1,1,0);
 
my @solution = solve_funcs $funcs, @guesses, 20, 1e-8;
 
say "Solution: ", @solution;
</syntaxhighlight>
{{out}}
<pre>
Solution: [0.8936282344764825 0.8945270103905782 -0.04008928615915281]
</pre>
 
=={{header|Wren}}==
{{trans|Kotlin}}
{{libheader|Wren-matrix}}
{{libheader|Wren-fmt}}
<syntaxhighlight lang="wren">import "./matrix" for Matrix
import "./fmt" for Fmt
 
var solve = Fn.new { |funcs, jacobian, guesses|
var size = funcs.count
var gu1 = []
var gu2 = guesses.toList
var jac = Matrix.new(size, size)
var tol = 1e-8
var maxIter = 12
var iter = 0
while (true) {
gu1 = gu2
var g = Matrix.new([gu1]).transpose
var v = List.filled(size, 0)
for (i in 0...size) v[i] = funcs[i].call(gu1)
var f = Matrix.new([v]).transpose
for (i in 0...size) {
for (j in 0...size) jac[i, j] = jacobian[i][j].call(gu1)
}
var g1 = g - jac.inverse * f
gu2 = List.filled(size, 0)
for (i in 0...size) gu2[i] = g1[i][0]
iter = iter + 1
if (iter == maxIter) break
if ((0...size).all { |i| (gu2[i] - gu1[i]).abs <= tol }) break
}
return gu2
}
 
/* solve the two non-linear equations:
y = -x^2 + x + 0.5
y + 5xy = x^2
given initial guesses of x = y = 1.2
 
Example taken from:
http://www.fixoncloud.com/Home/LoginValidate/OneProblemComplete_Detailed.php?problemid=286
 
Expected results: x = 1.23332, y = 0.2122
*/
var f1 = Fn.new { |x| -x[0] * x[0] + x[0] + 0.5 - x[1] }
var f2 = Fn.new { |x| x[1] + 5 * x[0] * x[1] - x[0] * x[0] }
var funcs = [f1, f2]
var jacobian = [
[ Fn.new { |x| - 2 * x[0] + 1 }, Fn.new { |x| -1 } ],
[ Fn.new { |x| 5 * x[1] - 2 * x[0] }, Fn.new { |x| 1 + 5 * x[0] } ]
]
var guesses = [1.2, 1.2]
var sols = solve.call(funcs, jacobian, guesses)
Fmt.print("Approximate solutions are x = $.7f, y = $.7f", sols[0], sols[1])
 
/* solve the three non-linear equations:
9x^2 + 36y^2 + 4z^2 - 36 = 0
x^2 - 2y^2 - 20z = 0
x^2 - y^2 + z^2 = 0
given initial guesses of x = y = 1.0 and z = 0.0
 
Example taken from:
http://mathfaculty.fullerton.edu/mathews/n2003/FixPointNewtonMod.html (exercise 5)
 
Expected results: x = 0.893628, y = 0.894527, z = -0.0400893
*/
System.print()
var f3 = Fn.new { |x| 9 * x[0] * x[0] + 36 * x[1] *x[1] + 4 * x[2] * x[2] - 36 }
var f4 = Fn.new { |x| x[0] * x[0] - 2 * x[1] * x[1] - 20 * x[2] }
var f5 = Fn.new { |x| x[0] * x[0] - x[1] * x[1] + x[2] * x[2] }
funcs = [f3, f4, f5]
jacobian = [
[ Fn.new { |x| 18 * x[0] }, Fn.new { |x| 72 * x[1] }, Fn.new { |x| 8 * x[2] } ],
[ Fn.new { |x| 2 * x[0] }, Fn.new { |x| -4 * x[1] }, Fn.new { |x| -20 } ],
[ Fn.new { |x| 2 * x[0] }, Fn.new { |x| -2 * x[1] }, Fn.new { |x| 2 * x[2] } ]
]
guesses = [1, 1, 0]
sols = solve.call(funcs, jacobian, guesses)
Fmt.print("Approximate solutions are x = $.7f, y = $.7f, z = $.7f", sols[0], sols[1], sols[2])</syntaxhighlight>
 
{{out}}
<pre>
Approximate solutions are x = 1.2333178, y = 0.2122450
 
Approximate solutions are x = 0.8936282, y = 0.8945270, z = -0.0400893
</pre>
 
=={{header|zkl}}==
This doesn't use Newton-Raphson (with derivatives) but a hybrid algorithm.
<langsyntaxhighlight lang="zkl">var [const] GSL=Import.lib("zklGSL"); // libGSL (GNU Scientific Library)
 
// two functions of two variables: f(x,y)=0
Line 306 ⟶ 1,039:
v.format(11,8).println(); // answer overwrites initial guess
 
xy:=fs.run(True,v.toList().xplode()).println(); // Vectordeltas tofrom Listzero</syntaxhighlight>
fs.apply('wrap(f){ f(xy.xplode()) }).println(); // deltas from zero</lang>
{{out}}
<pre>
Line 314 ⟶ 1,046:
</pre>
A condensed solver (for a different set of functions):
<langsyntaxhighlight lang="zkl">v:=GSL.VectorFromData(-10.0, -15.0);
GSL.multiroot_fsolver(T( fcn(x,y){ 1.0 - x }, fcn(x,y){ 10.0*(y - x*x) }),v)
.format().println(); // --> (1,1)</langsyntaxhighlight>
{{out}}
<pre>
Line 322 ⟶ 1,054:
</pre>
Another example:
<langsyntaxhighlight lang="zkl">v:=GSL.VectorFromData(1.0, 1.0, 0.0); // initial guess
fxyzs:=T(
fcn(x,y,z){ x*x*9 + y*y*36 + z*z*4 - 36 }, // 9x^2 + 36y^2 + 4z^2 - 36 = 0
Line 329 ⟶ 1,061:
(v=GSL.multiroot_fsolver(fxyzs,v)).format(12,8).println();
 
fxyzs.run(True,v.toList().xplode()).println(); // deltas from zero</syntaxhighlight>
xyz:=v.toList();
fxyzs.apply('wrap(f){ f(xyz.xplode()) }).println(); // deltas from zero</lang>
{{out}}
<pre>
PureFox
9,476

edits

Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/279861...353378"