Hello Rohit,

I have a follow-up question. When I plot eq1 around one of the solutions (the negative one), it seems the equation is not actually defined. Would you know why that point would be a solution if it's not on one of the equations? Please see my code and plot below:

Thank you again!

\[Beta] = 0.95
R = 1/\[Beta]
W = 10^6
A = 100
\[Phi] = 0.2
\[Lambda] = 0.1
eq1 = p - \[Phi]/R - ((1 - \[Phi]) \[Phi] (1 - \[Lambda]) \[Beta])/(((W - T p) R + \[Phi] T) R + \[Phi]*(1 - \[Phi]) T)/(R (\[Lambda]/((W - T p) R + \[Phi] T  ) + ((1 - \[Lambda]) \[Beta] R)/(((W - T p) R + \[Phi] T) R + \[Phi] (1 - \[Phi]) T)))

sol = Solve[eq1 == 0 && p == A/T, {p, T}]

Which gives 2 solutions: {p -> -0.0000199018, T -> -5.02466(10^6)}, {p -> 0.319279, T -> 313.206}. And then I plot the eq1 for T between -5.025(10^6) and 0 to see what it looks like around the negative solution T= -5.02466(10^6). The plot shows that there is no corresponding value for p when T= -5.02466(10^6).

ContourPlot[
 p - \[Phi]/R - ((1 - \[Phi]) \[Phi] (1 - \[Lambda]) \[Beta])/(((W - T p) R + \[Phi] T) R + \[Phi]*(1 - \[Phi]) T)/(R (\[Lambda]/((W - T p) R + \[Phi] T  ) + ((1 - \[Lambda]) \[Beta] R)/(((W - T p) R + \[Phi] T) R + \[Phi] (1 - \[Phi]) T))) == 0, {T, -5.025*(10^6), 0}, {p, -0.1, 0.4}]

T is the x-axis, p is the y-axis

Hello Rohit,

Thanks for the tip on Clear[], not clear[]! Now I don't need to restart the kernel each time. My main issue is that one of the solutions ({p -> -1.98998(10^-8), T -> -5.02518(10^9)}) I got from solving the system of equations is not even a point on equation 1. That is, when I plug T = -5.02518(10^9) into equation 1 in my main post, I get {p -> 0.325402}, instead of {p -> -1.9899810^-8} . Would you by any chance know how this can be a solution even when it only satisfies one equation?

Thank you!