# Obtaining both Solutions of Linear, Second-Order ODE using DSolve

Bug introduced in 10.4 or earlier and persisting through 11.3

Reported as CASE:3849226.

While attempting to provide a general answer for question 138321, I obtained a result,

Flatten@DSolve[{D[(1 - R/r) D[Φ[r], r], r] - ((n (n + 1))/r^2) Φ[r] == 0}, Φ[r], r,
Assumptions -> n ∈ Integers && n > 0]
{* {Φ[r] -> (r^2 C[1] Hypergeometric2F1[1 - n, 2 + n, 3, r/R])/R^2 +
C[2] MeijerG[{{}, {1 - n, 2 + n}}, {{0, 2}, {}}, r/R]} *)


which appears to be incorrect. For instance, plot the two ostensibly independent solutions for n == 2 (and x == r/R).

Plot[Evaluate[{x^2 Hypergeometric2F1[1 - n, 2 + n, 3, x],
-MeijerG[{{}, {1 - n, 2 + n}}, {{0, 2}, {}}, x]/3} /. n -> 2], {x, 0, 1.2},
PlotRange -> All, Exclusions -> None]


Two problems are evident. First, the two solutions are identical (up to a constant multiplier, here 1/3) for x < 1. Second, neither solution is singular at x == 1, even though it is not difficult to show that one solution should have a logarithmic singularity there. Incidentally,

FunctionExpand[MeijerG[{{}, {1 - n, 2 + n}}, {{0, 2}, {}}, r/R], Abs[r/R] < 1]
(* Hypergeometric2F1[-1 - n, n, 1, 1 - r/R] *)


which I am confident is proportional to r^2 Hypergeometric2F1[1 - n, 2 + n, 3, r/R])/R^2 for all positive integer n. although I have not searched the myriad Hypergeometric Function identities to prove it.

So, my questions are

1. Is this a bug?
2. Is there a DSolve work-around to obtain the second independent solution of the ODE?

Addendum: Correct n == 2 solution

Further insight can be gained by solving the ODE for a specific value of n, which DSolve can handle correctly. For instance, with n == 2,

s2 = With[{n = 2}, FullSimplify@Flatten@DSolveValue[{D[(1 - R/r) D[Φ[r], r], r] -
((l (l + 1))/r^2) Φ[r] == 0}, Φ[r], r]]
(* (3 r^2 (4 r - 3 R) R^5 C[1] - 8 R (-24 r^2 + 6 r R + R^2) C[2] +
48 r^2 (4 r - 3 R) C[2] (-Log[r] + Log[r - R]))/(12 R^5) *)


The C are arbitrary constants. For convenience in plotting, set them to C[1] -> - R^-3 and C[2] -> - R^2/(16 Pi).

s2s = FullSimplify[{-Coefficient[s2, C[1]]/R^3, -ReIm@Coefficient[s2, C[2]] R^2/(16 Pi)}
/. r -> x R, R > 0 && x > 0]
(* {((3 - 4*x)*x^2)/4, {(1 + 6*(1 - 4*x)*x + 6*(3 - 4*x)*x^2*Re[Log[(-1 + x)/x]])/(24*Pi),
Piecewise[{{((3 - 4*x)*x^2)/4, x < 1}}, 0]}} *)

Plot[Evaluate[s2s], {x, 0, 1.2}, PlotRange -> {-.4, .4}, Exclusions -> None]


Comparing this plot with the first plot in the question makes clear that the DSolve solution for arbitrary positive integer n is incomplete. (The first solution is blue, Re of the second solution orange, and Im of the second solution green.) Note that the green curve is not independent of the blue curve and could be eliminated by taking a linear combination of the first and second solutions for 0 < x < 1.

• I only got to look at this now. It seems to me that DSolve[] completely ignored your n ∈ Integers, since the two solutions are indeed independent for noninteger n, but collapse into a single one otherwise (reminiscent of the situation when solving Bessel's equation, where you need a special limiting operation to get $Y$ in addition to $J$). – J. M.'s torpor Oct 27 '17 at 2:29
• To add: it might be profitable to try to bring the DE you have to hypergeometric form, based on the form of the solution you have, so that the solutions pairs featured here might be applicable. – J. M.'s torpor Oct 27 '17 at 2:32

eqn = D[(1 - R/r) D[Φ[r], r],
r] - ((n (n + 1))/r^2) Φ[r] == 0;

assume = n ∈ Integers && n > 0;

sol = Assuming[assume, DSolve[eqn, Φ, r][[1]]]

(*  {Φ -> Function[{r},
(r^2*C[1]*Hypergeometric2F1[
1 - n, 2 + n, 3, r/R])/
R^2 + C[2]*MeijerG[
{{}, {1 - n, 2 + n}},
{{0, 2}, {}}, r/R]]}  *)


The solution is valid for Abs[r/R] < 1

(eqn /. sol) // FullSimplify[#, assume && Abs[r/R] < 1] &

(*  True  *)


Or for Abs[r/R] > 1

(eqn /. sol) // FullSimplify[#, assume && Abs[r/R] > 1] &

(*  True  *)


However, the general solution fails for Abs[r/R] == 1

(eqn /. sol) // FullSimplify[#, assume && Abs[r/R] == 1] &

(*  Indeterminate == 0  *)


EDIT: The components are separately solutions

(eqn /. (sol /. C[1] -> 0)) // FullSimplify[#, assume && Abs[r/R] < 1] &

(*  True  *)

(eqn /. (sol /. C[2] -> 0)) // FullSimplify[#, assume && Abs[r/R] < 1] &

(*  True  *)

(eqn /. (sol /. C[1] -> 0)) // FullSimplify[#, assume && Abs[r/R] > 1] &

(*  True  *)

(eqn /. (sol /. C[2] -> 0)) // FullSimplify[#, assume && Abs[r/R] > 1] &

(*  True  *)


The discontinuity comes from the second component when Abs[r/R] == 1

(eqn /. (sol /. C[1] -> 0)) // FullSimplify[#, assume && Abs[r/R] == 1] &

(*  Indeterminate == 0  *)


The first component remains valid for Abs[r/R] == 1

(eqn /. (sol /. C[2] -> 0)) // FullSimplify[#, assume && Abs[r/R] == 1] &

(*  True  *)


Plotting the two independent solutions

Manipulate[
Plot[Evaluate[
List @@ (Φ[r] /. sol1) /.
{C[1] -> c1, C[2] -> c2,
r -> x*R, n -> m}],
{x, -1.1, 1.1},
Exclusions -> {-1, 1},
PlotRange -> All,
AxesLabel -> {Style[ToString[x == r/R, TraditionalForm], 14, Bold],
None},
PlotLegends -> Placed["Expressions", Above]],
{{c1, 1}, -5, 5, Appearance -> "Labeled"},
{{c2, 1}, -5, 5, Appearance -> "Labeled"},
{{m, 2, "n"}, Range[5]}]


• This answer does not show that the two solutions are independent for Abs[x] < 1, and they are not, as I noted in my question. Neither is 0 an independent solution for Abs[x] > 1. – bbgodfrey Feb 23 '17 at 13:23
• @bbgodfrey - setting either C[1] or C[2] to zero will separate the solutions and remain valid solutions. – Bob Hanlon Feb 23 '17 at 16:46
• Please see the addendum to my question. The two solutions in your plot are not independent. – bbgodfrey Feb 23 '17 at 17:34