Equation solving

I am implementing an equation solver into Omicron.

If I have an equation on the form:

f(x) = 0

If quite easy to solve using Newton-Raphson's method. This part of the equation solver works fine.

This equation system is a bit trickier:

f(x,y) = 0
g(x,y) = 0

That should be possible to solve using a similar method. But as of yet, I haven't studied any multi-dimensional analysis, so could someone please help me with an algorithm?

And what about the general case:

f1(x1,x2,...,xn) = 0
f2(x1,x2,...,xn) = 0
...
fn(x1,x2,...,xn) = 0

I think your still able to apply newton's
method to higher dimensions. I'm not sure what
changes you will have to make but probably
there is some way of taking the partial derivatives and
using them.

Hmm, I might be able to figure it out.
Instead of calculating the intersection of lines, you calculate the intersection of planes.

I'll do some seaching.

Now I got it.
I found an algorithm by looking at these two pages:
http://www.economics.unimelb.edu.au/...n-raphson1.pdf
http://uranus.ee.auth.gr/lessons/1/10.html

This is the actual code:

Code:

//INPUT
//Equation to solve
// f(x,y) == 0
// g(x,y) == 0

//Test equations
function f(x,y) = x^2 - y
function g(x,y) = 2x + 3y - 16

//guess the solution
x = 8
y = 40
//END INPUT

h = 0.0000001
for count = 1...10
  //Get the derivatives
  dfx = (f(x+h,y) - f(x-h,y))/(2h)
  dgx = (g(x+h,y) - g(x-h,y))/(2h)
  dfy = (f(x,y+h) - f(x,y-h))/(2h)
  dgy = (g(x,y+h) - g(x,y-h))/(2h)

  //Create the Jacobian matrix
  J = matrix(2,2)
  J[1] = (dfx , dfy)
  J[2] = (dgx , dgy)
  //Create the matrix with the function values
  F = matrix(2,1)
  F[1,1] = f(x,y)
  F[2,1] = g(x,y)

  //The matrix with the previous solutions
  X = matrix(2,1)
  X[1,1] = x
  X[2,1] = y

  //The solution matrix is now calcualted
  S = X - invert(J) * F

  //Obtain the variables from the solution matrix
  x =  S[1,1]
  y =  S[2,1]

next

println "x="+x
println "y="+y
println "f(x,y)="+f(x,y)
println "g(x,y)="+g(x,y)