2
$\begingroup$

Consider the equation

eq = x == y - Log[1 + y];

I would like to solve this equation for y

yofx = y /. Solve[eq, y]

{-1 - ProductLog[-E^(-1 - x)]}

Unfortunately, this solution seems to be incorrect, as can be verified by inspecting the following two plots

{Plot[y - Log[1 + y], {y, 0, 10}], Plot[yofx, {x, 0, 10}]}

enter image description here

On the left, the plot correctly rises and asymptotes to a straight line for large y and x. On the right, the plot should have been a flipped image of the left plot, with x- and y-axes interchanged, but the curve looks nothing like it.

Is it possible to get the correct solution for this equation in Mathematica? If so, how can it be done?

Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ I get a warning from Solve that the answer may be wrong. Why not follow the advice in the message? $\endgroup$
    – Michael E2
    Apr 12 at 3:22
  • $\begingroup$ @MichaelE2 Using Reduce produces an output that still contains an unsolved equation for a constant specifying a branch cut. Instead of trying to solve that equation, I just found by trial and error that the correct branch constant value is -1, since the constant must be an integer. Thanks for the suggestion! $\endgroup$
    – Kagaratsch
    Apr 12 at 3:29
  • 1
    $\begingroup$ If you use Reduce[x == y - Log[1 + y] , y] you'll get a generalized solution. Since it seems you are interested in positive values of x and y, you could try Reduce[x == y - Log[1 + y] && x > 0 && y > 0, y], which narrows your solution to y==-1 - ProductLog[-1, -E^(-1 - x)]. If you then plot it using Plot[-1 - ProductLog[-1, -E^(-1 - x)], {x, 0, 10}], you'll get what looks like the axes-flipped version of your left plot. $\endgroup$
    – theorist
    Apr 12 at 3:30
  • $\begingroup$ @theorist Thank you, that's even better! $\endgroup$
    – Kagaratsch
    Apr 12 at 3:32
  • 2
    $\begingroup$ It's worth noting that the inverse function Mathematica provided is the inverse of $y - \ln(1 + y)$ restricted to the domain $(-1, 0]$ rather than $[0, \infty)$. So it's not wrong per se, it's just that Mathematica made an assumption different from yours. $\endgroup$
    – Michael Seifert
    Apr 12 at 15:49

2 Answers 2

Reset to default
4
$\begingroup$

To plot an equation use ContourPlot or, in this case, ParametricPlot

eq = x == y - Log[1 + y];

ContourPlot[Evaluate@eq, {x, 0, 10}, {y, 0, 10},
 FrameLabel -> (Style[#, 16] & /@ {x, y})]

enter image description here

ParametricPlot[{eq[[-1]], y}, {y, 0, 10},
 Frame -> True,
 PlotRange -> {{0, 10}, {0, 10}},
 PlotRangePadding -> Scaled[.02],
 FrameLabel -> (Style[#, 16] & /@ {x, y})]

(* same plot *)

Using your approach,

yofx = SolveValues[{eq, y > 0}, y, Reals][[1]] // Quiet

Plot[Evaluate@yofx, {x, 0, 10},
 Frame -> True,
 PlotRange -> {{0, 10}, {0, 10}},
 PlotRangePadding -> Scaled[.02],
 FrameLabel -> (Style[#, 16] & /@ {x, y}),
 AspectRatio -> 1]

(* same plot *)

Improve this answer
$\endgroup$
1
  • $\begingroup$ Adding a condition y > 0 and requirement Reals to Solve is a very nice way to get the result directly, without looking at Reduce! Thank you! $\endgroup$
    – Kagaratsch
    Apr 12 at 3:48
4
$\begingroup$

This is a typical case where conditional and branch cut check routines are not sophisticated enough to get result, although you give all neccesary conditions (x>0 and y>0 or y<0).

But Solve knows the underlying rules, as can be seen, when you insert numbers for variables (here EulerGamma for x). Of course you have to check validity of that trick, e.g. with a graphic or insert solution in eq as shown below.

eq = x == y - Log[1 + y];

ContourPlot[Evaluate@eq, {x, 0, 6}, {y, -3, 3}, 
   GridLines -> {{EulerGamma}, {-1 - 
    ProductLog[-E^(-1 - EulerGamma)], -1 - 
    ProductLog[-1, -E^(-1 - EulerGamma)]}}]

enter image description here

{sol1 = Solve[eq /. x -> EulerGamma, y, Reals], sol1 // N}

(*   {{{y -> -1 - ProductLog[-E^(-1 - EulerGamma)]}, {y -> -1 - 
     ProductLog[-1, -E^(-1 - EulerGamma)]}}, {{y -> -0.729215}, {y -> 
    1.48916}}}   *)

Get the two solutions

{y1[x_], y2[x_]} = y /. sol1 /. EulerGamma -> x

(*   {-1 - ProductLog[-E^(-1 - x)], -1 - ProductLog[-1, -E^(-1 - x)]}   *)

Test it applying Exp to both sides of equations

Map[Exp[#] &, {eq /. y -> y1[x], eq /. y -> y2[x]}, {2}] // 
 FullSimplify[#, x > 0] &

(*   {True, True}   *)

Edit

Most simple way to get solutions, is just apply Exp to both sides of equation.

neweq = Map[Exp[#] &, eq, {1}]

(*   E^x == E^y/(1 + y)   *)

Solve[neweq && x > 0, y, Reals]

(*   {{y -> ConditionalExpression[-1 - ProductLog[-E^(-1 - x)], 
    x > 0]}, {y -> 
   ConditionalExpression[-1 - ProductLog[-1, -E^(-1 - x)], x > 0]}}   *)
Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.