Using Newton's method for 2 equations for finding concentration at different points
1 visualizzazione (ultimi 30 giorni)
Mostra commenti meno recenti
Hi all, I'm currently trying to use Newton's method on two equations (I've attatched them) for finding the concentration of the lung and liver at three different R values 0.6, 0.8, and 0.9 from their intial guesses R = 0.6: Clung = 360.6 and Cliver = 284.6, R = 0.8: Clung = 312.25 and Cliver = 236.25, and R = 0.9: Clung = 215.5 and Cliver = 139.5. What I have so far is shown below, I'm hitting a wall and it wont allow me to create my J without giving an error if I use .*, etc. If I don't use .* in the J it gives me an error on dimensions instead. Any help would be greatly appreciated. Thank you
clearvars
clc
close all
%% Newton's Method to Systems
%% I have 2 equations, 3 sets of initials
maxiter = 100;
tol = 1e-6;
k = 0;
err(1) = 1;
%xold = [0.5; 0.25];
R_Old = [0.6; 0.8; 0.9];
C_Lung_Old = [360.6; 312.25; 215.5];
C_Liver_Old = [284.6; 236.25; 139.5];
Mat_Old = [R_Old,C_Lung_Old,C_Liver_Old];
while k<maxiter && err(k+1)>tol
k = k+1;
f1_1 = -780.*(1-R_Old(1))-R_Old(1).*(0.5.*C_Liver_Old(1)+1.5.*C_Lung_Old(1)) + (8.75.*C_Lung_Old(1)./(2.1 + C_Lung_Old(1))) + 2.*C_Lung_Old(1) ;
f1_2 = -780.*(1-R_Old(2))-R_Old(2).*(0.5.*C_Liver_Old(2)+1.5.*C_Lung_Old(2)) + (8.75.*C_Lung_Old(2)./(2.1 + C_Lung_Old(2))) + 2.*C_Lung_Old(2) ;
f1_3 = -780.*(1-R_Old(3))-R_Old(3).*(0.5.*C_Liver_Old(3)+1.5.*C_Lung_Old(3)) + (8.75.*C_Lung_Old(3)./(2.1 + C_Lung_Old(3))) + 2.*C_Lung_Old(3) ;
f2_1 = -0.5.*C_Lung_Old(1) + 0.322*(118*C_Liver_Old(1)./(7.0 +C_Liver_Old(1))) + 0.5*C_Liver_Old(1) ;
f2_2 = -0.5.*C_Lung_Old(2) + 0.322*(118*C_Liver_Old(2)./(7.0 +C_Liver_Old(2))) + 0.5*C_Liver_Old(2) ;
f2_3 = -0.5.*C_Lung_Old(3) + 0.322*(118*C_Liver_Old(3)./(7.0 +C_Liver_Old(3))) + 0.5*C_Liver_Old(3) ;
nF = [f1_1,f1_2,f1_3;f2_1,f2_2,f2_3];
J = [(C_Liver_Old+3*C_Lung_Old-1560)/2 (20*(3*R_Old-4)*C_Lung_Old*(5*C_Lung_Old+21)+1323*R_Old-5439)/(2*(10*C_Lung_Old+21).^2) 2*R_Old R_Old/2; ...
0 1/2 -((125*C_Liver_Old^2+1750*C_Liver_Old+72618)/(250*(C_Liver_Old+7)^2))];
xx = J\nF;
xnew = xx+Mat_Old;
f1 = -3*xnew(1).^2-xnew(2).^2 +10;
f2 = -xnew(1).^2-xnew(2)+1;
errtemp = [abs(f1) abs(f2)];
err(k+1) = max(errtemp);
Mat_Old = xnew;
end
semilogy(err)
return
0 Commenti
Risposta accettata
Torsten
il 7 Nov 2022
Modificato: Torsten
il 7 Nov 2022
Why don't you use "fsolve" ? Don't you have a license for the optimization toolbox ?
I don't understand J, f1 and f2 in your code.
maxiter = 100;
tol = 1e-6;
R_Old = [0.6; 0.8; 0.9];
C_Lung_Old = [360.6; 312.25; 215.5];
C_Liver_Old = [284.6; 236.25; 139.5];
sol = [];
for i = 1:3
f = @(x) [-780.*(1-R_Old(i))-R_Old(i)*(0.5*x(2)+1.5*x(1)) + (8.75*x(1)/(2.1 + x(1))) + 2*x(1) ; ...
-0.5*x(1) + 0.322*(118*x(2)/(7.0 +x(2))) + 0.5*x(2)] ;
Mat_Old = [C_Lung_Old(i);C_Liver_Old(i)];
err(1) = 1;
k = 0;
while k<maxiter && err(k+1)>tol
k = k+1;
fv = f(Mat_Old);
df = [(f(Mat_Old+1e-6*[1;0])-fv)/1e-6,(f(Mat_Old+1e-6*[0;1])-fv)/1e-6];
xx = df\fv;
xnew = Mat_Old-xx;
err(k+1) = norm(abs(fv));
Mat_Old = xnew;
end
sol = [sol,Mat_Old];
end
sol
0 Commenti
Più risposte (0)
Vedere anche
Categorie
Scopri di più su Communications Toolbox in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!