# How can I solve a system of non-linear equations for positive values only?

9 views (last 30 days)
Julia S on 14 Aug 2020
Commented: Julia S on 28 Sep 2020
Hello,
I wrote a code to solve a system of non-linear equations, but I am not getting the right answer when I run it. I thought the initial guess was the problem, so I used experimental values from a research paper, but the solution is not even close it. Additionally, I get a comment that the solver was stopped prematurely. I am not very experienced in Matlab. Could someone help me find what I am doing wrong?
function f = root(x)
%data
delta = [0.1]*1e-4; QA = 3.57; QB = 20; QC = 60; QD = 1000; ph = [9.86923]*76; pl = [0.986923]*76;
nF = 277778; yAF = 0.7841; yBF = 0.2084; yCF = 0.0003; yDF = 1-yAF-yBF-yCF; A = 2260000;
% Material Balance
f(1) = x(13)*x(1) + x(14)*x(5) - nF*yAF;
f(2) = x(13)*x(2) + x(14)*x(6) - nF*yBF;
f(3) = x(13)*x(3) + x(14)*x(7) - nF*yCF;
f(4) = x(13)*x(4) + x(14)*x(8) - nF*yDF;
% Material Balance in terms of flux
f(5) = x(9)*A - x(13)*x(1);
f(6) = x(10)*A - x(13)*x(2);
f(7) = x(11)*A - x(13)*x(3);
f(8) = x(12)*A - x(13)*x(4);
% Flux Equations
f(9) = x(9) - (QA/delta)*(ph*x(5) - pl*x(1));
f(10) = x(10) - (QB/delta)*(ph*x(6) - pl*x(2));
f(11) = x(11) - (QC/delta)*(ph*x(7) - pl*x(3));
f(12) = x(12) - (QD/delta)*(ph*x(8) - pl*x(4));
% Component Equations
f(13) = x(5) + x(6) + x(7) + x(8) - 1;
f(14) = x(1) + x(2) + x(3) + x(4) - 1;
end
followed by:
x0 = [0.42, 0.48, 0.05, 0.05, 0.825, 0.075, 0.05, 0.05, 0.00258, 0.00295, 0.00031, 00031, 13888.9, 263889]';
[x, fval] = fsolve(@root, x0)
The output I get is:
Solver stopped prematurely.
fsolve stopped because it exceeded the function evaluation limit,
options.MaxFunctionEvaluations = 1.400000e+03.
x =
1.0e+11 *
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
1.6223
-1.9271

Ayush Gupta on 3 Sep 2020
Edited: Ayush Gupta on 3 Sep 2020
fsolve doesn’t allow bounds but lsqnonlin does. Change the last piece of code to get answers for positive values only.
x = lsqnonlin(@root,x0,zeros(size(x0)))
##### 2 CommentsShow 1 older commentHide 1 older comment
Julia S on 28 Sep 2020
Thank you so much!