Using lsqcurvefit with fsolve
3 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Hi everyone,
I am trying to solve a nonlinear equation using fsolve. The equation is:
for i = 1:length(t)
sol(i) = fsolve(@(x) (Param + x)-(K5/x + 2.*(ydata3(i).*K3.*K4./(x^2 + x.*K3 + K3.*K4)) +...
(ydata3(i).*x.*K3./(x^2 + x.*K3 + K3.*K4))+ ...
2.*(Param.*K1.*K2./(x^2 + x.*K1 + K1.*K2)) +...
(ydata3(i).*x.*K1./(x^2 + x.*K1 + K1.*K2))+ ...
(ydata2(i).*K6)./(K6 + x)), x0);
end
However, the parameter, Param, is not known, I have to regress it from lsqcurve fitting using:
Param0 = 0.035;
Param = lsqcurvefit(@chargeBal,Param0,t,ydata1);
Where the YDATA will have the size of 4 x 1, whereas the Param, has a size of 1 x 1.
How can I make this to work?
The complete code is:
function pH_Example
global t pH K1 K2 K3 K4 K5 K6 ydata1 ydata2 ydata3 Param sol
pH = [11.48
6.86
6.72
6.25];
t = [0
30
60
100];
ydata3 = [2.58E-01
6.28E-02
5.76E-02
4.83E-02];
ydata2 = [ 1.25E-02
9.75E-03
7.93E-03
5.83E-03];
ydata1 = [3.31131E-12
1.38038E-07
1.38038E-07
1.38038E-07];
K1 = 1.09e-6;
K2 = 3.69e-10;
K3 = 6.2e-4;
K4 = 3.11e-8;
K5 = 3.02e-14;
K6 = (1.91e-2 + 5.38e-4)/2;
function H = chargeBal(Param,t)
x0 = 10^(-11.48);
for i = 1:length(t)
sol(i) = fsolve(@(x) (Param + x)-(K5/x + 2.*(ydata3(i).*K3.*K4./(x^2 + x.*K3 + K3.*K4)) +...
(ydata3(i).*x.*K3./(x^2 + x.*K3 + K3.*K4))+ ...
2.*(Param.*K1.*K2./(x^2 + x.*K1 + K1.*K2)) +...
(ydata3(i).*x.*K1./(x^2 + x.*K1 + K1.*K2))+ ...
(ydata2(i).*K6)./(K6 + x)), x0);
end
function sol = pHfun(t,x)
K1 = 1.09e-6;
K2 = 3.69e-10;
K3 = 6.2e-4;
K4 = 3.11e-8;
K5 = 3.02e-14;
K6 = (1.91e-2 + 5.38e-4)/2;
x0 = 10^(-11.48);
t = [0
30
60
100];
pH = [11.48
6.86
6.72
6.25];
sol = (Param + x)-(K5./x + 2.*(ydata3.*K3.*K4./(x^2 + x.*K3 + K3.*K4)) +...
(ydata3.*x.*K3./(x^2 + x.*K3 + K3.*K4))+ ...
2.*(Param.*K1.*K2./(x^2 + x.*K1 + K1.*K2)) +...
(ydata3.*x.*K1./(x^2 + x.*K1 + K1.*K2))+ ...
(ydata2.*K6)./(K6 + x));
end
H = sol(i);
end
Param0 = 0.035;
Param = lsqcurvefit(@chargeBal,Param0,t,ydata1);
fprintf(1,'\tParameters:\n')
for k1 = 1:length(Param)
fprintf(1, '\t\tParam(%d) = %8.5f\n', k1, Param(k1))
end
ph = 3 - log10(sol(i));
plot(t,ph,'r--')
hold on
plot (t,pH,'kd')
xlabel('Time (sec)')
ylabel('pH')
hold off
legend('pH-Mod','pH-Exp')
end
3 commentaires
Réponses (0)
Voir également
Catégories
En savoir plus sur Nonlinear Least Squares (Curve Fitting) dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!