Revision history [back]

Okay, I'll suppose this is not homework (i it is, please take the time to understand the solutions : it might be worth it...). It turns out exhibiting some interesting features (and non-features) of Sage and some packages interfaced to it, so let's dig it.

Sage (i. e. Maxima)

The first attempt:

var('x')
y = function('y')(x)
de = diff(y, x) == (4*x + y)^2
desolve(de, y, x)

fails with a suggestion:

    NotImplementedError: Maxima was unable to solve this
    ODE. Consider to set option contrib_ode to True.

which, when applied, gives a curious result:

Sol=desolve(de, y, ivar=x, contrib_ode=True)
Sol

[[x == _C - 1/2*arctan(1/2*sqrt(t)), y(x) == -4*x - sqrt(t)
 [x == _C + 1/2*arctan(1/2*sqrt(t)), y(x) == -4*x + sqrt(t)]]

This parametric form introduces non only an integration constant _C but also a parameter t ; both x and y are functions of t, which is only used through its square root. We can use this to solve sqrt(t) in x (hint : the justification of this solution is not trivial…) and get two solutions…

(_C, t)=var("_C, t")
SolSage=map(lambda u:u[1].subs(u[0].solve(sqrt(t))), Sol)
SolSage

… which turn out to be the same, and can be checked:

bool(SolSage[0].rhs().diff(x)==de.subs(SolSage[0]).rhs())

True

One notes that the constant _C is determined by the initial conditions but modulo $\pi$.

Sympy

By default, Sympy returns a Taylor series of the solution, which is not very interesting here. We'll ask for /all/ solutions of our equations:

D=sympy.dsolve(*map(sympy.sympify, [de, y(x)]), hint="all")
D

{'1st_power_series': Eq(y(x), x**3*(C1**2*(3*C1**2 + 4) + 12*C1**2 + 16)/3 + x**5*(8*C1**2*(3*C1**2 + 8) + 4*C1**2*(9*C1**2 + 8) + C1**2*(2*C1**2*(3*C1**2 + 8) + C1**2*(9*C1**2 + 8) + 36*C1**2 + 32) + 144*C1**2 + 128)/15 + C1 + C1*x**2*(C1**2 + 4) + C1*x**4*(C1**2*(3*C1**2 + 8) + 12*C1**2 + 32)/3 + C1**2*x + O(x**6)),
 'best': Eq(y(x), x**3*(C1**2*(3*C1**2 + 4) + 12*C1**2 + 16)/3 + x**5*(8*C1**2*(3*C1**2 + 8) + 4*C1**2*(9*C1**2 + 8) + C1**2*(2*C1**2*(3*C1**2 + 8) + C1**2*(9*C1**2 + 8) + 36*C1**2 + 32) + 144*C1**2 + 128)/15 + C1 + C1*x**2*(C1**2 + 4) + C1*x**4*(C1**2*(3*C1**2 + 8) + 12*C1**2 + 32)/3 + C1**2*x + O(x**6)),
 'best_hint': '1st_power_series',
  'default': '1st_power_series',
  'lie_group': Eq(y(x), -4*x + 2*tan(C1 + 2*x)),
  'order': 1}

The `lie_group solution can be expressed in Sage as:

D['lie_group']._sage_()

y(x) == -4*x + 2*tan(C1 + 2*x)

which is identical to the solution given by Sage up to the expression of the constant.

Giac

Unfortunately, Sage is not (yet) able to translate the differential equation in giac. But a "manual" translation allows us to call giac's solver, which gives:

from giacpy_sage import * SolGiac=libgiac.desolve(libgiac.diff(y(x),x)-de.rhs(),y(x))[0] SolGiac

2*tan(2*x-2*c_0)-4*x

i. e. again the same result up to the expression of the constant…

Fricas

I have been unable to use Fricas solver to solve this equation. But that's probably because I do not (yet) know the `fricas' interface in detail.

Mathematica

Again, Sage is currently unable to translate the differential equation. A manual call:

mathematica("DSolve[Derivative[1][y][x] == (4*x + y[x])^2, y[x], x]")

{{y[x] -> -2*I - 4*x + (-I/4 + E^((4*I)*x)*C[1])^(-1)}}

gives a result, manually translatable to Sage as -4*x + 4/(4*C1*e^(4*I*x) - I) - 2*I, totally different from those obtained so far.

A heavy pounding on this expression by deMoivre and algebraic transformations may reveal that this expresses the same relation between x and y as before. I haven't done this.