3
$\begingroup$

I have been given a project, I need to show the use of a version of Newton's method to solve these non-linear equations. The version of Newton's method I am required to use is: $ X_{n+1} = X_n - J^{(-1)} F(X_n) $. I have all the values required here, and this works to find the point $ X_1 $. However I need to find a code that inputs the following points $ X_2, X_3, ..., X_n $ automatically. Here are the given values:

f[x_, y_] := x^2 + y^2 - 5;
g[x_, y_] := x^3 - y^3 - 7;

x0 = 2.1;
y0 = 0.9;

f[x0, y0]
g[x0, y0]

0.22

1.532

M = {{2*x0, 2*y0}, {3*x0^2, -3*y0^2}}

{{4.2, 1.8}, {13.23, -2.43}}

J = Inverse[M]

{{0.0714286, 0.0529101}, {0.388889, -0.123457}}

F0 = {{f[x0, y0]}, {g[x0, y0]}}
X0 = {{x0}, {y0}}
X1 = X0 - J.F0

{{0.22}, {1.532}}
{{2.1}, {0.9}}
{{2.00323}, {1.00358}}

Thank you in Advance!

Improve this question
$\endgroup$
3
  • 3
    $\begingroup$ Have a look at FixedPoint and FixedPointList and their option SameTest. Nest(List) and NestWhile(List) also come to mind. You can also use FindRoot which has all this built-in. $\endgroup$
    – Henrik Schumacher
    Commented Dec 11, 2018 at 12:32
  • $\begingroup$ Ok Thank you, I will try some of the commands you have listed. Unfortunately, FindRoot does not help me in my situation as I need to show all the iterations leading up to the point, rather than just find the point itself. $\endgroup$
    – Andrew Bradley
    Commented Dec 11, 2018 at 13:02
  • $\begingroup$ Oh, this can also be done with FindRoot: Try Reap[FindRoot[Sin[x] == 0.2, {x, Pi/42}, EvaluationMonitor :> Sow[x]]]. Btw.: Using LinearSolve instead of Inverse should be faster for larger systems and should prevent certain problems with precision loss. $\endgroup$
    – Henrik Schumacher
    Commented Dec 11, 2018 at 14:20

2 Answers 2

Reset to default
2
$\begingroup$ 
jac = D[{f[x, y], g[x, y]}, {{x, y}, 1}];
Xlist = NestList[# - Inverse[jac /. Thread[{x, y} -> #]].{f @@ #, g @@ #} &, {x0, y0}, 5]

{{2.1, 0.9}, {2.0032275, 1.0035802}, {2.0000033, 1.0000051}, {2., 1.}, {2., 1.}, {2., 1.}}

You can get the same from FindRoot:

{res, {steps}} = Reap[FindRoot[{f[x, y], g[x, y]}, {{x, x0}, {y, y0}}, 
                       Method -> "Newton", StepMonitor :> Sow[{x, y}]]]

{{x -> 2., y -> 1.}, {{{2.0032275, 1.0035802}, {2.0000033, 1.0000051}, {2., 1.}, {2., 1.}}}}

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thank you so much, you are my hero! $\endgroup$
    – Andrew Bradley
    Commented Dec 11, 2018 at 13:48
1
$\begingroup$

I'll show it with FixedPointList

f[x_, y_] = {x^2 + y^2 - 5, x^3 - y^3 - 7};
j[x_, y_] = Grad[f[x, y], {x, y}];

with LinearSolve:

FixedPointList[(# - LinearSolve[j[Sequence @@ #], f[Sequence @@ #]]) &, {2.1, 0.9}, 3]
{{2.1, 0.9}, {2.00323, 1.00358}, {2., 1.00001}, {2., 1.}}

with Inverse:

FixedPointList[(# - Inverse[j[Sequence @@ #]].f[Sequence @@ #]) &, {2.1, 0.9}, 3]
{{2.1, 0.9}, {2.00323, 1.00358}, {2., 1.00001}, {2., 1.}}
Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.