[I'm using V13.0.1]
There are three parts to this answer:
- What
DSolve
does in this case.
- How to transform the Heun equation to the OP's ODE.
- How to get the general solution to the OP's ODE.
First answer: What DSolve
does
I can explain, literally, half of it:
ode = (1 - z) z Z''[z] - 2 z Derivative[1][Z][z] -
3/4 Z[z]/(z + a^2 (1 + z)) == 0;
coeff = Coefficient[-Z''[z] /. First@Solve[ode, Z''[z]], {Z'[z],
Z[z]}] // Simplify
(* {2/(-1 + z), 3/(4 (-1 + z) z (z + a^2 (1 + z)))} *)
DSolve`DSolveSecondOrderODE; (* autoload 2nd order fns *)
basis = DSolve`DSolveSecondOrderODEDump`DSolveSpecialOrder2Linear[
coeff, z]
(*
{Indeterminate,
z HeunG[-(a^2/(1 + a^2)),
2 - 3/(4 (1 + a^2)) + 2 (-1 - a^2/(1 + a^2)), 2, 1, 2, 2, z]}
*)
basis . {C[1], C[2]}
(* Indeterminate *)
ode /. Z ->
Function[z,
z HeunG[-(a^2/(1 + a^2)),
2 - 3/(4 (1 + a^2)) + 2 (-1 - a^2/(1 + a^2)), 2, 1, 2, 2, z]
] // FullSimplify
(* True *)
So what's going on? WRI improved DSolve
in V13 (see New in 13: Symbolic & Numeric Computation), and I believe that's what we're seeing. The particular solver finds half the basis for the solution space, but the other basis function comes out Indeterminate
. This usually happens if something like 0/0
happens, which could happen because things weren't simplified carefully. (If I have time, I'll see if I can track it down.) We can see from checking with FullSimplify
that the second basis function satisfies the ode
.
We can also see that the standard linear combination of the basis functions results in Indeterminate
(basis . {C[1], C[2]}
).
Consider reporting it as a bug. Since this is new functionality, they probably have not run into this edge case. It looks like it could be handled better.
How to transform the Heun equation to the OP's ODE.
@josh's interest prompted me to share the basis for one of my comments below. It is not much help to me, but maybe it will be for others. It shows only the relationship between the OP's equation and the Heun equation.
odeHeun = (-q + z α β) Z[z] +
((-1 + z) (-a + z) γ +
(-1 + z) z (1 + α + β - γ - δ) +
z (-a + z) δ) Z'[z] +
(-1 + z) z (-a + z) Z''[z] == 0;
(* transform Heun to OP's *)
xf = Block[{HeunG = List},
Thread[
HeunG[a, q, α, β, γ, δ, z] ->
HeunG[-(a^2/(1 + a^2)), -(3/(4 (1 + a^2))), 0, 1, 0, 2, z]]
]
(*
{a -> -(a^2/(1 + a^2)),
q -> -(3/(4 (1 + a^2))), α -> 0, β -> 1, γ -> 0, δ -> 2,
z -> z}
*)
(* Check xf transforms the Heun equation to the OP's *)
Solve[odeHeun /. xf, Z''[z]] // FullSimplify
Solve[ode, Z''[z]] // FullSimplify
% === %%
(* True *)
(* Extract the general basis and apply xf to get OP's case *)
First@DSolve[odeHeun, Z[z], z] /. (Thread[{C[1], C[2]} -> #] & /@
IdentityMatrix[2]) /. HeunG -> Inactive[HeunG]
% /. xf
% // Activate

We have γ = 0
in solution that evaluates to Indeterminate
. The series expansion of the general HeunG
function shows this and γ
a negative integer are exceptional cases (as implied in the DLMF):
Series[HeunG[a, q, α, β, γ, δ, z], {z, 0, 3}]

Possibly, the current state of knowledge means there is no way to express the complete solution, but I'm not an expert on Heun functions. Maybe the best place to seek an answer to the question whether there is a complete solution is on one of the other mathematics sites.
See below.
Update 2: How to get the general solution
At the singular point $z=0$, the OP's ODE has characteristic exponents (or indicial roots) $0$ and $1$.
We have one solution above $Z_1$ that solution corresponds to a series solution expanded about the singular point $z=0$ with characteristic exponent $1$. The difference in the characteristic exponents is an integer, so there is a Frobenius solution of the form
$$ Z_2 = \sum_{k=0}^\infty b_k z^k + a \log z \cdot Z_1 \,.$$
Whether this may be expressed in terms of higher transcendental functions is unclear to me. (I'd guess not, but I don't know.) However, AsymptoticDSolveValue
may be used to get a finite expansion of the Frobenius series. It takes a long time to get more than a few terms, though.
Since the ODE is linear, we may also try to solve it by reduction of order, setting $Z_2(z)=Y(z)Z_1(z)$. This produces a general solution in terms of an unevaluated integral (which is the sole basis for guessing that a solution in terms of higher transcendental functions is not possible).
(* z1sol from spelunking DSolve: *)
z1sol = Z ->
Function[z,
z HeunG[-(a^2/(1 + a^2)),
2 - 3/(4 (1 + a^2)) + 2 (-1 - a^2/(1 + a^2)), 2, 1, 2, 2, z]];
(* DSolve will do the last step in the reduction of order
* for us, so we need only make the substitution *)
odeY = odeHeun /. xf /. Z -> Function[z, Y[z] Z[z]] /. z1sol //
FullSimplify;
ysol = First@DSolve[{odeY,
(* irritating way to get the initial point for the
* integral different than 1, which is a singular point:
*)
dummy'[z] == 0, dummy[1/2] == 0},
{Y, dummy}, z];
gensol = Z -> Function @@ {z,
(Y[z] + C[2] /. ysol) (Z[z] /. z1sol) // FullSimplify}
(* You can omit the dummy variable and edit the integral
* with the following substitution Rule: *)
(*bettergensol = gensol //. Inactive[Integrate][f_,{x_,1,b_}] :>
RuleCondition[Inactive[Integrate][f,{x,1/2,b}],True]*)
(* Check: *)
odeHeun /. xf /. gensol // FullSimplify

On the dummy
variable.
Adding the dummy equation dummy'[z] == 0
makes the system one dimension larger. This can have consequences, sometimes ones you would wish to avoid like DSolve
thinking the system is too complicated to try to solve (a problem in previous version, which I assume is still present). However, starting the integral at a point where HeunG
is undefined (or Indeterminate
in this case) is rather unhelpful. The replacement Rule
method for bettergensol
is probably the better way to handle the integral. The dummy
method was to show that there is a way to get DSolve
to do it for you. There does not seem to be any other way to tell DSolve
where to start the integral.
Summary/synopsis
This is a bug. The fact that gensol
is easily obtained by standard methods shows that DSolve
could be improved in its handling of the Heun equation. Even if that were not possible, the Indeterminate
result for one of the basis functions should be handled better.
One may obtain the Frobenius series with order = 3; AsymptoticDSolveValue[ode, Z, {z, 0, order}]
, if that is desired.
HeunG
is the solution to the Heun equation corresponding to the characteristic exponent 0
at the singular point z == 0
for which HeunG[..., 0] = 1
. It is denoted $H\!\ell$ in DLMF; the solution in z1sol
corresponds to eq. (31.3.5). (I'm not an expert, so let me know if I made any blunders.)
HeunG[-(1/2), -(11/8), 2, 1, 2, 2, 2.]
is another bug, which I will report to WRI. It is obtained by setting a -> 1
, z -> 2.
in the HeunG
factor in z1sol
. It happens for most settings of a
as well. Note 0
, 1
, ∞
and the first argument of HeunG
are singular points of the Heun equation. When the first argument is equal to 0
or 1
(confluent case) or z
is equal to one of singular points, sometimes HeunG
does not evaluate or returns Indeterminate
.