Clear Filters
Clear Filters

Using lsqcurvefit with fsolve

4 views (last 30 days)
Dursman Mchabe
Dursman Mchabe on 4 Sep 2021
Edited: Matt J on 6 Sep 2021
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 Comments
Dursman Mchabe
Dursman Mchabe on 6 Sep 2021
Edited: Dursman Mchabe on 6 Sep 2021
Hi Matt,
Thanks a lot for the comment.
Param is of size 1 x 1, whereas x is of size 1 x 4.
Thanks
D
Matt J
Matt J on 6 Sep 2021
Edited: Matt J on 6 Sep 2021
And? How does tha invalidate my suggestion?

Sign in to comment.

Answers (0)

Community Treasure Hunt

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

Start Hunting!