102 4 Non-linear Equations
together with range reduction , is used in the intrisic computational function which computes
square roots.
Newton’s method can be generalized to systems of several non-linear equations and vari-
ables. Consider the case with two equations
f 1 (x 1 ,x 2 ) = 0
f 2 (x 1 ,x 2 ) = 0, (4.32)
which we Taylor expand to obtain
0 =f 1 (x 1 +h 1 ,x 2 +h 2 ) =f 1 (x 1 ,x 2 )+h 1 ∂f 1 /∂x 1 +h 2 ∂f 1 /∂x 2 +...
0 =f 2 (x 1 +h 1 ,x 2 +h 2 ) =f 2 (x 1 ,x 2 )+h 1 ∂f 2 /∂x 1 +h 2 ∂f 2 /∂x 2 +.... (4.33)Defining the Jacobian matrixJˆwe have
Jˆ=(
∂f 1 /∂x 1 ∂f 1 /∂x 2
∂f 2 /∂x 1 ∂f 2 /∂x 2)
, (4.34)
we can rephrase Newton’s method as
(
xn 1 +^1
xn 2 +^1
)
=
(
xn 1
xn 2)
+
(
hn 1
hn 2)
, (4.35)
where we have defined (
hn 1
hn 2
)
=−Jˆ−^1
(
f 1 (xn 1 ,xn 2 )
f 2 (xn 1 ,xn 2 ))
. (4.36)
We need thus to compute the inverse of the Jacobian matrix andit is to understand that
difficulties may arise in caseJˆis nearly singular.
It is rather straightforward to extend the above scheme to systems of more than two non-
linear equations.
The code for Newton-Raphson’s method can look like this
/*
This function
calculates a root between x1 and x2 of a function pointed to
*by (funcd) using the Newton-Raphson method. The user-defined
function funcd() returns both the function value and its first
derivative at the point x,
*The root is returned with an accuracy of +- xacc.
/
doublenewtonraphson(void(funcd)(double,double,double),doublex1,doublex2,
doublexacc)
{
int j;
doubledf, dx, f, rtn;
rtn = 0.5(x1 + x2); // initial guess
for(j = 0; j < max_iterations; j++){
(funcd)(rtn, &f, &df);
dx = f/df;
rtn -= dx;
if((x1 - rtn)(rtn - x2) < 0.0){
cout <<"\n\nError in function newtonraphson:"<< endl ;
cout <<"Jump out of interval bracket"<< endl;
}
if(fabs(dx) < xacc)returnrtn;
}