Differential Equations - Maple Help
For the best experience, we recommend viewing online help using Google Chrome or Microsoft Edge.

Online Help

All Products    Maple    MapleSim

Home : Support : Online Help : System : Information : Updates : Maple 13 : Differential Equations

Updates to Differential Equation (DE) Solvers in Maple 13



Exact Solutions

Numerical Solutions



The theme for exact, symbolic Ordinary Differential Equation solving in Maple 13 is the exploration of alternative methods for finding solutions to DE problems that are beyond the scope of standard methods or for which standard methods give overly complicated results. For example, important nonlinear ODE families, including some which do not have point symmetries, can be solved by constructing and then solving related linear ODEs of higher order. Similarly, by exploring symmetries in non-standard ways, the solution for various ODE initial value problems can be expressed in simple form even when the ODE general solution cannot be computed. These abilities are now built in to the standard DE solver (dsolve), and the techniques are also available separately through a variety of new tools.


For Partial Differential Equations, the functionality of most of the symmetry commands of PDEtools has been significantly expanded in many respects and four new commands were added. Among the novelties there is one new command for the fast computation of polynomial solutions to linear and nonlinear PDE systems, and the internal PDEtools Library, with 45 specialized routines for programming purposes, is now available.


For numerical solutions of ODEs, improvements have been made to event handling and to accuracy and error control for ODE and index-1 DAE problems.

Exact Solutions

Ordinary Differential Equations (ODEs)


The dsolve command can now solve varied problems formerly out of reach, as well as compute remarkably simpler solutions for various ODE families for which former solutions were too complicated.

New solutions for nonlinear ODEs, from 1st to 5th order: linearization by differentiation


In some cases, a nonlinear ODE that cannot be solved using standard methods can still be multiplied by an appropriate factor so that after differentiation (one or two times), the resulting expression contains a linear ODE of higher order as a factor. Solving this linear ODE permits expressing the solution to the original problem explicitly and in rather compact manner.



PDEtools:-declare(y(x), prime=x);

yxwill now be displayed asy

derivatives with respect toxof functions of one variable will now be displayed with '


The following shows a first order ODE nonlinear in y', as an arbitrary function G of two arguments, and its corresponding new solution.

ode[1] := G((2*y(x)-2+diff(y(x),x)*x)*x^3/(3*x^4+2*x^3),-x^3*(-y(x)+1+(x+1)*diff(y(x),x))/(3*x+2)) = 0;



sol[1] := dsolve(ode[1], implicit);



So, for any form of G in ode[1], the same G in the output above solves the problem. ODE solutions can be verified as usual with odetest.

odetest(sol[1], ode[1]);



Here are two nonlinear ODE examples, of third and fourth order, extracted from scientific publications related to applications, and their new solutions computed quickly and in compact manner using linearization by differentiation.

ode[2] := diff(y(x),x,x,x) = 1/8*(-8*x*diff(y(x),x,x)*diff(y(x),x) + 4*x^2*diff(y(x),x,x)^2 + diff(y(x),x)^2 + 12*y(x))/x^2/diff(y(x),x);






ode[3] := diff(y(x),x,x,x,x) = 1/2*(2*diff(y(x),x,x,x)*diff(y(x),x)-diff(y(x),x,x)^2+c^4*y(x)^2)/y(x);






New ODE-IVP solutions from the ODE symmetries


For various nonlinear ODE initial value problems (ODE-IVP), a rational solution exists, but the general ODE solution involves elliptic or special functions, making it difficult to compute the simpler form when adjusting the integration constants to match the initial conditions. In Maple 13, dsolve automatically explores the ODE symmetries to compute an ODE-IVP solution in its simpler rational form directly, even if the general solution is entirely out of reach.



ode[4] := diff(y(x), x, x) + 3*(diff(y(x), x))/x - 8*y(x)^3/x^6 = 0;



ic[4] := y(1) = -1/2, D(y)(1) = -1/2;



The general solution of ode[4] involves the inverse of an elliptic integral, represented by the elliptic function JacobiSN.




It is not easy to adjust the integration constants C1, C2 so that the preceding solution satisfies the conditions y(1) = -1/2, D(y)(1) = -1/2. However, by exploring the ODE symmetries,




dsolve can compute directly the following simple solution to the whole ODE initial value problem.

sol_ic[4] := dsolve([ode[4], ic[4]]);



You can use odetest to test solutions against the initial conditions.

odetest(sol_ic[4], [ode[4], ic[4]]);



ode[5] := diff(y(x), x,x,x) - diff(y(x), x,x)*y(x) + diff(y(x),x)^2 = 0;



ic[5] := y(0) = 1, D(y)(0) = 1/6;



For ode[5] alone dsolve obtains only a reduction of order, which contains an ODE of Abel type for _g(_f) whose solution is not known.




This result is the best one achievable at present but there is no useful way to adjust it to match the conditions y(0) = 1, D(y)(0) = 1/6. In Maple 13, however, dsolve first detects the symmetries of ode[5] and instead of just failing to compute the general solution, it explores these symmetries, solving the problem as a whole.

sol_ic[5] := dsolve([ode[5], ic[5]]);



map(odetest, [sol_ic[5]], [ode[5], ic[5]]);



New solutions for linear ODEs with non-rational coefficients


The existing routines in dsolve for this problem were optimized to solve more complicated problems and in simpler form, or using simpler special functions.



ode[6] := diff(y(x),x,x) + tan(x)*diff(y(x),x) + a^2*cos(x)^2*sin(x)^(2*n-2)*y(x) = 0;



For the ode[6], the solution in Maple 12 involves hypergeometric functions. In Maple 13, the solution is computed in terms of simpler Bessel functions.




New solutions for nonlinear ODEs exploring dynamical symmetries


Maple's dsolve started exploring dynamical symmetries in 1997. The related solving approach, however, is algebraically complicated. For Maple 13, corresponding routines were enhanced substantially to handle these kinds of problems, taking advantage of the properties of subfamilies of these symmetries (for example, self-invariant contact symmetries). In this way, new or concretely simpler solutions can now be computed.



ode[7] := diff(y(x), x,x) = diff(y(x), x)+(ln(diff(y(x), x))-x)*diff(y(x), x)^2/(y(x)-diff(y(x), x));



This ODE admits only one point symmetry, but in addition it has a computable self-invariant contact (dynamical) symmetry that involves y'.




Remarkably, the point symmetry (first one in the preceding output) is of no use to integrate the equation; it leads to an intractable first order ODE. The dynamical symmetry, on the other hand, results in the following new implicit solution, involving the Ei (exponential integral) composed with the LambertW function.

dsolve(ode[7], implicit);



New parametric solutions for nonlinear ODEs


For ODEs of all differential orders, dsolve can now compute parametric solutions directly from their symmetries. This is useful in various cases where the standard solution, of the form yx=..., involves the algebraic inversion of special functions or integrals, while the parametric solution is simpler.



ode[8] := diff(y(x),x,x) = y(x)+_F1(diff(y(x),x)-y(x));






This solution is the best that can be achieved but involves composite integrals where the inner one is algebraically solved with respect to the integration variable using RootOf. These solutions can frequently be re-expressed in parametric form, free of RootOf, which can now be computed directly by dsolve from the ODE symmetries.

dsolve(ode[8], parametric);



DEtools[IVPsol], a new command for specializing general solutions and a related new option in dsolve


The standard process of solving an ODE-IVP consists of first computing a solution that is more general, which depends on arbitrary constants, and then adjusting these constants to match the initial or boundary conditions. The initial more general solution, however, can frequently be written in different ways, which can either facilitate or make more complicated the posterior adjustment of the solution.


To take advantage of this situation, a new command, IVPsol, has been added to the DEtools package, with the purpose of adjusting an ODE solution to given initial or boundary conditions. Also, related new flexibility has been added to dsolve: it now accepts a given solution to be used directly as a departure point to match the initial/boundary conditions.



ode[9] := x*diff(y(x),x,x)+4*y(x)^2 = 0;



A general solution for ode[9] is not known. So, in principle, it would not be possible to solve an initial value problem involving it, such as in the following example with these initial conditions.

ic[9] := y(1) = -1/2, D(y)(1) = 1/2;



But suppose you know a particular solution for ode[9] that you might have computed by any other means, for instance:

sol[9] := y(x) = _C1/(2*x);



You can now pass this solution directly to dsolve so that it uses it to solve the problem instead of trying and failing to compute one itself.

dsolve([ode[9], ic[9]], usesolution = sol[9]);



Analogously, you can pass the initial conditions and the solution directly to DEtools[IVPsol].

DEtools[IVPsol]([ic[9]], sol[9]);



DEtools: extended functionality of parametricsol and particularsol


The DEtools commands particularsol and parametricsol can now solve a new range of higher order and nonlinear ODE problems by combining standard methods with alternative ways of using symmetries.



ode[10] := diff(y(x), x,x,x) - diff(y(x), x,x)*y(x) + diff(y(x),x)^2 = 0;



No general solution is known for ode[10]. A particular solution for it, however, depending on one arbitrary constant, is now computable using the following.




For ode[8], introduced in a previous paragraph, a parametric solution is now computable using the following.




This functionality also has been added as an option to dsolve; see previous paragraphs.

Partial Differential Equations (PDEs)


A large number of improvements and additions were made to the PDEtools package, setting a new benchmark for the state-of-the-art in symbolic computation and Partial Differential Equations.

The new commands in PDEtools


Four new commands were added to PDEtools: CharacteristicQInvariants, ConsistencyTest, PolynomialSolutions, and SymmetrySolutions. As their names suggest, the new commands can perform fast consistency tests, compute polynomial solutions for PDE systems, transform given solutions into other different solutions exploring the symmetries of the problem, and can compute invariants of differential equation systems using a simpler approach, exploring the CharacteristicQ function of the system's symmetries.






Polynomial solutions for PDE systems can now be computed using the new PolynomialSolutions command.

PDEtools:-declare((f, g)(x,y));

fx,ywill now be displayed asf

gx,ywill now be displayed asg


sys[1] := diff(f(x,y),x)*diff(g(x,y),x) + diff(f(x,y),y)*diff(g(x,y),y) + g(x,y)*(diff(f(x,y),x,x) + diff(f(x,y),y,y)) = -1;



psol[1] := PolynomialSolutions(sys[1]);



map(pdetest, [psol[1]], sys[1]);



New different solutions for PDE systems can now be computed from old known solutions exploring their symmetries; consider, for example, the wave equation in two dimensions:

pde[2] := diff(u(x, t), x,x) - diff(u(x, t), t,t) = 0;



Here is a solution separable by product using pdsolve.

sol[2] := pdsolve(pde[2], HINT = `*`, build);



Using a new feature, symmetries of polynomial type:

S[2] := Infinitesimals(pde[2], typeofsymmetry = polynomial)[-1..-1];



A new solution is obtained transforming sol[2] and exploring the symmetry S[2].

newsol[2] := SymmetrySolutions(sol[2], S[2]);



This new solution is verified using pdetest.

pdetest(newsol[2], pde[2]);



Significant functionality added to compute solutions using symmetries


Most of the symmetry commands of PDEtools have had their functionality significantly expanded, providing enough flexibility to get a solution, or one of the necessary form, even for tough problems. New features that stand out are:


The main solver of the symmetry approach, InvariantSolutions has been packed with options to facilitate obtaining solutions of the most varied forms;


All of InvariantSolutions, SimilaritySolutions, and Infinitesimals can now solve problems faster, exploring symmetries of polynomial type or with restricted dependency, making solutions easier to compute;


InvariantSolutions and SimilaritySolutions automatically specialize arbitrary functions entering the infinitesimals of PDE systems, including an option to avoid this specialization when desired;


The key command of the symmetry approach, DeterminingPDE, can now compute the determining system also for integrating factors and conserved currents.




In many cases, a particular solution for a PDE system, which can also be computed faster than a more general solution, is enough to solve the problem. The new features permit searching for these particular solutions in varied manners. Recall, for example, the wave equation in two dimensions introduced in a previous section:




The following specifies that a solution depending only on t would suffice.

InvariantSolutions(pde[2], numberofsolutions = 1, dependency = t);



This other input specifies to construct a solution using infinitesimals only depending on one variable, either x or t, which are easy to compute.

InvariantSolutions(pde[2], dependencyofinfinitesimals = {x, t});



With this input, you specify the use of infinitesimals of polynomial type, linear in the jet variables (almost the simplest ones that can be computed, requiring only inverting a matrix); note the display option to show the corresponding invariants and solution.

InvariantSolutions(pde[2], degreeofinfinitesimals = 1, display);















These invariants shown are a fingerprint of each solution; they can now be copied directly from the screen and pasted in the call to InvariantSolutions to reproduce precisely the corresponding solution.

InvariantSolutions(pde[2], invariants = [t/x, u(x,t)]);



The general form of Infinitesimals for pde[2] involves arbitrary functions _Fnt+x and _Fmxt.




It is not possible to compute a SimilarityTransformation or, for that matter, any symmetry related transformation with them because of the arbitrary functions. In order to use these infinitesimals, the arbitrary functions are now automatically specialized by all the symmetry commands (except for Infinitesimals), as you could see by interactively using the new specialize_Fn option of Infinitesimals.

S := Infinitesimals(pde[2], specialize_Fn, displayfunctionality = false);



In this way, all these commands can deal with the symmetries of pde[2] properly.

SymmetryTransformation(S[5], u(x, t), v(r, s));



Invariant solutions can then be computed automatically because this specialization of arbitrary functions is now automatic.




The internal PDEtools Library of routines is now available for programming purposes


The PDEtools internal library for manipulating and programming with differential equations is now available, starting in Maple 13 with 45 subroutines.




The PDEtools Library of routines is made available mainly for programming purposes. Some of these routines, however, are valuable for mathematical reasons. For example, the result previously obtained using the new PolynomialSolutions command for sys[1] is computed by first estimating upper bounds for the degree of the solution.





# Parameters
# * PDESYS - a PDE system; it can contain PDEs, ODEs, and also non-differential equations
# * F - a list of the unknown functions of PDESYS
# * increasefornonautonomous - optional, can be true (default) or false, to increase the upper bounds by 1 in the presence of nonautonomous equations in PDESYS
# Description: UpperBounds returns a sequence of sets of equations, where the lhs is a function of F and the rhs is an upper bound for the degree of a solution polynomial in the independent variables.
# Example
# > sys := diff(f(x,y),x)*diff(g(y),y) + diff(f(x,y),y)*diff(g(y),y) + f(x,y)*(diff(f(x,y),x,x) + diff(f(x,y),y,y)) = x*y:
# > PDEtools:-Library:-UpperBounds(sys, [f(x,y), g(y)]);
# {f(x, y) = 3, g(y) = 2}
# > PDEtools:-Library:-UpperBounds(sys, [f(x,y), g(y)], increasefornonautonomous = false);
# {f(x, y) = 2, g(y) = 1}
UpperBounds( PDESYS, F, X::list(name) := _GetIndepVars(F),
             { increasefornonautonomous::truefalse := true } )

PDEtools:-Library:-UpperBounds(sys[1], [f, g](x,y));



Relevant extension in dpolyform and casesplit


The casesplit and dpolyform commands can now compute with arbitrary functions of algebraic expressions (not symbol variables). In other words, given an expression for some unknowns (that is, non-rational in arbitrary functions whose arguments are non-rational in the independent variables), compute a differential polynomial representation (PDE system), solved by the given expression, where the PDE system is rational in the unknowns and its derivatives, and has rational coefficients.


Motivation: arbitrary functions of algebraic expressions (typically differential invariants) appear frequently when solving intermediary PDE problems that enter the formulation of applications or larger PDE problems. The new ability permits doing differential elimination with these objects and so uncoupling nonlinear PDE systems that involve them, extending in significant ways what can now be done with PDE systems and their solutions.




Consider the following expression for some gx,t that involves the arbitrary functions f__1t+x+f__2xt.

sol[3] := g(x,t) = _F1(x+t) + _F2(x-t);



What is the PDE satisfied by this expression? It can now be computed: it is the wave equation in two dimensions.

pde[3] := dpolyform(sol[3], no_Fn);



The ability to compute differential polynomial representations for some mapping of variables x,t,... is old, based on the knowledge of the mathematical properties of the mapping applied to the variables (for example, sint+x). The ability to compute these PDE representations with arbitrary mappings (for example, f__1t+x) that have no particular mathematical properties is new in Maple 13. Derivatives of arbitrary mappings evaluated at expressions can also be represented in rational PDE form.

sol[4] := g(x,t) = D(_F1)(ln(x*t) + _F2(1/exp(t)));



pde[4] := dpolyform(sol[4], no_Fn);



These results can be verified as usual using pdetest.

pdetest(sol[4], pde[4]);



The casesplit command also accepts a new optional argument, caseplot, so that in addition to the splitting into cases, it also produces a plot with a graphical representation of the splitting.

ode[11] := diff(y(t),t,t)^2 + 2*diff(y(t),t,t)*y(t)^3*diff(y(t),t) - 4*y(t)^2*diff(y(t),t)^3;



casesplit(ode[11], caseplot);

========= Pivots Legend =========







Extended functionality in pdetest regarding boundary conditions


The pdetest command can now test solutions for correctness regarding boundary conditions.




Consider the following PDE, boundary condition, and solution.

pde[3] := diff(u(x,t),t) = k*diff(diff(u(x,t),x),x)+Q;



bc[3] := u(0,t) = 2*exp(k*t)-1/k*Q;



sol[3] := u(x,t) = _C1^2*exp(x+k*t)-(_C1^2-2)*exp(-x+k*t)-1/2/k*Q*x^2+1/_C1^2/k*Q*(_C1^2-2)*x-1/k*Q;



You can test whether the sol[1] solves pde[1] using pdetest; the new functionality is that you can now test whether it solves the boundary condition bc[1].

pdetest(sol[3], [pde[3], bc[3]]);



The boundary conditions can involve derivatives:

bc[3.1] := D[1,1](u)(0,t) = 2*exp(k*t)-1/k*Q;



pdetest(sol[3], [pde[3], bc[3.1]]);



Numerical Solutions

There have been improvements in several aspects of numerical ODE solution for Maple 13.


The event handling features introduced in Maple 12 have been extended with the ability to control initial triggering of events; the addition of new event programming constructs (tobegin, toend, breakiteration, delayhalt); the addition of the side construct for round-off control; the ability to restrict a root-finding trigger to an increasing or decreasing region of the trigger expression; and the ability to interactively query when events have fired (eventfired). See dsolve,Events for more details.


There is improved robustness of core solvers with respect to accuracy and error control for pure ODE problems and index-1 DAE variables. This provides better detection of solution discontinuities (for more information see examples,dsolve_numeric_NewErrorControl).


The error control improvements mean that Maple now gives a correct answer to the following problem:





A gap in the explicit Runga-Kutta technique error control mechanism has been fixed, and Maple now gives a correct answer to the following problem:





Maple can also perform more sophisticated error detection:


Error, (in dsn2) cannot evaluate the solution further right of .99999999, probably a singularity


This returns an error because the ODE is singular at t=1.


There is improved performance for DAE problems having relatively few degrees of freedom using an improved projection strategy.


There is more effective optimization with optimize=true option, and the addition of the new compile option for use of the Compiler with hardware precision evalhf-capable problems using the default rkf45 and rosenbrock solvers.


There is also the addition of the interactive numfun query option to provide the number of evaluations of the ODE performed in computation of a solution value for an IVP problem (supported for all IVP solvers that support maxfun).

See Also

Index of New Maple 13 Features