I want to find the five unknown of a system of five nonlinear equations. Thank you !

Here the code:

clc;

clear;

close all;

clearvars;

syms ET i1 i2 i3 i4 i5 x1 x2 x3 x4 x5 % constants defined below

syms k1 k2 k3 k4 k7 % unknown variables

ET=1

i1=0.2

i2=0.32

i3=0.43

i4=0.3

i5=0.15

x1=100

x2=200

x3=400

x4=600

x5=800

eq1 = 0==-i1+k7*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1))*(ET-(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)))/(1+2*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)));

eq2 = 0==-i2+k7*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2))*(ET-(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)))/(1+2*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)));

eq3 = 0==-i3+k7*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3))*(ET-(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)))/(1+2*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)));

eq4 = 0==-i4+k7*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4))*(ET-(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)))/(1+2*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)));

eq5 = 0==-i5+k7*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5))*(ET-(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)))/(1+2*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)));

eq1=subs(eq1); eq2=subs(eq2); eq3=subs(eq3); eq4=subs(eq4); eq5=subs(eq5);

[k1,k2,k3,k4,k7]=fsolve( [eq1, eq2, eq3, eq4, eq5], [k1,k2,k3,k4,k7])

Answer 1

[Torsten]

function main
x0 = rand(5,1);
x = fsolve(@fun,x0)
end
function res = fun(x)
k1 = x(1);
k2 = x(2);
k3 = x(3);
k4 = x(4);
k7 = x(5);
ET=1;
i1=0.2;
i2=0.32;
i3=0.43;
i4=0.3;
i5=0.15;
x1=100;
x2=200;
x3=400;
x4=600;
x5=800;
res(1) = -i1+k7*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3 /(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1))*(ET-(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)))/(1+2*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)));
res(2) = -i2+k7*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2))*(ET-(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)))/(1+2*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)));
res(3) = -i3+k7*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3))*(ET-(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)))/(1+2*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)));
res(4) = -i4+k7*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4))*(ET-(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)))/(1+2*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)));
res(5) = -i5+k7*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5))*(ET-(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)))/(1+2*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)));
end

[Thierry Rebetez]

Dear Torsten,

Thank you very much for your answer.

I tried your code, it is much better than mine, but I still have a problem. I added this line of code:

options = optimoptions(@fsolve, 'MaxFunctionEvaluations', 1000000, 'MaxIterations', 1000000);

x = fsolve(@fun,x0, options);

I get the following message:

fsolve stopped because the relative size of the current step is less than the

default value of the step size tolerance squared, but the vector of function values

is not near zero as measured by the default value of the function tolerance.

Do you have an idea how to get rid of this issue? Thanks a lot for your invaluable help!

[Torsten]

What do you get when you include

res = fun(x)

after

x = fsolve(@fun,x0)

?

If all elements of res are small, then fsolve succeeded to find a solution.

Take a look here:

https://de.mathworks.com/matlabcentral/answers/164331-nonlinear-equation-solver-issue

[Thierry Rebetez]

Hi Torsten,

Thank you for your quick answer.

I put new values for i1 to i5 and for x1 to x5

(ET=0.01

i1=0.02;

i2=0.032;

i3=0.043;

i4=0.03;

i5=0.015;

x1=0.01;

x2=0.02;

x3=0.04;

x4=0.06;

x5=0.08;)

The equation is now solved, but the solution differ from time to time when I run the code again and the solutions are sometimes negative and/or imaginary.

In fact, I need only positive and real solutions for the k's. How can I implement this ?

Another question: I want to solve:

i1=k7*k3/(k4+k7)*(-(k2+k1*x1))+.....

or, either formulated:

0 == -i1+k7*k3/(k4+k7)*(-(k2+k1*x1))+....

and not:

res(1) = -i1+k7*k3/(k4+k7)*(-(k2+k1*x1

Thank you a lot for your invaluable help.

[Torsten]

Since you take sqrt's in your expressions, it's not guaranteed that the solution remains real-valued.

If it remains real-valued, you can get positive values by using this code:

function main
x0 = rand(5,1);
x = fsolve(@fun,x0);
k1 = x(1)^2
k2 = x(2)^2
k3 = x(3)^2
k4 = x(4)^2
k7 = x(5)^2
end
function res = fun(x)
k1 = x(1)^2;
k2 = x(2)^2;
k3 = x(3)^2;
k4 = x(4)^2;
k7 = x(5)^2;
ET=1;
i1=0.2;
i2=0.32;
i3=0.43;
i4=0.3;
i5=0.15;
x1=100;
x2=200;
x3=400;
x4=600;
x5=800;
res(1) = -i1+k7*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3 /(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1))*(ET-(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)))/(1+2*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)));
res(2) = -i2+k7*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2))*(ET-(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)))/(1+2*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)));
res(3) = -i3+k7*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3))*(ET-(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)))/(1+2*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)));
res(4) = -i4+k7*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4))*(ET-(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)))/(1+2*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)));
res(5) = -i5+k7*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5))*(ET-(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)))/(1+2*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)));
end

fsolve tries to adjust k1,k2,k3,k4 and k7 such that res(1),res(2),res(3),res(4) and res(5) become 0.

So res(i) = ... this is the correct syntax for fsolve.

Maybe you should read the fsolve-documentation, at least what is written under

https://de.mathworks.com/help/optim/ug/fsolve.html

[Thierry Rebetez]

Dear Torsten,

Thank you again for your answer!

I think now I am a step further. Here the code updated with different initial values closed to the reality of the study:

*********************************************************************************************************

function main

x0 = rand(5,1);

options = optimoptions(@fsolve, 'MaxFunctionEvaluations', 1000000, 'MaxIterations', 1000000,'OptimalityTolerance', 1e-07);

x = fsolve(@fun,x0,options)

res = fun(x)

k1 = x(1)^2;

k2 = x(2)^2;

k3 = x(3)^2;

k4 = x(4)^2;

k7 = x(5)^2;

end

function res = fun(x)

k1 = x(1)^2;

k2 = x(2)^2;

k3 = x(3)^2;

k4 = x(4)^2;

k7 = x(5)^2;

ET=1E-6

i1=2E-7;

i2=3.2E-7;

i3=4.3E-7;

i4=3E-7;

i5=1.5E-7;

x1=0.0001;

x2=0.0002;

x3=0.0004;

x4=0.0006;

x5=0.0008;

res(1) = -i1+k7*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3 /(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1))*(ET-(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)))/(1+2*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)));

res(2) = -i2+k7*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2))*(ET-(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)))/(1+2*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)));

res(3) = -i3+k7*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3))*(ET-(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)))/(1+2*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)));

res(4) = -i4+k7*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4))*(ET-(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)))/(1+2*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)));

res(5) = -i5+k7*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5))*(ET-(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)))/(1+2*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)));

end

I get for instance this result of values for k's:

*******************************************************************************************

Equation solved at initial point.

fsolve completed because the vector of function values at the initial point

is near zero as measured by the selected value of the function tolerance, and

the problem appears regular as measured by the gradient.

<stopping criteria details>

x =

0.3111

0.9234

0.4302

0.1848

0.9049

*******************************************************************************

Now I tried to put back the constants k's solved by Matlab in the equations, but the result is wrong. Why?

Another question: why, by each new resolution of the k's, Matlab get another solution ? Should it not be only one solution ?

Thank you very much for your help !!

[Torsten]

Now I tried to put back the constants k's solved by Matlab in the equations, but the result is wrong. Why?

It's nearly impossible that the initial values are solutions to the problem. So it is not surprising that the x values don't solve your equations. What do you get back with the command res = fun(x) ? Is the res-vector small ?

Another question: why, by each new resolution of the k's, Matlab get another solution ? Should it not be only one solution ?

A polynomial of degree n can have n solutions. Can you exclude such thing for the case you are working on ?

[Thierry Rebetez]

Dear Torsten,

Thank you very much for your invaluable help and your quick answer. I have added the line

res = fun(x)

after

x = fsolve(@fun,x0,options)

It resolves the problem partially, but I get the following message as in my last post. How to get rid of "Equation solved at initial point"?:

Equation solved at initial point.

fsolve completed because the vector of function values at the initial point

is near zero as measured by the selected value of the function tolerance, and

the problem appears regular as measured by the gradient.

<stopping criteria details>

x =

0.4387

0.3816

0.7655

0.7952

0.1869

ET =

1.0000e-06

res =

1.0e-06 *

-0.2000 -0.3200 -0.4300 -0.3500 -0.2000

Equation solved. The final point is the initial point.

The sum of squared function values, r = 4.898000e-13, is less than sqrt(options.FunctionTolerance) = 1.000000e-05.

The relative norm of the gradient of r, 1.527667e-21, is less than options.OptimalityTolerance = 1.000000e-10.

Optimization Metric Options

relative norm(grad r) = 1.53e-21 FunctionTolerance = 1e-10 (selected)

r = 4.90e-13 sqrt(OptimalityTolerance) = 1.0e-05 (selected)

******************************************************************

Here again the last code I used:

function main

x0 = rand(5,1);

options = optimoptions(@fsolve, 'MaxFunctionEvaluations', 1000000, 'MaxIterations', 1000000,'OptimalityTolerance', 1e-10);

x = fsolve(@fun,x0,options)

res = fun(x)

k1 = x(1)^2;

k2 = x(2)^2;

k3 = x(3)^2;

k4 = x(4)^2;

k7 = x(5)^2;

end

function res = fun(x)

k1 = x(1)^2;

k2 = x(2)^2;

k3 = x(3)^2;

k4 = x(4)^2;

k7 = x(5)^2;

ET=1E-6

i1=2E-7;

i2=3.2E-7;

i3=4.3E-7;

i4=3.5E-7;

i5=2E-7;

x1=0.0001;

x2=0.0002;

x3=0.0004;

x4=0.0006;

x5=0.0008;

res(1) = -i1+k7*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3 /(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1))*(ET-(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)))/(1+2*k3/(k4+k7)*(-(k2+k1*x1-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x1-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x1+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x1)));

res(2) = -i2+k7*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2))*(ET-(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)))/(1+2*k3/(k4+k7)*(-(k2+k1*x2-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x2-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x2+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x2)));

res(3) = -i3+k7*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3))*(ET-(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)))/(1+2*k3/(k4+k7)*(-(k2+k1*x3-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x3-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x3+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x3)));

res(4) = -i4+k7*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4))*(ET-(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)))/(1+2*k3/(k4+k7)*(-(k2+k1*x4-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x4-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x4+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x4)));

res(5) = -i5+k7*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5))*(ET-(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)))/(1+2*k3/(k4+k7)*(-(k2+k1*x5-k7*k3/(k4+k7)*ET)+sqrt((k2+k1*x5-k7*k3/(k4+k7)*ET)^2+4*k2*ET*(2*k1*k3/(k4+k7)*x5+k7*k3/(k4+k7))))/(2*k3/(k4+k7)*(k7+2*k1*x5)));

end

*************************************************************************************

Thank you very much!!

[Torsten]

Whatever value for x0 I choose, my solver gives the warning that the Jacobian matrix in the starting point is singular. Are you sure your equations are correct ? Where do they stem from ?

[Thierry Rebetez]

Do you want that I give you my email so we could talk more easily and I could send then my study and the way I come to these equations? I think they are correct however. Thank you for replying!

[Thierry Rebetez]

Hi Torsten

I think it is better to continue to chat on this platform if you want :) I think indeed that these equations are correct. The complete equation we have tried to solve is a combination of these three equations here:

h1 = (-(k2+k1.*x-k7.*K1.*ET)+sqrt((k2+k1.*x-k7.*K1.*ET).^2+4.*k2.*ET.*(2.*k1.*K1.*x+k7.*K1)))./(2.*K1.*(k7+2.*k1.*x));

i1=k7.*K1.*h1.*(ET-h1)./(1+2.*K1.*h1);

K1 = k3./(k4+k7);

with all k's positive or null.

I've tried with other values of ET, i1 to i5 and x1 to x5. It gives the following result:

ET =

0.0100

Equation solved, inaccuracy possible.

The vector of function values is near zero, as measured by the default value

of the function tolerance. However, the last step was ineffective.

<stopping criteria details>

x =

4.8110

1.9126

136.9516

-0.0000

43.0017

ET =

0.0100

res =

-0.0058 0.0065 -0.0034 0.0014 -0.0003

It is better but not perfect. What do you think about it?

Another question: the x which are given above (x =4.8110...) must be squared to get the k's, right?

Thank you very much for your help!