Combining my newton's method function for nonlinear equations and my lu decomposition function.
4 views (last 30 days)
Show older comments
Hi, I'm trying to combine my newton's method function for nonlinear equations with my LU decomposition function.
This is my function for newton's method:
function x = newtonSysGB(J, fun, x0, xtol, ftol, max_iter)
F_value = fun(x0);
for i = 1:max_iter
delta = J(x0)\-F_value;
Delta_norm = norm(delta);
x0 = x0 + delta;
F_value = fun(x0);
if cond(fun(x0))>10^5
fprintf('Warning!')
end
if cond(fun(x0))>10^15
error('Condition number too high!')
end
F_norm = norm(F_value);
if F_norm < ftol || Delta_norm < xtol
x = x0;
break
end
if F_norm > ftol || Delta_norm > xtol
disp ('Need more iterations!')
end
end
J and fun are two other m-files which are:
function F_vector = fun(x)
F_vector = [x(1)*exp(0.1*x(2)) + 0.1*x(3) - 100; x(1)*exp(0.2*x(2)) + 0.2*x(3) - 120; x(1)*exp(0.3*x(2)) + 0.3*x(3) - 150];
end
function J_matrix = J(x)
J_matrix = zeros(3,3);
J_matrix(1,1) = exp(0.1*x(2)); % df1x1
J_matrix(1,2) = 0.1*exp(0.1*x(2))*x(1); % df1x2
J_matrix(1,3) = 0.1; % df1x3
J_matrix(2,1) = exp(0.2*x(2)); % df2x1
J_matrix(2,2) = 0.2*exp(0.2*x(2))*x(1); % df2x2;
J_matrix(2,3) = 0.2; % df2x3;
J_matrix(3,1) = exp(0.3*x(2)); % df3x1
J_matrix(3,2) = 0.3*exp(0.3*x(2))*x(1); % df3x2
J_matrix(3,3) = 0.3; % df3x3
end
I want to change the bolded part of it, which is just like a system of linear equations to solving it using this LU decomposition function I have.
This is the LU decomposition function that I have:
function X_sol = LUgauss(A,b)
n = length(A);
L = zeros(n);
U = zeros(n);
P = eye(n);
for k=1:n
[~,r] = max(abs(A(k:end,k)));
r = n-(n-k+1)+r;
A([k r],:) = A([r k],:);
P([k r],:) = P([r k],:);
L([k r],:) = L([r k],:);
L(k:n,k) = A(k:n,k) / A(k,k);
U(k,1:n) = A(k,1:n);
A(k+1:n,1:n) = A(k+1:n,1:n) - L(k+1:n,k)*A(k,1:n);
end
U(:,end) = A(:,end);
B = P*b;
n = length(B);
y(1,1) = B(1)/L(1,1);
for i = 2:n
y(i,1) = (B(i)-L(i,1:i-1)*y(1:i-1,1))./L(i,i);
end
n = length(y);
X_sol(n,1) = y(n)/U(n,n);
for i = n-1:-1:1
X_sol(i,1) = (y(i)-U(i,i+1:n)*X_sol(i+1:n,1))./U(i,i);
end
end
THANK YOU SOOO MUCH!
0 Comments
Answers (1)
Fabio Freschi
on 8 Nov 2019
Edited: Fabio Freschi
on 8 Nov 2019
Is this what you are asking for?
delta = LUgauss(J(x0),-F_value);
Note that the buil-in lu Matlab function can do the work for you
2 Comments
See Also
Categories
Find more on Calculus in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!