Need help getting the intercept of two lines in a ContourPlot

Question

I am trying to find a set of parameters $\{\epsilon,\epsilon_b\}$, so that the two local minima have exactly the difference of 15.1 (see code below) and the left sided local minimum and the local maximum have a distance of 2.5 (see code below). I am able to solve this problem in form of two ContourPlots and Show, but I am not able to find the intercept of the two ContourLines with NSolve. Maybe someone can help me with that and explain why I can solve this problem graphically, but not numerically.

The code:

U = Function[{r, ϵ, ϵb},E^(-((2(r-2σ)^2)/[Sigma]^2)) ϵb + 4 ϵ (-(σ^6/r^6) + σ^12/r^12) + (lc Log[(lc (1 - r^2/lc^2) Sinh[(r (3 - r^2/lc^2))/(lc (1 - r^2/lc^2))])/(r (3 - r^2/lc^2))])/(β σ)];

T = 298.15;
kB = 8.3145*10^-3;
β = 1/(kB T);
lc = 65;
σ = 2.5;

xfCoord[ϵ_?NumberQ, ϵb_?NumberQ] := FindMinimum[{U[r, ϵ, ϵb],0 < r < 2 σ}, {r, σ}]
xf[ϵ_?NumberQ, ϵb_?NumberQ] := r /. FindMinimum[{U[r, ϵ, ϵb],0 < r < 2 σ}, {r, σ}][[2]]
Uf[ϵ_?NumberQ, ϵb_?NumberQ] := FindMinimum[{U[r, ϵ, ϵb],0 < r < 2 σ}, {r, σ}][[1]]
xuCoord[ϵ_?NumberQ, ϵb_?NumberQ] :=FindMinimum[{U[r, ϵ, ϵb],3 σ < r < lc}, {r,3 σ}]
xu[ϵ_?NumberQ, ϵb_?NumberQ] :=r /. FindMinimum[{U[r, ϵ, ϵb],3 σ < r < lc}, {r, 3 σ}][[2]]
Uu[ϵ_?NumberQ, ϵb_?NumberQ] := FindMinimum[{U[r, ϵ, ϵb], 3 σ < r < lc}, {r,3 σ}][[1]]
xbCoord[ϵ_?NumberQ, ϵb_?NumberQ] :=FindMaximum[{U[r,ϵ, ϵb], σ < r <0.9 lc}, {r, 2 σ}]
xb[ϵ_?NumberQ, ϵb_?NumberQ] :=r /. FindMaximum[{U[r,ϵ, ϵb], σ < r <0.9 lc}, {r, 2 σ}][[2]]
Ub[ϵ_?NumberQ, ϵb_?NumberQ] :=FindMaximum[{U[r, ϵ, ϵb], σ < r <0.9 lc}, {r, 2 σ}][[1]]

Δx[ϵ_?NumberQ, ϵb_?NumberQ] := xb[ϵ, ϵb] - xf[ϵ, ϵb]
ΔG[ϵ_?NumberQ, ϵb_?NumberQ] := Uu[ϵ, ϵb] - Uf[ϵ, ϵb]

p1 = ContourPlot[Δx[ϵ, ϵb] == 2.5, {ϵ, 2, 20}, {ϵb, 1, 20},FrameLabel -> {"ϵ","ϵb"}];
p2 = ContourPlot[ΔG[ϵ, ϵb] == 15.1, {ϵ, 2, 20}, {ϵb, 1, 20},FrameLabel -> {"ϵ","ϵb"},ContourStyle -> Green];
Show[p1, p2]

This now does not work:

NSolve[{Δx[ϵ, ϵb] == 2.5, ΔG[ϵ, ϵb] == 15.1}, {ϵ, ϵb}]

It just gives back the same expression, meaning it can not evaluate it, I guess.

Answer 1

That's because NSolve[ ] isn't the right tool for the job as it deals primarily with linear and polynomial equations.

Try FindRoot[ ] instead:

FindRoot[{Δx[ϵ, ϵb] == 2.5, ΔG[ϵ, ϵb] == 15.1}, {{ϵ, 5}, {ϵb, 5}}]

(* {ϵ -> 14.6863, ϵb -> 5.28574} *)

Answer 2

U = Function[{r, ϵ, ϵb}, 
   E^(-((2 (r - 2 σ)^2)/σ^2)) ϵb + 
    4 ϵ (-(σ^6/r^6) + σ^12/
        r^12) + (lc Log[(lc (1 - 
             r^2/lc^2) Sinh[(r (3 - r^2/lc^2))/(lc (1 - r^2/lc^2))])/(r (3 - 
             r^2/lc^2))])/(β σ)];

T = 298.15;
kB = 8.3145*10^-3;
β = 1/(kB T);
lc = 65;
σ = 2.5;

To restrict a function's arguments to numeric values use NumericQ rather than NumberQ. NumberQ will not accept numeric symbols (e.g., π, E, GoldenRatio).

NumberQ /@ {π, E, GoldenRatio}

(*  {False, False, False}  *)

xfCoord[ϵ_?NumericQ, ϵb_?NumericQ] :=
  FindMinimum[
   {U[r, ϵ, ϵb], 0 < r < 2 σ},
   {r, σ}];
xf[ϵ_?NumericQ, ϵb_?NumericQ] :=
  r /. FindMinimum[{U[r, ϵ, ϵb], 
      0 < r < 2 σ}, {r, σ}][[2]];
Uf[ϵ_?NumericQ, ϵb_?NumericQ] :=
  FindMinimum[
    {U[r, ϵ, ϵb], 0 < r < 2 σ},
    {r, σ}][[1]];
xuCoord[ϵ_?NumericQ, ϵb_?NumericQ] :=
  FindMinimum[{U[r, ϵ, ϵb], 3 σ < r < lc},
   {r, 3 σ}];
xu[ϵ_?NumericQ, ϵb_?NumericQ] :=
  r /. FindMinimum[
     {U[r, ϵ, ϵb], 3 σ < r < lc},
     {r, 3 σ}][[2]];
Uu[ϵ_?NumericQ, ϵb_?NumericQ] :=
  FindMinimum[
    {U[r, ϵ, ϵb], 3 σ < r < lc},
    {r, 3 σ}][[1]];
xbCoord[ϵ_?NumericQ, ϵb_?NumericQ] :=
  FindMaximum[
   {U[r, ϵ, ϵb], σ < r < 0.9 lc},
   {r, 2 σ}];
xb[ϵ_?NumericQ, ϵb_?NumericQ] :=
  r /. FindMaximum[
     {U[r, ϵ, ϵb], σ < r < 0.9 lc},
     {r, 2 σ}][[2]];
Ub[ϵ_?NumericQ, ϵb_?NumericQ] :=
  FindMaximum[
    {U[r, ϵ, ϵb], σ < r < 0.9 lc},
    {r, 2 σ}][[1]];

Δx[ϵ_?NumericQ, ϵb_?NumericQ] :=
  xb[ϵ, ϵb] - xf[ϵ, ϵb];
ΔG[ϵ_?NumericQ, ϵb_?NumericQ] :=
  Uu[ϵ, ϵb] - Uf[ϵ, ϵb];

ContourPlot[
 {Δx[ϵ, ϵb] == 
   2.5, ΔG[ϵ, ϵb] == 15.1},
 {ϵ, 2, 20}, {ϵb, 1, 20},
 FrameLabel -> {"ϵ", "ϵb"}]

Use FindRoot. The initial estimates can be taken from the plot

FindRoot[{Δx[ϵ, ϵb] == 
   2.5, ΔG[ϵ, ϵb] == 15.1},
 {{ϵ, 15}, {ϵb, 5}}]

(*  {ϵ -> 14.6863, ϵb -> 5.28574}  *)