Symbolic differential equation solving deserves better syntax

[Jason Merrill]

I thought I would take my first shot at contributing to Sage by making
desolve return a SymbolicExpression instead of a string. Doing so
turned out to be pretty easy, but I got hung up trying to get the
other methods in desolvers.py to do the same. Anyway, after getting
my hands in this code a little bit, I think it needs big changes--so
much so that incremental improvements don't seem so worthwhile.

In particular, I don't like that desolve_laplace is a different
function from desolve is a different function from desolve_system, and
I don't like the rather cryptic way of specifying boundary conditions.

I guess Mathematica is the leader on solving differential equations
symbolically, and pending other great ideas, I think their syntax is
worth copying. Here's an example of the DSolve syntax in Mathematica:

DSolve[{y''[x] + x^2 y[x] == 0 , y[0] == 0, y'[0] == 1}, y, x]

The arguments are a list of equations, the dependent variable, which
can also be a list if there is more than one dependant variable, and
finally the independent variable (or variables for PDEs). What I
really like about this is that boundary conditions are specified as
equations, and not as a list of symbols separated from their meaning.
I also like that there is only one function and one syntax, regardless
of what order the equation is, or whether 0 or more boundary values
are given, and whether the boundary values specify a Boundary Value
Problem or an Initial Value Problem.

This would make more work for Sage behind the scenes, since it has to
pick the equation apart and figure out what routine to send it off to
based on what kind of boundary values are given, and what the order of
the equation is, etc. But the point of a CAS is to do the algebra so
I don't have to!

I mostly hope to provoke discussion about syntax, rather than
implementation, but I'll mention that one (of several) stumbling
blocks to implementing the above is that it's not easy in Sage to ask
of a differential equation: "What order are you." This is in the same
vein as W5 at http://wiki.sagemath.org/symbench, and I hope it will
get easier soon with some more pattern matching type abilities.

Regards,

JM

[Tim Lahey]

On Aug 21, 2008, at 8:52 PM, Jason Merrill wrote:

>
> I guess Mathematica is the leader on solving differential equations
> symbolically, and pending other great ideas, I think their syntax is
> worth copying. Here's an example of the DSolve syntax in Mathematica:
>
> DSolve[{y''[x] + x^2 y[x] == 0 , y[0] == 0, y'[0] == 1}, y, x]
>
> The arguments are a list of equations, the dependent variable, which
> can also be a list if there is more than one dependant variable, and
> finally the independent variable (or variables for PDEs). What I
> really like about this is that boundary conditions are specified as
> equations, and not as a list of symbols separated from their meaning.
> I also like that there is only one function and one syntax, regardless
> of what order the equation is, or whether 0 or more boundary values
> are given, and whether the boundary values specify a Boundary Value
> Problem or an Initial Value Problem.

I'd recommend that the inputs be four lists. The first is the
differential equations, the second are the boundary/initial conditions
and the third are the dependent variables, and the fourth the
independent
variables. That way, the conditions are cleanly separated from the
equations which may be useful in the solution process.

Cheers,

Tim.

[Alec Mihailovs]

> I guess Mathematica is the leader on solving differential equations
> symbolically, and pending other great ideas, I think their syntax is
> worth copying. Here's an example of the DSolve syntax in Mathematica:

I think, Maple is better at that, especially for partial differential
equations. In particular, Lie symmetries and Heun functions are explored
much more in Maple than in Mathematica. It still has some bugs, but they
usually come from other sources, such as extremely buggy int, solve, or
simplify.

Alec

[Jason Merrill]

That sounds good too, as long as boundary conditions are input in the
form of equations rather than grunts. I like it a little less in the
case that you don't want to supply any boundary conditions--then you'd
have to supply an empty list to avoid resorting to keywords, right?
On the other hand, it might make it a little easier if you wanted to
try the same equation with a bunch of different sets of boundary
conditions or something like that. And as for "useful in the solution
process", I don't mind if the parser has to work hard on my behalf ;-)

JM

[Tim Lahey]

On Aug 21, 2008, at 10:22 PM, Jason Merrill wrote:
> That sounds good too, as long as boundary conditions are input in the
> form of equations rather than grunts. I like it a little less in the
> case that you don't want to supply any boundary conditions--then you'd
> have to supply an empty list to avoid resorting to keywords, right?
> On the other hand, it might make it a little easier if you wanted to
> try the same equation with a bunch of different sets of boundary
> conditions or something like that. And as for "useful in the solution
> process", I don't mind if the parser has to work hard on my behalf ;-)
>

Yes, I've found that since you're often getting the boundary conditions
and the equations separately, it's easier to have them as two separate
lists. If I remember correctly, this is what Maple does. Yes, if there
are no boundary conditions, you'd have to provide an empty list. I
don't think that's really much of a problem.

Cheers,

Tim.

[Jason Merrill]

I've never used Maple, but the Wikipedia example shows equation and
boundary conditions in a single list:
http://en.wikipedia.org/wiki/Maple_(software)#Solution_of_linear_differential_equations,
which is the same as Mathematica and my earlier suggestion. The thing
I don't like about an empty list is that it doesn't really represent
anything mathematical, and things like that which only serve to please
a parser are harder to learn and remember. If you are given the
equations and boundary conditions as separate lists, it's pretty easy
to combine them:

dsolve(eqns+bcs,y,x)

though the concatenation operator being + is unfortunate in a
mathematical context. Also, it's annoying if one of eqns or bcs is
actually a single element and not a list.

JM

[Alec Mihailovs]

Yes, Maple puts both ODE and initial conditions in one set, as

dsolve({ODE, ICs}, y(x), options)

Alec

[Martin]

> I guess Mathematica is the leader on solving differential equations
> symbolically, and pending other great ideas, I think their syntax is
> worth copying. Here's an example of the DSolve syntax in Mathematica:
>
> DSolve[{y''[x] + x^2 y[x] == 0 , y[0] == 0, y'[0] == 1}, y, x]

FriCAS / Axiom is supposed to be very good at linear differential
equations and differential equations of the form y'=f(x, y) - the code
is by Manuel Bronstein. It seems to be rather weak for others, it
cannot solve the equation above for example. I must admit, however,
that I do not know much about it. Furthermore, it does not attempt to
solve for the function, as the first example shows (sorry about the
ASCII-art, do you prefer LaTeX on the list?)

What concerns the syntax, it's as follows:

-- first you tell FriCAS that y is a function ("operator" in FriCAS
language)
y := operator 'y
-- then you input solve(eq, name of function, point to develop around,
initial values)
solve(D(y x, x, 2) + x^2* y x = 0, y, x=0, [0,1])

-- of course, initial values are optional
solve(x*(1-x^2)*D(y x, x) + (2*x^2 -1)*y x - x^3*(y x)^3 = 0, y, x)

5 2 4 2
(- 2x + 4)y(x) + 5x - 5x
----------------------------
2
5y(x)

Here is an example (taken from http://www.loria.fr/~zimmerma/ComputerAlgebra/ode_comp.ps.gz)
that maxima seems to be unable to do - at least for me in SAGE 3.1.
(The same report shows, that maxima is much better for nonlinear
equations)

solve(D(y x, x, 4)-4/x^2*D(y x, x, 2) + 8/x^3*D(y x, x)-8/x^4*y x, y,
x)

[particular= 0,
5 3 2 5 3 2 5 3 2
5 3 2
x - 5x + 5x - 1 x + 5x - 10x + 4 x - 10x + 5x + 4 x
- 10x + 20x + 4
basis=
[------------------,-------------------,-------------------,--------------------]]
x x
x x


Martin

[root]

>FriCAS / Axiom is supposed to be very good at linear differential
>equations and differential equations of the form y'=f(x, y) - the code
>is by Manuel Bronstein. It seems to be rather weak for others, it
>cannot solve the equation above for example. I must admit, however,
>that I do not know much about it. Furthermore, it does not attempt to
>solve for the function, as the first example shows (sorry about the
>ASCII-art, do you prefer LaTeX on the list?)

The Axiom project has created a test suite based on the Kamke
(E. Kamke "Differential Equations, Methods of Solutions and
Solutions" Leipzig (1956)) which, I believe, is the same test
suite used by Maple.

Portions of the test suite, including answers where Axiom can produce them,
is available at
<http://axiom.svn.sourceforge.net/viewvc/axiom/trunk/axiom/src/input/kamke0.input.pamphlet>
(also kamke1, kamke2, kamke3, kamke4, kamke5, kamke6, and kamke7)

Maple clearly dominates this area and has, by far, the best
implementation of an ODE solver. Like all of the math projects,
they dominate because they have the best person for the job.

Tim