# Why Can't DSolve Find a Solution for this ODE?

I wanted to find a basis for the set of solutions of the following ODE.

$$y^{''}+\frac{1}{x^2+1}y^{'}(x)+\left[-1-\frac{1}{x^2+1}\right]y(x)=0$$

But when I try to use DSolve as follows

ODE =
Derivative[2][y][x] + 1/(1 + x^2) Derivative[1][y][x] + (-1 - 1/(1 + x^2)) y[x] == 0
DSolve[ODE, y, x]


I get an endless evaluation!

As far as I know, there is not any analysis issues which may cause difficulties to find a solution. So I am wondering what is wrong?

Then I tried Maple and it just gave the solution as

• Does the Maple result count as a solution since it's given in terms of an integral? – QuantumDot Jul 9 '16 at 21:19
• @QuantumDot: I think it is better than nothing! ;) – Hosein Rahnama Jul 9 '16 at 22:01
• Getting stuck in a loop is something that can happen with any mathematical software when integration is involved, some software just automatically times out though. – Feyre Jul 9 '16 at 22:01
• The FullSimplify[] version of this equation in the positive $x$ domain is (2 + x^2) y[x] == y'[x] + (1 + x^2) y''[x] which is a Sturm‐Liouville equation. Potentially the methods explained here will help. – Young Jul 9 '16 at 22:55
• @Quantum, to be fair, Mathematica is not above giving integral solutions as well: DSolve[{y'[x] == Sin[Sin[x]], y[0] == 0}, y, x] – J. M.'s discontentment Jul 10 '16 at 14:34

When DSolve runs a long time, it's probably because either Integrate or Solve is chewing over a tough problem. Here's a way to make Integrate give up quicker, so that you can see if that is the problem.

Note/warning: The code returns Inactive[Integrate][...] instead of an unevaluated Integrate[...]. The latter seems preferable, but I was unable to figure out how to do it. Luckily, the inactive version seems to work.

Base code: The function withTimedIntegrate runs code with Integrate under a time constraint of tc seconds.

ClearAll[withTimedIntegrate];
SetAttributes[withTimedIntegrate, HoldFirst];
withTimedIntegrate[code_, tc_] := Module[{$in}, InternalInheritedBlock[{Integrate}, Unprotect[Integrate]; i : Integrate[___] /; ! TrueQ[$in] :=
Block[{$in = True}, TimeConstrained[i, tc, Inactivate[i, Integrate]] ]; Protect[Integrate]; code ] ];  OP's example: withTimedIntegrate[{dsol} = DSolve[ode == 0, y, x], 1]; // AbsoluteTiming dsol  OP's example with an IVP: withTimedIntegrate[ {dsolIVP} = DSolve[{ode == 0, y[0] == 1, y'[0] == 0}, y, x], 1] // AbsoluteTiming y[1] /. dsolIVP // Activate // N  Check: NDSolveValue[{ode == 0, y[0] == 1, y'[0] == 0}, y[1], {x, 0, 1}] (* 1.83793 *)  • (+1) Perfect illustration! :) – Hosein Rahnama Jul 13 '16 at 19:58 • A slight shortening: TimeConstrained[i, tc] /.$Aborted :> Inactivate[i, Integrate] can be replaced with TimeConstrained[i, tc, Inactivate[i, Integrate]]. – J. M.'s discontentment Aug 5 '17 at 7:58
• @J.M. Thanks! . – Michael E2 Aug 5 '17 at 11:05
• @ooo Re Q1: It's now called the Villegas-Gayley trick. It's used several times on the site, which you can search for. – Michael E2 Jan 30 '18 at 12:47
• Re Q2: 1. Match Integrate[___] /; ! TrueQ[$in]; $in is undefined (or should be), so TrueQ[$in] is False and ! TrueQ[$in] is True; so this definition matches any Integrate[..]. 2. Then Block is execute which sets $in = True. 3. Then TimeConstrained, which executes the original Integrate command (represented by i) under a time constraint. 4. Now $in is True, so this def. does not match and the regular Integrate is executed. 5 If Integrate finishes in time, the result is returned; if not, the inactivated Integrate command is returned. – Michael E2 Jan 30 '18 at 12:52

With

\$Version
(* 10.4.1 for Microsoft Windows (64-bit) (April 11, 2016) *)


the code in the question produces the desired answer, although slowly,

(* {{y -> Function[{x}, E^x*C[1] + E^x*C[2]*Integrate[E^(-ArcTan[K[1]] - 2*K[1]),
{K[1], 1, x}]]}} *)


Computation time, as measured by AbsoluteTiming, is about 40 minutes on my PC.

As is so often the case, DSolve performs much better when given some help. Begin with the substitution, y[x] -> Exp[x] z[x].

Unevaluated[D[y[x], {x, 2}] + D[y[x], {x, 1}]/(1 + x^2) -
y[x] (1 + 1/(1 + x^2))] /. y[x] -> Exp[x] z[x];
Simplify[% Exp[-x]] // Apart
(* ((3 + 2*x^2)*Derivative[1][z][x])/(1 + x^2) + Derivative[2][z][x] *)


One might think that DSolve could solve this greatly simplified ode in seconds, but in fact it takes 36 minutes! (Perhaps, DSolve is searching for a solution that does not involve Integrate.)

DSolve[% == 0, z[x], x]
(* {{z[x] -> C[2] + Integrate[E^(-ArcTan[K[1]] - 2*K[1])*C[1], {K[1], 1, x}]}} *)


The obvious substitution z'[x] -> w[x] finally allows DSolve to proceed quickly.

%% /. {z''[x] -> w'[x], z'[x] -> w[x]};
DSolve[% == 0, w[x], x]
(* {{w[x] -> E^(-2 x - ArcTan[x]) C[1]}} *)


Back substitution and an additional integration then yield the desired result.

• So it seems to be a system issue? I was using 10.4.0 for Microsoft Windows (32-bit) – Hosein Rahnama Jul 10 '16 at 8:15
• @H.R. I just added the computation time on my PC to the answer, almost 40 minutes. Did you wait that long? If so, it must be that the newest version is more capable. I do know that a lot of improvements were made to DSolve in version 10. – bbgodfrey Jul 10 '16 at 12:29
• I really appreciate your detailed answer and attention. Many thanks for that. It seems that DSolve needs some improvement in the cases where solution can be represented in an integral form! :) – Hosein Rahnama Jul 10 '16 at 15:31
• It seems that you just did the reduction of order method with the help of Mathematica for computations. – Hosein Rahnama Jul 10 '16 at 17:01
• @H.R. Yes, I wanted to obtain some insight into why DSolve` took so long. – bbgodfrey Jul 10 '16 at 17:06