# Transforming implicit solutions of an ODE involving InverseFunction to an explicit form

When I DSolve a nonlinear differential equation:

q1[t_] := 3 (y'[t]/y[t])^2 - (1/(y[t]^3)) - (1/(y[t]^4)) - 1
q2[t_] := -2 y''[t]/y[t] - (y'[t]/y[t])^2 - (1/(3 y[t]^4)) + 1


By:

DSolve[{ q1[t] == q2[t]}, y, t]


I’d like to ask how to interpret these solutions in a comprehensive way. I mean how to get rid of InverseFunction and the intergration in terms of the symbols K and have a simplified solution? Should the integration and the InverseFunction to be solved again in a separate step?

The initial conditions are arbitrary for y[t] or y'[t].

Any help is appreciated.

• btw you have typo DSolve[{q1[t]==q2[t]},a,t] the a should be y. "how to get rid of InverseFunction and the intergration" Mathematica does this when it can't solve intermediate steps. The integral looks like some version of elliptic integral which it can't solve analytically. So it returns the result as shown. i.e. this is the best it can do. Commented Feb 2, 2022 at 10:22
• For arbitrary initial conditions you cannot do better. But if initial condition is such that $c_1$ is zero then there is some hope. Commented Feb 2, 2022 at 10:32
• Would have been interesting to first solve numerically mySol = NDSolveValue[{q1[t] == q2[t], y[1] == 1, y'[1] == 1}, y, {t, 1, 2}] then fit the inverse expression to that solution. That's a bit of a challenge though. :)
– josh
Commented Feb 3, 2022 at 15:22
• fyi, this is the answer Maple gives q1 := 3*(diff(y(t),t)/y(t))^2 - (1/(y(t)^3)) - (1/(y(t)^4)) - 1; q2 := -2*diff(y(t),t\$2)/y(t) - (diff(y(t),t)/y(t))^2 - (1/(3*y(t)^4)) + 1; dsolve(q1=q2,y(t)) screen shot !Mathematica graphics which is similar to Mathematica but the inverse function is implicit in the Maple answer since the limit of integration has y on it. Commented Feb 4, 2022 at 9:51
• @Dr.phy: In order to interpret these solutions in a comprehensive way as you requested, I suggest you solve it by hand (yes, the hard way). Let p=y' and make it exact by computing an integration factor (y^2). This way you can better understand the inverse construct returned by Mathematica.
– josh
Commented Feb 4, 2022 at 11:59

StreamPlot shows a graphical solution

temp={y'[t], y''[t]} /. Solve[q1[t] == q2[t], y''[t]][[1, 1]] /. {y[t] -> y, y'[t] -> ys}
StreamPlot[tmp, {y, 0.001, 2}, {ys, -2, 2}, FrameLabel -> {y[t], y'[t]}]


Fortunately the ode can be transformed to first order and solved, because the ode doesn't depend explicitely on time t.

ode = q1[t] - q2[t] /. {y[t] -> y, y'[t] -> ys[y],y''[t] -> ys'[y] ys[y] }


It's sufficient for now to consider initiel conditions ys[y0]==0 (see graphical solution StreamPlot)

sol = Values@DSolve[{ode == 0, ys[y0] == 0}, ys, y] // Flatten

Plot[Table[Evaluate[{sol[[1]][y], sol[[2]][y]}], {y0, 0, 2, .1}], {y,0, 3}, AxesLabel -> {y[t], y'[t]}]


• @UlrichMeumann. Hello, thanks for your answer. May I ask do you plot yand y’ or y’ and y’’ ?because in the definition of temp and after Solve[] I don't get why do you mention y’’. Also please let me ask why the StreamPlot requires that y’[0]=0? and can we use ParametricPlot instead? Commented Feb 2, 2022 at 13:17
• Excuse me, I'm trying to understand the code you have used cause I'm not an expert in MA. I just wonder how DSlove[l1[t]==l2[t],y, t] returned a complicated solution, while Solve returned a solution could be plotted. Also what dose [[1,1]] mean afterSolve? Commented Feb 2, 2022 at 13:34
• @Dr.phy StreamPlot shows a graphical solution of the ode. The arrow in a point {y[t],ys[t]} show the direction {y'[t],ys'[t]} of the unknown solution. Commented Feb 2, 2022 at 16:36
• From this way of solution, can we use ParametricPlot for {y[t],y’’[t]}? Commented Feb 2, 2022 at 17:17
• @UlrichNeumann You should not write "fortunately this ode can be solved" since this kind of implicit solution is the same as the output of the original direct approach in the question. This does not provide any deeper insight toward an explicit exact solution. Does it? Commented Feb 4, 2022 at 2:25

Mathematica cannot get rid of InverseFunction in this case at least up to version 13.

If $$s^2(t)$$ is a cubic or quartic polynomial in $$t$$ and $$r(s,t)$$ is a rational function of $$s$$ and $$t$$ containing at least one odd power of $$s$$ then $$\int r(s,t)\; dt\;$$ is called an elliptic integral. If $$s^2(t)$$ is higher than $$4$$-th order polynomial then $$\int r(s,t)\; dt\;$$ is called hyperelliptic integral (see e.g. Hyper-elliptic integral) as is the case here.

Inverse functions to elliptic integrals are elliptic functions (i.e. doubly periodic meromorphic functions in the complex plane like e.g. Jacobi elliptic function $$sn$$ or Weierstrass elliptic function $$\wp$$, in Mathematica JacobiSN and WeierstrassP respectively), while inverse functions to hyperelliptic integrals are hyperelliptic functions. Even in case of elliptic functions and integrals Mathematica cannot provide a seamless way of inversing elliptic integrals to elliptic functions and vice versa. Compare this post How to solve a nonlinear second order ODE with another answer to the same question. In order to provide an explicit solution we had to transform appropriately original equation by changing the dependent variable, while another answer based on an automatic approach yields an implicit solution. This problem has been discussed on this forum mere extensively here 1, see also e.g. 2 or 3

The theory of hyperelliptic functions is not developed comarably to the theory elliptic functions and it still expects for someone like Weirerstrass. Consequently there is a lack of appropriate chapters in e.g. Digital Library of Mathematical Functions and of course one should not expect that Mathematica will overtake mathematical theories, nonetheless one can guess that several hyperelliptic integrals and functions would be expressible in terms of Heun functions and its inverses.

Here's a substitution that helps DSolve get an explicit solution, although -- note well -- the solution is still in terms of the equivalent of the hyperelliptic integrals that @Artes notes cannot be solved at present.

y2u = u[t] == E^(-Sqrt[3] t) y[t]^3;
ueq = q1[t] == q2[t] /.
DSolveDSolveToPureFunction@First@Solve[y2u, y[t]] // Simplify
usol = DSolve[ueq, u, t];
(ysol = Flatten[Solve[#, y[t]] & /@ (y2u /. usol), 1]) // Short

(*
2 E^(Sqrt[3] t) (2 Sqrt[3] u'[t] + u''[t]) ==
3 + (2 E^(-(t/Sqrt[3]))) / u[t]^(1/3)

{{y[t] -> -E^((t/Sqrt[3])) (<<1>> / (C[1] <<1>> <<1>>) + <<1>>)^(1/3)},
<<4>>,
{y[t] -> -(-1)^(2/3) <<1>> <<1>>^<<1>>}}
*)


The six solutions (2 x three cube-roots) are nearly 2MB:

ysol // ByteCount
(*  1995056  *)


Here's a related substitution that yields a nice-looking answer (to me), but no closer to being solved than any other solution so far:

y2v = v[t] == y[t]^3;
veq = q1[t] - q2[t] // Together // Numerator;
veq = veq == 0 /. DSolveDSolveToPureFunction@First@Solve[y2v, y[t]] //
Simplify;
vsol = DSolve[veq, v, t] /. Verbatim[Solve][eq_, __] :> eq /.
K[1] :> \[Eta] /. First@Solve[y2v, v[t]]