# Solving system of non-linear equations

I'm trying to plot the concentration of chemical reaction products (PL,PL,LPL) vs. molar ratio of initial components concentration (r). All concentrations are positive, rational and could not exceed initial concentration of components. The following code works fine:

sol = Solve[{K1*P*L == PL, K2*P*L == LP, K3*PL*L == LPL,
P0 == P + PL + LP + LPL, r*P0 == L + PL + LP + 2*LPL}, {P, L, LP,
PL, LPL}][[1]];
val = LPL /. sol;
val1 = PL /. sol;
val2 = LP /. sol;

complex[rVal_, K1Val_, K2Val_, K3Val_, P0Val_, DifInd_: 0] :=
Re@With[{r = rVal, K1 = K1Val, K2 = K2Val, K3 = K3Val, P0 = P0Val,
índex = DifInd}, (Evaluate@D[val, {r, índex}])];
complex1[rVal_, K1Val_, K2Val_, K3Val_, P0Val_, DifInd_: 0] :=
Re@With[{r = rVal, K1 = K1Val, K2 = K2Val, K3 = K3Val, P0 = P0Val,
índex = DifInd}, (Evaluate@D[val1, {r, índex}])];
complex2[rVal_, K1Val_, K2Val_, K3Val_, P0Val_, DifInd_: 0] :=
Re@With[{r = rVal, K1 = K1Val, K2 = K2Val, K3 = K3Val, P0 = P0Val,
índex = DifInd}, (Evaluate@D[val2, {r, índex}])];

Manipulate[Plot[
Evaluate[{
complex[r, K1, K2, K3, P0, 0], rcomplex1[r, K1, K2, K3, P0, 0],
rcomplex2[r, K1, K2, K3, P0, 0]
}],
{r, 0.1, 40},
ImageSize -> 500, Frame -> True,
PlotStyle -> {Thick, {Red, Thick, Dashed}, Green}],
Style["Constants", 12, Bold],

{{K1, 100000, Subscript[Style["K", Italic], 1]}, 1000, 1000000, 10,
Appearance -> "Labeled", ImageSize -> Tiny},
{{K2, 100001, Subscript[Style["K", Italic], 2]}, 1000, 1000000, 10,
Appearance -> "Labeled", ImageSize -> Tiny},
{{K3, 100000, Subscript[Style["K", Italic], 3]}, 1000, 1000000, 10,
Appearance -> "Labeled", ImageSize -> Tiny},

Delimiter, Style["Concentration", 12, Bold],
{{P0, 0.00001, Subscript[Style["P", Italic], 0]}, 0.00001, 0.001,
0.00001, Appearance -> "Labeled", ImageSize -> Tiny},
ControlPlacement -> Right]


But for another reaction I could not plot the same curves. Any idea why it could happen?

sol = Solve[{K1*T*T == TT, K2*TT*S == TST, T0 == T + 2*TT + 2*TST,
r*T0 == S + TST}, {T, S, TT, TST}];

val = TT /. sol;
val1 = TST /. sol;

complex[rVal_, K1Val_, K2Val_, T0Val_, DifInd_: 0] :=
Re@With[{r = rVal, K1 = K1Val, K2 = K2Val, T0 = T0Val,
índex = DifInd}, (Evaluate@D[val, {r, índex}])];
complex1[rVal_, K1Val_, K2Val_, T0Val_, DifInd_: 0] :=
Re@With[{r = rVal, K1 = K1Val, K2 = K2Val, T0 = T0Val,
índex = DifInd}, (Evaluate@D[val1, {r, índex}])];

Manipulate[Plot[
Evaluate[{complex[r, K1, K2, T0, 0], complex1[r, K1, K2, T0, 0]}],
{r, 0.1, 40},
ImageSize -> 500, Frame -> True,
PlotStyle -> {Thick, {Red, Thick, Dashed}, Green}],

Style["Constants", 12, Bold],
{{K1, 1000, Subscript[Style["K", Italic], 1]}, 1000, 1000000, 10,
Appearance -> "Labeled", ImageSize -> Tiny},
{{K2, 100001, Subscript[Style["K", Italic], 2]}, 1000, 1000000, 10,
Appearance -> "Labeled", ImageSize -> Tiny},

Delimiter, Style["Concentration", 12, Bold],
{{T0, 0.00001, Subscript[Style["P", Italic], 0]}, 0.00001, 0.001,
0.00001, Appearance -> "Labeled", ImageSize -> Tiny},

ControlPlacement -> Right]


-
You have 4 solutions, two of which are complex. It looks like you are getting crossings of solution branches. Might try adding Cubics -> False and Quartics -> False to Solve, and using either of the first two solutions, to see if that improves things. – Daniel Lichtblau Nov 25 '12 at 20:14