4 Nonlinear equations solve for 4 variables

%script [ error- fsolve stopped because the problem appears to be locally singular.]

%solving E1,E2,E3,Noi

%Initial guesses

%E1o = 0.00009;

%E2o = 0.00009;

%E3o = 0.0009;

%Noio = 0.1;

X0 = [E1o; E2o; E3o; Noio];

%call fsolve function

%Call fsolve function

fhandle = @fsolve_function;

opts=optimoptions('fsolve','Display','iter','TolFun',1e-30,'TolX',1e-30);

[X,fval,exitflag] = fsolve (fhandle, X0,opts);

disp ('E1, E2, E3, Noi')

%equilibrium composition- function

function F=fsolve_function(X)

%unpack variables

E1=X(1);

E2=X(2);

E3=X(3);

Noi=X(4);%initial value of oxygen

K1=4.669E20;

K2=1.29;

K3=3.44E19;

yo2eq=0.00906;%atm

%Initial guesses

E1o = 0.00009;

E2o = 0.00009;

E3o = 0.0009;

Noio = 0.1;

%input equations

eq1=((2*E1-E2)^2/((Noi-E1)*(Noi+E1-E2+2*E3)))-K1;

eq2=(E2*(E2-2*E3)/(2*E1-E2)/(2*E3-E2))-K2;

eq3=((2*E3-E2)^2*(Noi+E1-E2+2*E3)/(E2-2*E3)^2/(Noi-E1))-K3;

eq4=((Noi-E1-E3)/(Noi+E1-E2+2*E3))-yo2eq;

F=[eq1; eq2; eq3; eq4; ];

end

Answer 1

[Torsten]

%Initial guesses
E1o = 0.00009; 
E2o = 0.00009; 
E3o = 0.0009; 
Noio = 0.1; 
X0 = [E1o; E2o; E3o; Noio];

%call fsolve function
%Call fsolve function
fhandle = @fsolve_function;
opts=optimoptions('fsolve','Display','iter','TolFun',1e-30,'TolX',1e-30);
[X,fval,exitflag] = fsolve (fhandle, X0,opts)

                                         Norm of      First-order   Trust-region
 Iteration  Func-count     f(x)          step         optimality    radius
     0          5     2.25508e+37                      4.47e+38               1
     1         10      1.3183e+36       0.102144       5.42e+37               1
     2         15     8.23905e+34      0.0517217       6.78e+36               1
     3         20     5.14896e+33       0.025735       8.48e+35               1
     4         25     3.21746e+32      0.0127713       1.06e+35               1
     5         30     2.00999e+31     0.00631321       1.33e+34               1
     6         35     1.25491e+30     0.00310237       1.67e+33               1
     7         40     7.82371e+28     0.00151075       2.11e+32               1
     8         45     4.86093e+27    0.000725366       2.66e+31               1
     9         50     2.99428e+26    0.000340595        3.4e+30               1
    10         55     1.80325e+25     0.00015428        4.4e+29               1
    11         60     1.01923e+24    6.58533e-05       5.82e+28               1
    12         65     4.74074e+22    2.54535e-05       7.81e+27               1
    13         70     1.15438e+21    8.24765e-06       9.18e+26               1
    14         75     3.58404e+18    1.72945e-06       4.71e+25               1
    15         80     2.66772e+11    5.38742e-08       1.28e+22               1
    16         85     3.34614e+06    9.76111e-08       4.53e+19               1
    17         90         1.03933    1.96575e-08       2.52e+16               1
    18         95     3.45209e-07    1.98259e-08       1.46e+13               1
    19        100     8.51577e-18    1.67752e-20       7.23e+07               1

Equation solved, solver stalled.

fsolve stopped because the relative size of the current step is less than the
value of the step size tolerance squared and the vector of function values
is near zero as measured by the value of the function tolerance.

X = 4×1

1.0e-04 * -0.0783 -0.6969 -0.0048 -0.0783

fval = 4×1

1.0e-08 * 0.2918 0.0000 0.0000 0.0000

exitflag = 2

disp ('E1, E2, E3, Noi')

E1, E2, E3, Noi

%equilibrium composition- function
function F=fsolve_function(X)
%unpack variables
E1=X(1);
E2=X(2);
E3=X(3);
Noi=X(4);%initial value of oxygen
K1=4.669E20;
K2=1.29;
K3=3.44E19;
yo2eq=0.00906;%atm

%input equations
eq1=(2*E1-E2)^2-K1*(Noi-E1)*(Noi+E1-E2+2*E3);
eq2=E2*(E2-2*E3)-K2*(2*E1-E2)*(2*E3-E2);
eq3=(2*E3-E2)^2*(Noi+E1-E2+2*E3)-K3*(E2-2*E3)^2*(Noi-E1);
eq4=(Noi-E1-E3)-yo2eq*(Noi+E1-E2+2*E3);
F=[eq1; eq2; eq3; eq4; ];
end

[Mahe Rukh]

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

value of the step size tolerance squared and the vector of function values

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

How to solve this? and all the variables should have poisitve values.

[Torsten]

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

value of the step size tolerance squared and the vector of function values

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

This is a message of success - fsolve found a solution.

How to solve this? and all the variables should have poisitve values.

If you need positive solutions, try different initial values or straight away "MultiStart" with the following code.

%Initial guesses
E1o = 0.00009; 
E2o = 0.00009; 
E3o = 0.0009; 
Noio = 0.1; 
X0 = [E1o; E2o; E3o; Noio];
X0 = sqrt(X0);
%Call fsolve function
fhandle = @fsolve_function;
opts=optimoptions('fsolve','Display','iter','TolFun',1e-30,'TolX',1e-30);
[X,fval,exitflag] = fsolve (fhandle, X0,opts)

                                         Norm of      First-order   Trust-region
 Iteration  Func-count     f(x)          step         optimality    radius
     0          5     2.25508e+37                      2.83e+38               1
     1          6     2.25508e+37              1       2.83e+38               1
     2         11      1.9304e+34           0.25       1.86e+36            0.25
     3         12      1.9304e+34       0.225686       1.86e+36           0.625
     4         17     1.21718e+34      0.0564215       1.23e+36          0.0564
     5         22     6.75648e+33      0.0564215       6.41e+35          0.0564
     6         27      6.9413e+32      0.0479351       1.16e+35          0.0564
     7         32      6.2024e+31      0.0551061       1.91e+34          0.0564
     8         37     5.43061e+30      0.0470679       3.14e+33           0.138
     9         42     4.68477e+29      0.0389483        5.2e+32           0.138
    10         47     3.99873e+28      0.0321637       8.64e+31           0.138
    11         52     3.38917e+27      0.0266262       1.44e+31           0.138
    12         57     2.85981e+26      0.0220939       2.41e+30           0.138
    13         62     2.40677e+25       0.018363       4.02e+29           0.138
    14         67     2.02263e+24      0.0152779       6.73e+28           0.138
    15         72     1.69894e+23      0.0127192       1.13e+28           0.138
    16         77     1.42749e+22      0.0105931       1.89e+27           0.138
    17         82     1.20097e+21     0.00882445       3.16e+26           0.138
    18         87     1.01325e+20     0.00735215       5.31e+25           0.138
    19         92     8.59611e+18     0.00612601       8.93e+24           0.138
    20         97     7.37005e+17     0.00510462       1.51e+24           0.138
    21        102     6.44621e+16     0.00425365       2.58e+23           0.138
    22        107     5.85144e+15     0.00354461       4.49e+22           0.138
    23        112     5.67861e+14      0.0029538       8.07e+21           0.138
    24        117     6.16912e+13     0.00246147       1.54e+21           0.138
    25        122     7.95278e+12     0.00205121        3.2e+20           0.138
    26        127     1.27873e+12     0.00170934       7.47e+19           0.138
    27        132      2.5927e+11     0.00142445       1.95e+19           0.138
    28        137     6.38719e+10     0.00118704       5.63e+18           0.138
    29        142     1.80017e+10    0.000989197       1.73e+18           0.138
    30        147     5.50032e+09     0.00082433       5.54e+17           0.138
    31        152     1.75696e+09    0.000686942       1.81e+17           0.138
    32        157     5.74453e+08    0.000572451          6e+16           0.138
    33        162     1.90077e+08    0.000477043          2e+16           0.138
    34        167     6.32736e+07    0.000397536       6.67e+15           0.138
    35        172     2.11268e+07     0.00033128       2.23e+15           0.138
    36        177     7.06494e+06    0.000276066       7.46e+14           0.138
    37        182     2.36437e+06    0.000230055        2.5e+14           0.138
    38        187          791578    0.000191713       8.37e+13           0.138
    39        192          265070    0.000159761        2.8e+13           0.138
    40        197         88771.7    0.000133134       9.39e+12           0.138
    41        202         29731.6    0.000110945       3.14e+12           0.138
    42        207         9958.26    9.24541e-05       1.05e+12           0.138
    43        212         3335.54    7.70451e-05       3.53e+11           0.138
    44        217         1117.29    6.42042e-05       1.18e+11           0.138
    45        222          374.27    5.35035e-05       3.96e+10           0.138
    46        227         125.379    4.45863e-05       1.33e+10           0.138
    47        232         42.0043    3.71552e-05       4.44e+09           0.138
    48        237         14.0732    3.09627e-05       1.49e+09           0.138
    49        242          4.7155    2.58022e-05       4.99e+08           0.138
    50        247         1.58019    2.15019e-05       1.67e+08           0.138
    51        252     2.03083e-07    6.05235e-09          6e+04           0.138
    52        253     2.03083e-07    1.79133e-05          6e+04           0.138
    53        254     2.03083e-07    4.47833e-06          6e+04        4.48e-06
    54        255     2.03083e-07    1.11958e-06          6e+04        1.12e-06
    55        256     2.03083e-07    2.79896e-07          6e+04         2.8e-07
    56        257     2.03083e-07     6.9974e-08          6e+04           7e-08
    57        258     2.03083e-07    1.74935e-08          6e+04        1.75e-08
    58        263      2.3587e-08    4.37337e-09       2.04e+04        4.37e-09
    59        264      2.3587e-08    4.37337e-09       2.04e+04        4.37e-09
    60        269     1.47466e-09    1.09334e-09       5.11e+03        1.09e-09
    61        270     1.47466e-09    2.73336e-09       5.11e+03        2.73e-09
    62        275     5.75459e-10     6.8334e-10       3.19e+03        6.83e-10
    63        276     5.75459e-10     6.8334e-10       3.19e+03        6.83e-10
    64        281     3.63573e-11    1.70835e-10            798        1.71e-10
    65        282     3.63573e-11    4.27087e-10            798        4.27e-10
    66        287     1.44114e-11    1.06772e-10            499        1.07e-10
    67        288     1.44114e-11    1.06772e-10            499        1.07e-10
    68        293     1.25046e-12     2.6693e-11            125        2.67e-11
    69        294     1.25046e-12    6.67324e-11            125        6.67e-11
    70        299     7.14699e-13    1.66831e-11           77.9        1.67e-11
    71        300     7.14699e-13    1.66831e-11           77.9        1.67e-11
    72        305     3.93389e-13    4.17077e-12           19.5        4.17e-12
    73        306     3.93389e-13    1.04269e-11           19.5        1.04e-11
    74        311     3.80309e-13    2.60673e-12           12.2        2.61e-12
    75        312     3.80309e-13    2.60673e-12           12.2        2.61e-12
    76        317     3.72465e-13    6.51683e-13           3.04        6.52e-13
    77        318     3.72465e-13    1.62921e-12           3.04        1.63e-12
    78        323     3.72145e-13    4.07302e-13            1.9        4.07e-13
    79        324     3.72145e-13    4.07302e-13            1.9        4.07e-13
    80        329     3.71954e-13    1.01826e-13          0.476        1.02e-13

No solution found.

fsolve stopped because the last step was ineffective. However, the vector of function
values is not near zero, as measured by the value of the function tolerance. 

X = 4×1

1.0e-04 * 0.2733 0.8154 0.4145 0.4965

fval = 4×1

1.0e-06 * -0.0036 -0.0000 -0.6099 0.0000

exitflag = -2

norm(fsolve_function(X))

ans = 6.0988e-07

X = X.^2

X = 4×1

1.0e-08 * 0.0747 0.6648 0.1718 0.2465

%equilibrium composition- function
function F=fsolve_function(X)
%unpack variables
E1=X(1)^2;
E2=X(2)^2;
E3=X(3)^2;
Noi=X(4)^2;%initial value of oxygen
K1=4.669E20;
K2=1.29;
K3=3.44E19;
yo2eq=0.00906;%atm
%input equations
eq1=(2*E1-E2)^2-K1*(Noi-E1)*(Noi+E1-E2+2*E3);
eq2=E2*(E2-2*E3)-K2*(2*E1-E2)*(2*E3-E2);
eq3=(2*E3-E2)^2*(Noi+E1-E2+2*E3)-K3*(E2-2*E3)^2*(Noi-E1);
eq4=(Noi-E1-E3)-yo2eq*(Noi+E1-E2+2*E3);
F=[eq1; eq2; eq3; eq4; ];
end