Newton's method for 2 dimension vectors

12 views (last 30 days)
Kevser Cifci on 18 Sep 2021
Commented: Alan Stevens on 20 Sep 2021
I want to write a code which gives me the roots of a vector function with sing the Newton's method.
that's what I wrote:
%%file newton2d.m
function [out] = newton2d (f,Jf,x,y, eps)
k=0;
Jf = jacobian(f, [x, y]) %error
h = f(x,y)./Jf(x,y); %error
while abs(h)< eps && (k < 1000)
h = f(x,y)./Jf(x,y);
[x,y] = [x,y] - h;
k = k+1;
end
out = [x,y];
end
But there is a problem with the jacobian: i don't know where to insert it?
What can I do to make newton2d function work?
Thank you !
Kevser Cifci on 18 Sep 2021
I changed a little bit my code:
%%file newton2d.m
function [out] = newton2d(x, y, f, Jf, eps)
k = 0;
h = Jf(x,y)\f(x,y)';
while abs(h)>=eps && k<1000
h = Jf(x,y)\f(x,y)'; %% error
[x,y] = [x,y]- h';
k = k+1;
end
out = [x,y];
end
Why it's still not working?
Then I create a vector function with it's Jacobian:
%%file f4.m
function out = f4(x,y)
out = [(x^3)-(3*x*(y^2))-1 , (3*(x^2)*y)-(y^3)]
end
%%file Jf4.m
function out = Jf4(x,y)
out = jacobian(f, [x, y]);
end
How to call f4 in my jacobian function?
newton2d(1,2,@f4, @Jf4, 10^(-8))
it's not working when I call my functions...
Thank you for your help !

Alan Stevens on 19 Sep 2021
LIke this?
% Functions
f = @(XY) [XY(1).^3 - 3*XY(1).*XY(2).^2 - 1;
3*XY(1).^2.*XY(2) - XY(2).^3];
J = @(XY) [3*XY(1).^2 - 3*XY(2).^2, -6*XY(1)*XY(2);
6*XY(1).*XY(2), 3*XY(1).^2 - 3*XY(2).^2];
% Initial guesses
XY = [2; 1];
tol = 10^-6;
err = 1;
maxits = 100;
its = 0;
while (err>tol) && (its<maxits)
XYold = XY;
XY = XYold - J(XYold)\f(XYold);
err = norm(XY-XYold);
its = its+1;
end
x = XY(1); y = XY(2);
disp('parameter values:')
parameter values:
disp(['x = ',num2str(x), ' y = ', num2str(y)])
x = 1 y = -1.0223e-17
disp('function values:' )
function values:
disp(f(XY))
1.0e-16 * 0 -0.3067
Alan Stevens on 20 Sep 2021
More like this then:
% Functions
f = @(xy) [xy(1).^3 - 3*xy(1).*xy(2).^2 - 1;
3*xy(1).^2.*xy(2) - xy(2).^3];
% Jacobian of f
J = @(xy) [3*xy(1).^2 - 3*xy(2).^2, -6*xy(1)*xy(2);
6*xy(1).*xy(2), 3*xy(1).^2 - 3*xy(2).^2];
% Initial guesses
xy = [2; 1];
xy = newton2d(f,J,xy);
x = xy(1); y = xy(2);
disp('parameter values:')
parameter values:
disp(['x = ',num2str(x), ' y = ', num2str(y)])
x = 1 y = -1.0223e-17
disp('function values:' )
function values:
disp(f(xy))
1.0e-16 * 0 -0.3067
function xy = newton2d(f,J,xy)
tol = 10^-6;
err = 1;
maxits = 100;
its = 0;
while (err>tol) && (its<maxits)
xyold = xy;
xy = xyold - J(xyold)\f(xyold);
err = norm(xy-xyold);
its = its+1;
end
end

R2019a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!