# Why is first substituting for the known variable and then solving for an unknown variable, using Solve[], correct whereas doing the reverse is wrong?

Consider the following code:

(*Initializations*)
c=10;
f=Log[1+c]-c/(1+c)


I am trying to figure out whether I should first substitute for \[Psi] and then solve for x or do it the reverse way.

Method 1:

N[Solve[\[Psi] == 1/f Log[1+cx]/x /. \[Psi]->1, x][[1,1]] ]


gives the result:

x -> 2.06561


whereas Method 2:

N[Solve[\[Psi] == 1/f Log[1+cx]/x, x][[1,1]] /. \[Psi]->1 ]


gives the result:

x -> -0.174588-0.019561 I


When I plot the function:

Plot[{1, 1/f Log[1+cx]/x}, {x,0,3}]


I get:

So, the answer of Method 1 is correct. But why does first substituting for \[Psi] and then solving for x give the correct answer as opposed to first solving for x and then substituting for \[Psi]?

• Compare N[Solve[\[Psi] == 1/f Log[1 + c*x]/x /. \[Psi] -> 1, x, Reals]][[1]] with N[Solve[\[Psi] == 1/f Log[1 + c*x]/x, x, Reals] /. \[Psi] -> 1][[1]] Commented Aug 3, 2022 at 18:01
• Thanks for the response. The two expressions you mention give the correct answer. However, the real problem I am trying to solve involves finding the quantity $\frac{\mathrm{d}^2f(x)}{\mathrm{d}\psi^2}$ expressed as a function of $\psi$. In order to express x as a function of \[Psi], I would have to use Solve[\[Psi] == 1/f Log[1 + c*x]/x, x, Reals]. However, Solve[] is not able to do this if I use the Reals keyword. I cannot directly substitute the value of \[Psi] here because I need to do the differentiation with respect to \[Psi] first. Commented Aug 3, 2022 at 18:19

Mathematica will automatically assume that all quantities are complex, unless explicitly told otherwise. In your case, if you specify that you want the real solutions to your equation, the two approaches lead to the same result:

c = 10;
f = Log[1 + c] - c/(1 + c);

NSolve[ψ == 1/f Log[1 + c x]/x /. ψ -> 1, x, Reals]  (* Out: {{x -> 2.06561}} *)
NSolve[ψ == 1/f Log[1 + c x]/x, x, Reals] /. ψ -> 1  (* Out: {{x -> 2.06561}} *)


Note that I replaced N[Solve[...]] with a direct call to NSolve.

• Thanks for your reply. I posted a comment under my question saying that I may not be able to use the Reals keyword because in the real problem I am trying to solve, I need to find an expression for x as a function of ψ. Then, I need to differentiate this function. Only then can I substitute the value for ψ. However, Solve[] is not able to do the following: Solve[\[Psi] == 1/f Log[1 + c*x]/x, x, Reals]. However, if I don't use the Reals keyword, I get a few different solutions. But none of them evaluate to the correct answer when I substitute for ψ. Commented Aug 3, 2022 at 18:35

In comments, OP mentioned that the actual problem is to find the second derivative of the function $$x=f(\psi)$$, and solve cannot find an explicit analytical expression for that function. This, however, can be approached numerically.

We write a function that solves the equation of interest for each numerical value of $$\psi$$ provided:

ClearAll[sol]
sol[ψ_?NumericQ] := x /. FindRoot[ψ == 1/f Log[1 + c x]/x, {x, 1/ψ}]


We then:

1. Use the adaptive sampling capabilities of the Plot functions to get well-distributed values of the function over a range of interest
2. Extract the coordinates of those points using Cases
3. Construct an interpolating function int using Interpolate
int =
Interpolation@
First@
Cases[
Plot[sol[psi], {psi, 0, 2}, PlotRange ->All],
Line[p_] :> p, All
];


We can then numerically differentiate the interpolating function (int'') and plot the result, together with the original function. I decided to use a log plot because the values of the function vary over a very wide range:

ListLogPlot[
{int, int''},
PlotLegends -> {"f", "f'", "f''"}, Joined->True
]