Getting error Solve::nsmet

Question

I'm trying to solve the following system of equations and I can't seem to manipulate Mathematica to be get an answer. Any and all help would be appreciated.

Some notes: XNil and XNis should both be between 0 and 1 when solved. The value of 1700 in the equation should be able to range from 1356 to 1726.

HmCu = 13055;
TmCu = 1356;
HmNi = 17472;
TmNi = 1726;

GmCu[T_] = HmCu - T*(HmCu/TmCu)
GmNi[T_] = HmNi - T*(HmNi/TmNi)

R = 8.3144;
Lol = 11100;
Los = 10435;

In[30]= Rationalize[{GmCu[1700] + R*1700*Log[(1 - XNil)/(1 - XNis)] + 
Lol*XNil^2 - Los*XNis^2 == 0 && 
GmNi[1700] + R*1700*Log[(XNil)/(XNis)] + Lol*(1 - XNil)^2 - 
Los*(1 - XNis)^2 == 0}]

Out[30]= {-(1122730/339) + 11100 XNil^2 - 10435 XNis^2 + 
353362/25 Log[(1 - XNil)/(1 - XNis)] == 0 && 
227136/863 + 11100 (1 - XNil)^2 - 10435 (1 - XNis)^2 + 
20786/25 Log[XNil/XNis] == 0}

Solve[Out[30], {XNil, XNis}, Reals]

Solve::nsmet: This system cannot be solved with the methods available to Solve.

EDIT: I was able to find a solution using FindRoot instead of Solve. However, both solutions are above 1 which does not make sense within the context of the problem.

Answer 1

Your problem arrises from the numerical difficulty of the solving your equations. To see the problem, let's plot the zero contours of the left hand side of your equations.

HmCu = 13055;
TmCu = 1356;
HmNi = 17472;
TmNi = 1726;
GmCu[T_] = HmCu - T*(HmCu/TmCu);
GmNi[T_] = HmNi - T*(HmNi/TmNi);
R = 8.3144;
Lol = 11100;
Los = 10435;

eqs = 
  Rationalize @
    {GmCu[1700] + R*1700*Log[(1 - XNil)/(1 - XNis)] + Lol*XNil^2 - Los*XNis^2 == 0, 
     GmNi[1700] + R*1700*Log[(XNil)/(XNis)] + Lol*(1 - XNil)^2 - Los*(1 - XNis)^2 == 0}

Quiet @ 
  ContourPlot[Evaluate @ eqs, {XNil, 0, 4}, {XNis, 0, 4}, 
    ContourStyle -> {Red, Blue}]

There appears to be a zero near the upper right corner of the unit square, but the two traces are so nearly parallel that FindRoot is not going to do well when it uses Newton's method (the default method) for this problem. We will tell it to use the secant method by giving it a bracket of starting values rather than a single starting point. Further, we want the bracketing to fit rather closely around the root, and it will be necessary to work with arbitrary precision arithmetic.

To get good starting values, we zoom the plot in on the desired zero.

Quiet @ 
  ContourPlot[Evaluate @ eqs, {XNil, 0, 4}, {XNis, 0, 4},
    PlotRange -> {{.89, .91}, {.91, .93}},
    ContourStyle -> {Red, Blue}]

Now we can use FindRoot and get a good result.

FindRoot[eqs, {{XNil, .89`20, .91`20}, {XNis, .91`20, .93`20}}, 
  WorkingPrecision -> 40]

{XNil -> 0.9010823440563772363797513599485915326375, 
 XNis -> 0.9208285280785948177856956468479419307813}

Answer 2

@m_goldberg's comment got me to wondering if I accidentally hit upon a good starting point with the following:

FindRoot[eqs, {{XNil, 0.5}, {XNis, 0.51}}]
(*  {XNil -> 0.901082, XNis -> 0.920829}  *)

It turns out the the OP's problem has an pretty infeasible region (red means the solution did not converge to the solution above):

With[{n = 1024, (* n = 256 is faster and looks ok; 128 takes 4+ sec. for me *)
  s0 = {XNil, XNis} /. FindRoot[eqs, {{XNil, 0.5}, {XNis, 0.51}}]},
 Graphics[{PointSize[0.0001], Red,
   Point[
     Flatten[
      ParallelTable[
       Block[{s = {XNil, XNis} /. Quiet@FindRoot[eqs, {{XNil, x}, {XNis, y}}]},
        If[EuclideanDistance[s, s0] > 0.001 // TrueQ, {x, y}, Nothing]
        ],
       {x, 1./n, 1. - 1./n, 1./n}, {y, 1./n, 1. - 1./n, 1./n}],
      1]
    ]
   },
  Frame -> True]
 ]