[[_high_order, _with_linear_symmetries]]
Book solution method
TO DO
Mathematica ✓
cpu = 0.839093 (sec), leaf count = 310
Maple ✓
cpu = 0.215 (sec), leaf count = 148
DSolve[-(a^4*x^3*y[x]) - x*y''[x] + 2*x^2*y'''[x] + x^3*y''''[x] == 0,y[x],x]
Mathematica raw output
{{y[x] -> C[1]*HypergeometricPFQ[{}, {3/4, 5/8 - Sqrt[5]/8, 5/8 + Sqrt[5]/8}, (a
^4*x^4)/256] + 2^(-3 - Sqrt[5])*((1 + I)*2^(1/2 + Sqrt[5])*a*x*C[2]*Hypergeometr
icPFQ[{}, {5/4, 7/8 - Sqrt[5]/8, 7/8 + Sqrt[5]/8}, (a^4*x^4)/256] + (-1)^((3 - S
qrt[5])/8)*a^(3/2 - Sqrt[5]/2)*x^(3/2 - Sqrt[5]/2)*(4^Sqrt[5]*C[3]*Hypergeometri
cPFQ[{}, {1 - Sqrt[5]/4, 9/8 - Sqrt[5]/8, 11/8 - Sqrt[5]/8}, (a^4*x^4)/256] + (-
1)^(Sqrt[5]/4)*a^Sqrt[5]*x^Sqrt[5]*C[4]*HypergeometricPFQ[{}, {9/8 + Sqrt[5]/8,
11/8 + Sqrt[5]/8, 1 + Sqrt[5]/4}, (a^4*x^4)/256]))}}
Maple raw input
dsolve(x^3*diff(diff(diff(diff(y(x),x),x),x),x)+2*x^2*diff(diff(diff(y(x),x),x),x)-x*diff(diff(y(x),x),x)-a^4*x^3*y(x) = 0, y(x),'implicit')
Maple raw output
y(x) = _C1*hypergeom([],[3/4, 5/8+1/8*5^(1/2), 5/8-1/8*5^(1/2)],1/256*a^4*x^4)+_
C2*x*hypergeom([],[5/4, 7/8+1/8*5^(1/2), 7/8-1/8*5^(1/2)],1/256*a^4*x^4)+_C3*x^(
3/2-1/2*5^(1/2))*hypergeom([],[1-1/4*5^(1/2), 9/8-1/8*5^(1/2), 11/8-1/8*5^(1/2)]
,1/256*a^4*x^4)+_C4*x^(3/2+1/2*5^(1/2))*hypergeom([],[1+1/4*5^(1/2), 9/8+1/8*5^(
1/2), 11/8+1/8*5^(1/2)],1/256*a^4*x^4)