Updates to Differential Equation Solvers in Maple 8

Differential equations solving capabilities have been greatly enhanced in Maple 8 regarding exact and numeric DE solutions.

•

New implementation of integrated factor methods added to dsolve

•

Symmetry algorithms refined

•

New routines for solving two parameterized Abel ODE classes added, further enhancing dsolve

•

Improvements made in computation of hypergeometric solutions

•

Algorithms for higher order ODEs enhanced

•

Maple can numerically solve timebased linear PDEs over rectangular domains.

•

Maple numerical ODE solution capabilities include numerous additions and improvements


Also in Maple 8, the DEtools and PDEtools packages contain enhancements and new tool commands.


Exact Solutions: Ordinary Differential Equations


•

The Maple dsolve command now includes a new implementation of integrating factor methods for systematically determining (when they exist) integrating factors of varied types for Ordinary Differential Equations. The integrating factor method is comparable with the symmetry method: both consist of transforming the nonlinear ODE problem into a overdetermined linear PDE system, the former searching for integrating factors  the cornerstone for reductions of order  and the latter searching for symmetry generators, with which reductions of order through changes of variables can be performed. The integrating factor approach implemented in Maple 8 can also perform multiple reductions of order in one step when more than one integrating factor is found. For nonlinear ODE problems, this combination of symmetry and integrating factor methods makes Maple's ODE solver unique among its peers.


The integrating factor approach leads to special function solutions for various nonlinear ODE families.

>

PDEtools[declare](y(x), prime=x);

 (1) 
>

ode := diff(y(x), x$3) = (1/y(x)*diff(y(x),x)2/x)*diff(y(x),x,x)+1/y(x)/x*diff(y(x),x)^2+1/x^3*y(x);

 (2) 
 (3) 
 (4) 

An example where the linearization process implicit in this approach is more transparent:

>

ode := diff(y(x),x,x,x) = 13/8/x^3*y(x)  1/8*(20*diff(y(x),x,x)*x^2 + 26*diff(y(x),x)*x)/x^3  1/8*(24*diff(y(x),x)*x^3*diff(y(x),x,x)  20*diff(y(x),x)^2*x^2)/x^3/y(x);

 (5) 

This nonlinear ODE example has three integrating factors, which are determined solving a related linear ODE. The integrating factors are then used simultaneously by the new Maple 8 ODE routines to reduce the ODE order by three in one step. To follow the solving process, you can turn ON the userinfo facility for dsolve.

>

infolevel[dsolve] := 4;

 (6) 
Methods for third order ODEs:
 Trying classification methods 
trying 3rd order ODE linearizable_by_differentiation
differential order: 3; trying a linearization to 4th order
trying differential order: 3; missing variables
trying differential order: 3; exact nonlinear
trying 3rd order, integrating factor of the form mu(y'') for some mu
Trying the formal computation of integrating factors depending on any 2 of [x, y, y', y'']
*** Sublevel 2 ***
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
> Calling odsolve with the ODE, 4*diff(diff(diff(y(x),x),x),x)*x^3+7*y(x)+10*diff(diff(y(x),x),x)*x^27*diff(y(x),x)*x, y(x)
*** Sublevel 2 ***
Methods for third order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
< LODE of Euler type successful
Successful computation of 3 integrating factors: [y(x)*x, x^(1/4+1/4*29^(1/2))*y(x), x^(1/41/4*29^(1/2))*y(x)]
Attempting computing related first integrals...
< high order, formal computation of integrating factors successful
 
 (7) 
•

For second order nonlinear ODEs, the symmetry algorithms introduced in Maple 7 have been refined, now leading to simpler solutions in varied cases.


The following ODE is from Kamke's book. Two nontrivial symmetries for it are found by dsolve routines leading to a solution in terms of Jacobi elliptic functions (note the use of the new dsolve option useint).

>

PDEtools[declare](f(x));

 (8) 
>

ode := 2*diff(y(x),x,x)*y(x)  diff(y(x),x)^2  6*diff(f(x),x,x)/diff(f(x),x)*y(x)*diff(y(x),x) + 2*(2*diff(f(x),x,x,x)/diff(f(x),x)+6*diff(f(x),x,x)^2/diff(f(x),x)^2)*y(x)^28*y(x)^3 = 0;

 (9) 
>

sol := dsolve(ode, y(x), useint);

Methods for second order ODEs: Trying classification methods 
trying 2nd order Liouville
trying 2nd order WeierstrassP
trying 2nd order JacobiSN
differential order: 2; trying a linearization to 3rd order
trying 2nd order ODE linearizable_by_differentiation
trying 2nd order, 2 integrating factors of the form mu(x,y)
trying differential order: 2; missing variables
> trying 2nd order, dynamical_symmetries, fully reducible to Abel through one integrating factor of the form G(x,y)/(1+H(x,y)*y')^2
 trying a change of variables {x > y(x), y(x) > x} and reentering methods for dynamical symmetries 
> trying 2nd order, dynamical_symmetries, fully reducible to Abel through one integrating factor of the form G(x,y)/(1+H(x,y)*y')^2
trying 2nd order, integrating factors of the form mu(x,y)/(y')^n, only the singular cases
trying symmetries linear in x and y(x)
trying differential order: 2; exact nonlinear
trying 2nd order, integrating factor of the form mu(y')
trying 2nd order, integrating factor of the form mu(x,y')
trying 2nd order, integrating factor of the form mu(x,y)/(y')^n, only the general case
trying 2nd order, integrating factor of the form mu(y,y')
trying differential order: 2; mu polynomial in y'
trying 2nd order, integrating factor of the form mu(x,y)
differential order: 2; looking for linear symmetries
trying 2nd order, integrating factor of the form mu(x,y)/(y')^n, only the general case
> trying 2nd order, dynamical_symmetries, only a reduction of order through one integrating factor of the form G(x,y)/(1+H(x,y)*y')^2
 trying a change of variables {x > y(x), y(x) > x} and reentering methods for dynamical symmetries 
> trying 2nd order, dynamical_symmetries, only a reduction of order through one integrating factor of the form G(x,y)/(1+H(x,y)*y')^2
 Trying Lie symmetry methods, 2nd order 
> Computing symmetries using: way = 3
> Computing symmetries using: way = 5
 
< successful computation of symmetries.2nd order, symmetry integrators to be used are: [can, can2, gon]
2nd order, trying integration using method can
> Computing canonical coordinates for the symmetry [1/diff(f(x),x), 2/diff(f(x),x)^2*diff(diff(f(x),x),x)*y]
> Calling odsolve with the ODE, diff(y(x),x) = 2*diff(diff(f(x),x),x)/diff(f(x),x)*y(x), y(x)
*** Sublevel 2 ***
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
> Rewriting ODE in canonical coordinates by means of differential invariants
Try computing 1 more symmetries for ODE written in canonical coordinates
> Computing symmetries using: way = 3
Found another symmetry:
 
Computing a convenient ordering to use the 2 symmetries availableTry integration with the canonical coordinates of the symmetry [1/diff(f(x),x), 2/diff(f(x),x)^2*diff(diff(f(x),x),x)*y]
> Calling odsolve with the ODE, diff(_b(_a),_a) = 1/2*_b(_a)*(1+8*_a^3*_b(_a)^2)/_a, _b(_a), explicit, HINT = [[_a, 3/2*_b]]
*** Sublevel 2 ***
symmetry methods on request
1st order, trying reduction of order with given symmetries:
 
1st order, trying the canonical coordinates of the invariance group> Computing canonical coordinates for the symmetry [_a, 3/2*_b]
> Calling odsolve with the ODE, diff(y(x),x) = 3/2*y(x)/x, y(x)
*** Sublevel 3 ***
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
> The canonical coordinates may not have unique inverse.
Trying gauging the symmetry to the form [0, eta(x,y)]
> Computing canonical coordinates for the symmetry [0, 4*_b^3*_a^3_b]
< 1st order, canonical coordinates successful
< differential order: 2; canonical coordinates successful
 
 (10) 
>

map(odetest, {sol}, ode);

 (11) 
•

New routines for solving two parameterized Abel ODE classes were added. This further enhances the toplevel performance of Maple 7 dsolve in tackling ODEs of this type, that is, of the form . The two integrable classes are based on recent work.


New Abel solvable class with solutions expressed using hypergeometric 0F1 (Bessel) functions:

>

infolevel[dsolve] := 0:

>

ode := diff(y(x),x) = 1/2*(x^2x1)/(x+1)^2/x^2*y(x)^31/2*(3*x+2)/x^2/(x+1)*y(x)^2;

 (12) 
>

sol := dsolve(ode,implicit);

 (13) 
 (14) 

Another new Abel solvable class, with solutions expressed using hypergeometric 1F1 (Whittaker) functions:

>

ode := diff(y(x),x) = (2*x+1)/(2*x^21)^2*y(x)^3+(1+2*x+6*x^2)/(2*x^21)^2*y(x)^2;

 (15) 
>

sol := dsolve(ode,implicit);

 (16) 
 (17) 
•

In addition, for Abel type ODEs, all routines for solving equivalence problems with respect to any of the 27 classes known by Maple dsolve can now handle fractional or abstract powers in a systematic manner (see the new Maple 8 DEtools[power_equivalent] command).


The following Abel ODE is from the suite found in Polyanin Zaitsev Handbook on ODEs (see dsolve,references).

>

ode := diff(y(x),x) = 1+(12*x^(7/2)+a)/x^(5/2)/y(x);

 (18) 
>

sol := dsolve(ode, implicit, useInt);

 (19) 
 (20) 
•

For second order linear ODEs there have been radical improvements in the computation of hypergeometric solutions. Concretely, Maple 8 dsolve can resolve the equivalence between a given ODE and one having 2F1, 1F1 or 0F1 hypergeometric solutions, whenever that equivalence involves power composed with Mobius transformations (see dsolve,hypergeometric). These algorithms new in Maple 8 largely supercede the related ones found in previous releases.


2F1 hypergeometric type of solution

>

infolevel[dsolve] := 4;

 (21) 
>

ode := diff(y(x),x,x) = 3/16*(x^2+10)*(4*x^6+16*x^4+23*x^210)/(x^2+2)^2/(x^2+5)^2/x^2*y(x);

 (22) 
Methods for second order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying a Liouvillian solution using Kovacic's algorithm
< No Liouvillian solutions exists
> Trying a solution in terms of special functions:
> Bessel
> elliptic
> Legendre
> Whittaker
> hyper3: Equivalence to 1F1 under 'a power @ Moebius'
> hypergeometric
> heuristic approach
> hyper3: Equivalence to 2F1, 1F1 or 0F1 under 'a power @ Moebius'
< hyper3 successful: received ODE is equivalent to the 2F1 ODE
pFq equation: _a*_b*y(x)+((_b+_a+1)*x_c)*diff(y(x),x)+(x^2x)*diff(y(x),x,x) = 0
Transformation resolving the equivalence: {x = 3/2*t^2/(t^2+5), y(x) = 1/t^(1/4)/(t^2+5)^(1/8)/(t^2+2)^(1/2)*u(t)}
Values of hypergeometric parameters: [_a = 1/4, _b = 1/2], [_c = 3/4]
< hypergeometric successful
< special function solution successful
 
 (23) 
 (24) 

1F1 hypergeometric type of solution

>

ode := diff(y(x),x,x) = 27*x*y(x)/(2*x+1)/(x1)^4;

 (25) 
Methods for second order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying a Liouvillian solution using Kovacic's algorithm
< No Liouvillian solutions exists
> Trying a solution in terms of special functions:
> Bessel
> elliptic
> Legendre
> Whittaker
> hyper3: Equivalence to 1F1 under 'a power @ Moebius'
< hyper3 successful: received ODE is equivalent to the 1F1 ODE
< Whittaker successful
< special function solution successful
 
 (26) 
 (27) 
>

convert(sol, StandardFunctions);

 (28) 

0F1 hypergeometric type of solution

>

ode := diff(y(x),x,x) = 1/20/x^2*(405*x^65670*x^4+58604*x^2+13720)/(3*x^214)^3*y(x);

 (29) 
Methods for second order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying a Liouvillian solution using Kovacic's algorithm
< No Liouvillian solutions exists
> Trying a solution in terms of special functions:
> Bessel
> elliptic
> Legendre
> Whittaker
> hyper3: Equivalence to 1F1 under 'a power @ Moebius'
> hypergeometric
> heuristic approach
> hyper3: Equivalence to 2F1, 1F1 or 0F1 under 'a power @ Moebius'
< hyper3 successful: received ODE is equivalent to the 0F1 ODE
pFq equation: y(x)x*diff(y(x),x,x)diff(y(x),x)*_c = 0
Transformation resolving the equivalence: {x = 14/5*t^2/(3*t^214), y(x) = 1/(3*t^314*t)^(1/2)*u(t)}
Values of hypergeometric parameters: [], [_c = 1]
< hypergeometric successful
< special function solution successful
 
 (30) 
 (31) 
 (32) 

An example with fractional power coefficients, whose solution can be expressed in elementary form (see new in Maple 8 convert,elementary):

>

ode := diff(y(x),x,x) = 1/4*(2/3*1/(x^(5/6))+2/x^(2/3))/x*y(x) + 1/2*(1+x^(1/6))/x*diff(y(x),x);

 (33) 
Methods for second order ODEs:
 Trying classification methods 
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying an equivalence, under noninteger power transformations,
to LODEs admitting Liouvillian solutions.
> Trying a Liouvillian solution using Kovacic's algorithm
A Liouvillian solution exists
Reducible group (found an exponential solution)
Reducible group (found another exponential solution)
< Kovacic's algorithm successful
< Equivalence, under noninteger power transformations successful
 
 (34) 
 (35) 
>

convert(sol, elementary);

 (36) 
•

The Maple implementation of Kovacic's algorithm has been extended to handle fractional and abstract powers.

>

ode := diff(y(x),x,x) = 1/2*(x+a*x^(1/2))/x^2*diff(y(x),x)+1/4*b*(b*x+a)/x^2*y(x);

 (37) 
Methods for second order ODEs:
 Trying classification methods 
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying an equivalence, under noninteger power transformations,
to LODEs admitting Liouvillian solutions.
> Trying a Liouvillian solution using Kovacic's algorithm
A Liouvillian solution exists
Reducible group (found an exponential solution)
Group is reducible, not completely reducible
Solution has integrals. Trying a special function solution free of integrals...
> Trying a solution in terms of special functions:
> Bessel
> elliptic
> Legendre
> Whittaker
> hyper3: Equivalence to 1F1 under 'a power @ Moebius'
> hypergeometric
> heuristic approach
> hyper3: Equivalence to 2F1, 1F1 or 0F1 under 'a power @ Moebius'
> trying extended approach, based on solving the ODE satisfied by the "derivative" of y(x)
> ODE succesfully computed; trying to solve: diff(diff(y(x),x),x) = 4*a*(2*b*x^23*x+a)/x/(4*b^2*x^4+4*b*x^2*a4*a*x+a^2)*diff(y(x),x)1/4*y(x)*(4*b^2*x^44*b*x^2*a+4*a*xa^2)/x^4
< extended approach, failed in solving the ODE for the "derivative" of y(x)
> Mathieu
> Equivalence to the rational form of Mathieu ODE under 'a power @ Moebius'
> trying extended approach, based on solving the ODE satisfied by the "derivative" of y(x)
> ODE succesfully computed; trying to solve: diff(diff(y(x),x),x) = 4*a*(2*b*x^23*x+a)/x/(4*b^2*x^4+4*b*x^2*a4*a*x+a^2)*diff(y(x),x)1/4*y(x)*(4*b^2*x^44*b*x^2*a+4*a*xa^2)/x^4
< extended approach, failed in solving the ODE for the "derivative" of y(x)
> Heun: Equivalence to the GHE or one of its 4 confluent cases under 'a power @ Moebius'
< Heun successful: received ODE is equivalent to the "HeunD" ODE, case "c = 0"
Heun equation: diff(y(x),x,x)(2*x^5+4*x^3+x^4*_a2*x_a)/(x1)^3/(x+1)^3*diff(y(x),x)(_b*x^2+(_c2*_a)*x_d)/(x1)^3/(x+1)^3*y(x)
Transformation resolving the equivalence: {x = (2^(1/2)*(b*a)^(1/2)*ta)/(2^(1/2)*(b*a)^(1/2)*t+a), y(x) = 1/exp(1/2*(2^(1/2)*(b*a)^(1/2)*ta)*(2^(1/2)*(b*a)^(1/2)*t+a)/a/t)/t^(1/2)*u(t)}
Values of Heun parameters: [_a = 4*2^(1/2)*(b*a)^(1/2), _b = 14*2^(1/2)*(b*a)^(1/2), _c = 8*2^(1/2)*(b*a)^(1/2), _d = 4*2^(1/2)*(b*a)^(1/2)+1]
< Kovacic's algorithm successful
< Equivalence, under noninteger power transformations successful
 
 (38) 
•

The algorithms for higher order linear ODEs have been enhanced, incorporating three new algorithms, all working in the presence of nonrational objects like arbitrary functions or fractional powers.

>

ode := diff(y(x),x$3) = (2*diff(y(x),x,x) + 2*y(x))/(2*diff(f(x),x)^2  2*diff(f(x),x,x) + 2)*diff(f(x),x,x,x) + diff(y(x),x) + ((2*diff(y(x),x,x) + 2*y(x))*diff(f(x),x)^3 + ((6*diff(f(x),x,x) + 2)*diff(y(x),x,x) + 6*y(x)*diff(f(x),x,x) 2*y(x))*diff(f(x),x))/(2*diff(f(x),x)^22*diff(f(x),x,x) + 2);

 (39) 

For arbitrary values of , the following ODE has exponential solutions with simple dependence, from where a general solution can be constructed:

>

sol := dsolve(ode, y(x));

Methods for third order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying high order exact linear fully integrable
trying to convert to a linear ODE with constant coefficients
trying differential order: 3; missing the dependent variable
> pFq: Equivalence to the 3F2 or one of its 3 confluent cases under 'a power @ Moebius'
trying a solution in terms of MeijerG functions
> Try computing a Rational Normal Form for the given ODE...
< unable to resolve the Equivalence to a Rational Normal Form
checking substitution methods
checking for doubly periodic methods
trying reduction of order using simple exponentials
> Calling odsolve with the ODE, diff(_b(_a),_a) = (diff(diff(diff(f(_a),_a),_a),_a)+diff(f(_a),_a)^3+3*diff(f(_a),_a)*diff(diff(f(_a),_a),_a)diff(f(_a),_a)+diff(f(_a),_a)^2+diff(diff(f(_a),_a),_a)1)*_b(_a)/(diff(f(_a),_a)^2+diff(diff(f(_a),_a),_a)1), _b(_a)
*** Sublevel 2 ***
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
< reduction of order using simple exponentials successful
 
 (40) 
 (41) 

In the presence of fractional or abstract power, dsolve now attempts a Rational Normal Form to remove the nonrational elements (see the new in Maple 8 DEtools[power_equivalent] command). In the displayed userinfo, you can see how the varied methods compose between themselves to result in a special function solution:

>

infolevel[dsolve] := 4;

 (42) 
>

ode := diff(diff(diff(y(x),x),x),x) = 1/5*(4*x*n+10*x+54*n)/x/(x1)*diff(diff(y(x),x),x)  1/25*(9+4*a*(x1)^2*x^(4/5) + 4*n*(x1)*(3+13*x) + x*(68+91*x))/(x1)^2/x^2*diff(y(x),x) + 1/125*(4*a*(x1)^2*(4+9*x) + 5*x^(1/5)*(9+4*n*(x1)*(3+13*x) + x*(68+91*x)))/(x1)^3/x^(11/5)*y(x);

 (43) 
>

sol := dsolve(ode, y(x));

Methods for third order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying high order exact linear fully integrable
trying to convert to a linear ODE with constant coefficients
trying differential order: 3; missing the dependent variable
> pFq: Equivalence to the 3F2 or one of its 3 confluent cases under 'a power @ Moebius'
trying a solution in terms of MeijerG functions
> Try computing a Rational Normal Form for the given ODE...
RNF transformation found: {t = x^(4/5), u(t) = 1/x^(1/5)*y(x)*exp(Int(1/15*(4*x*n+10*x+54*n)/x/(x1),x))}
> Calling odsolve with the ODE, diff(diff(diff(u(t),t),t),t) = ((1/12*a*n+7/24*a)/t^2+(5/18*n+13/544/9*n^2+2/27*n^3)/t^3)*u(t)+(1/4*a/t+(1/3*n+1/3*n^2+1/12)/t^2)*diff(u(t),t), u(t)
*** Sublevel 2 ***
Methods for third order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying high order exact linear fully integrable
trying to convert to a linear ODE with constant coefficients
trying differential order: 3; missing the dependent variable
Equation is the LCLM of (1/6*(2*n1)/t1/(2*(2*n1)/a+t))*y(t)+diff(y(t),t), 1/36*(4*n8*n^2+4+9*a*t)/t^2*y(t)1/6*(2*n1)/t*diff(y(t),t)+diff(diff(y(t),t),t)
trying differential order: 1; missing the dependent variable
checking if the LODE is of Euler type
> Attemtping a differential factorization
trying exponential solutions
checking if the LODE is of Euler type
1, exponential solutions found
exponential solutions successful
< differential factorization successful
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying a Liouvillian solution using Kovacic's algorithm
< No Liouvillian solutions exists
> Trying a solution in terms of special functions:
> Bessel
< Bessel successful
< special function solution successful
< solving the LCLM ode successful
< solving the Rational Normal Form successful
 
 (44) 
 (45) 


Enhancements in the DEtools and PDEtools Packages



Maple 8 includes enhancements and new tool commands for both linear and nonlinear differential equations

•

A new command, DEtools[power_equivalent], receives a linear ODE of arbitrary order  or a nonlinear first order ODE of Abel type  with nonrational coefficients involving fractional or abstract powers. It returns an equivalent  linear or Abel  ODE, with rational coefficients, followed by the transformation mapping it into the given nonrational ODE. The returned ODE is always in normal form and the coefficients' degrees in the independent variable are the lowest which could be achieved through power transformations.


A nonlinear Abel type ODE example:

>

PDEtools[declare](y(x), prime=x);

 (46) 
>

ode := diff(y(x),x) = 1+(12*x^(7/2)+a)/x^(5/2)/y(x);

 (47) 
>

DEtools[power_equivalent]( ode );

 (48) 
>

ode := diff(y(x),x,x) = 1/2*(x+a*x^(1/2))/x^2*diff(y(x),x)+1/4*b*(b*x+a)/x^2*y(x);

 (49) 
>

DEtools[power_equivalent]( ode );

 (50) 

DEtools[power_equivalent] works by resolving an Equivalence problem between the received ODE and a minimal rational one under power transformations composed with transformations of the dependent variable; that is, of the form , where a is an arbitrary constant (possibly an algebraic expression) and P(x), Q(x) are arbitrary functions  for linear ODEs . This command is useful when tackling some typical type of nonrational ODEs using techniques for rational ones and is in use by dsolve in Maple 8 to solve varied new ODE families not solved in previous Maple releases.

•

The DEtools[intfactor] command can now systematically compute  when they exist  integrating factors for nonlinear ODEs (all orders) where the integrating factor can depend on as many variables {x, y, y',...} as the differential order of DE minus one. The technique explores Maple's ability to tackle PDE systems and performing differential elimination when the system contains nonpolynomial coefficients (see dpolyform).

>

ode := diff(y(x),x,x,x) = (diff(y(x),x)*x+y(x))*diff(y(x),x,x)/(y(x)*x);

 (51) 
>

DEtools[intfactor]( ode );

Methods for first order ODEs: Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
> Calling odsolve with the ODE, diff(y(x),x)*x+y(x), y(x)
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
 
 (52) 

This new capability for computing integrating factors is in use by dsolve in Maple 8 when solving nonlinear ODEs. Both DEtools[intfactor] and dsolve can accept extra arguments allowing the user to specify a particular form of the integrating factor. This is useful when the default approach fails in solving the problem or for other reasons.

>

ode := diff(y(x),x,x) = (y(x) + x*diff(y(x),x) + 2*x*y(x)*diff(y(x),x) + x^2*diff(y(x),x)^2)*diff(y(x),x)/x/y(x);

 (53) 
>

DEtools[intfactor]( ode, _mu = x*f(y,_y1) );

Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
> Calling odsolve with the ODE, diff(diff(y(x),x),x)*x^22*(1/2+x)*(diff(y(x),x)*xy(x)), y(x)
Methods for second order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying a Liouvillian solution using Kovacic's algorithm
A Liouvillian solution exists
Reducible group (found an exponential solution)
Group is reducible, not completely reducible
< Kovacic's algorithm successful
 
 (54) 

When an integrating factor is found, you can use DEtools[firint] to compute a related first integral (reduction of order; see also dsolve,education) or directly pass that integrating factor as a hint for dsolve.

>

ode := diff(y(x), x$3) = (1/y(x)*diff(y(x),x)2/x)*diff(y(x),x,x)+1/y(x)/x*diff(y(x),x)^2+1/x^3*y(x);

 (55) 
>

Mu := DEtools[intfactor]( ode );

Methods for first order ODEs:
 Trying classification methods 
trying a quadrature
trying 1st order linear
< 1st order linear successful
 
 (56) 
Integrating factor hint being investigated...
 
< differential order: 3; exact nonlinear successful
 
 (57) 
>

diff(y(x),x,x) = (11/4  x^(2*n)  2*x^n)*n^2/(x^n+1)^2*x^(2*n2)*y(x) + (n1)/x*diff(y(x),x);

 (58) 
Methods for second order ODEs:
 Trying classification methods 
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying an equivalence, under noninteger power transformations,
to LODEs admitting Liouvillian solutions.
> Trying a Liouvillian solution using Kovacic's algorithm
< No Liouvillian solutions exists
> Trying a solution in terms of special functions:
> Bessel
> elliptic
> Legendre
> Whittaker
> hyper3: Equivalence to 1F1 under 'a power @ Moebius'
< hyper3 successful: received ODE is equivalent to the 1F1 ODE
< Whittaker successful
< special function solution successful
 
 (59) 
>

diff(y(x),x,x) = (x+1)/(x1)^5*y(x);

 (60) 
Methods for second order ODEs:
 Trying classification methods 
trying a quadrature
checking if the LODE has constant coefficients
checking if the LODE is of Euler type
trying a symmetry of the form [xi=0, eta=F(x)]
checking if the LODE is missing 'y'
> Trying a Liouvillian solution using Kovacic's algorithm
< No Liouvillian solutions exists
> Trying a solution in terms of special functions:
> Bessel
> elliptic
> Legendre
> Whittaker
> hyper3: Equivalence to 1F1 under 'a power @ Moebius'
> hypergeometric
> heuristic approach
> hyper3: Equivalence to 2F1, 1F1 or 0F1 under 'a power @ Moebius'
< hyper3 successful: indirect Equivalence to 0F1 under ``^ @ Moebius`` is resolved
pFq equation: y(x)x*diff(y(x),x,x)diff(y(x),x)*_c = 0
Transformation resolving the equivalence: {x = 1/36*(t+1)^3/(t1)^3, y(x) = 1/(t+1)*u(t)}
Values of hypergeometric parameters: [], [_c = 4/3]
< hypergeometric successful
< special function solution successful
 
 (61) 
>

infolevel[dsolve] := 0;

 (62) 
•

The DEtools[redode] command can now compute the ODE, admitting possibly many different integrating factors at the same time (formerly, the command could compute it admitting only one integrating factor). This command is useful for studying the inverse relation between conservation laws and Hamiltonian systems, as well as understanding the integrating factor approach as an ODE solving technique.


Consider the following three integrating factors of some (unknown) third order ODE:

>

mu := [y(x)*x, x^(1/4+1/4*29^(1/2))*y(x), x^(1/41/4*29^(1/2))*y(x)];

 (63) 

The most general third order ODE, simultaneously admitting all these integrating factors, depends on an arbitrary function of and can be computed using redode as follows:

>

DEtools[redode]( mu, 3, y(x) );

 (64) 
>

map(DEtools[mutest], mu, (64));

 (65) 

The test performed above verifies  returning zero  that the given integrating factors (the list mu) indeed turn the returned ODE exact  see mutest.


The number of given integrating factors can be different than the differential order. For example, taking just the first two integrating factors of the previous example, the most general third order ODE turned exact by them actually depends  in addition  on an arbitrary function of two variables, :

 (66) 
>

DEtools[redode]( (66), 3, y(x) );

 (67) 
>

map(DEtools[mutest], mu[1..2], (67));

 (68) 
•

The DEtools[reduceOrder]  originally written for linear ODEs  was fixed and merged into DEtools[reduce_order] (introduced in Maple 7) for nonlinear ODEs, resulting in a single command with varied and extended capabilities. This renewed routine is now one of the main Maple tools for studying the processes by means of which ODEs of order higher than one  linear and nonlinear  are solved.

•

The capabilities of the conversion facility for ODEs (see convert,ODEs) have been extended. It is now possible to convert to normal form ODEs of type: linear (arbitrary order), Riccati and Abel. For Abel ODEs there is a new conversion to "Rational Normal Form". These tools are useful to study the solving properties of these relevant classes of equations. Also a conversion applicable to ODEs in general, called here 'y_x', to switch the roles between the dependent (y) and independent (x) variables, has been implemented. This is a useful conversion tool in different frameworks since this transformation enters many standard and alternative interactive solving approaches.

•

The PDEtools[dpolyform] command, introduced in Maple 7, now admits an optional extra argument to request a nonlinear  of lower differential order  representation for nonpolynomial objects, when that is possible, instead of a linear representation of higher differential order.



Numerical Solution of PDE



New to Maple 8, Maple has the ability to numerically solve timebased linear PDE over rectangular domains. The solution procedure uses the finite difference approach to obtain these numerical solutions. There are a number of aspects of this new capability.


Timebased systems



The default solver uses a centered implicit scheme, and is capable of finding numerical solutions of high order linear PDE and/or linear PDE systems. Note: The PDE systems must be sufficiently close to a standard form for the method to find the numerical solution.

>

PDE := diff(u(x,t),t,t)=1/10*diff(u(x,t),x,x,x,x);

 (69) 
>

IBC := {u(x,0)=sin(Pi*x)^2,D[2](u)(x,0)=0,
u(1,t)=0,D[1](u)(1,t)=0,
u(1,t)=0,D[1](u)(1,t)=0};

 (70) 
>

pds := pdsolve(PDE,IBC,numeric,spacestep=1/40);

 (71) 
>

pds:plot3d(t=0.5,color=[u,0,0],axes=boxed,
orientation=[45,45]);

>

pds:plot3d(t=0.5,color=[u,0,0],axes=boxed,
orientation=[45,45]);



Custom components



A call to pdsolve to obtain a numerical solution returns a module that has methods for producing 2 or 3D plots, animations, and evaluation procedures that can be used to compute solution values. In addition, the capability to change step sizes for the numerical computation is provided.




Numerical Solution of ODE



Additions and improvements have been made to Maple's capabilities for numerical ODE solution.


Additions


•

The default stiff and nonstiff methods can now produce a numerical solution over a range in terms of piecewise polynomial functions. These functions are constructed from the method interpolant at each step of the integration. For example,

>

dsys := {diff(y(x),x)=y(x),y(0)=1};

 (72) 
>

pw := dsolve(dsys,numeric,output=piecewise,range=0..0.2);

 (73) 
•

The plots[odeplot] command has been made more flexible for plotting multiple solution curves, and is now capable of producing animations of ODE system solutions.

>

dsys := {diff(y(x),x,x)=y(x),y(0)=0,D(y)(0)=1};

 (74) 
>

dsol := dsolve(dsys,numeric);

 (75) 
>

plots[odeplot](dsol,[[x,y(x),color=blue],
[x,diff(y(x),x),style=POINT],
[x,y(x)^2+diff(y(x),x)^2,linestyle=DASH]],
0..5);


produces a plot of the curves with the indicated options.


Addition of the option frames=50 after the range produces a 50 frame animation instead of a static plot.

•

It is now possible to obtain numerical solutions for ODE problems where the initial point has a removable singularity. For example, Bessel's equation can now be solved to compute BesselJ(1,x):

>

bsys := {x^2*diff(y(x),x,x)+x*diff(y(x),x)+(x^21)*y(x),
y(0)=0,D(y)(0)=1/2}:

>

Bes := subs(dsolve(bsys,numeric,output=listprocedure),y(x));

 (76) 
>

Bes(0.5),evalf(BesselJ(1,0.5));

 (77) 
>

Bes(1.0),evalf(BesselJ(1,1.0));

 (78) 

The removable singularity places restrictions on the initial data. If these restrictions are not met, it is reported:

>

dsolve({x^2*diff(y(x),x,x)+x*diff(y(x),x)+(x^21)*y(x),
y(0)=1,D(y)(0)=1},numeric);


This change is available for the rkf45, dverk78, gear, and classical methods, and is limited to singularities where all derivatives, up to the order of the system, are finite.

•

The higher accuracy dverk78 method now has a choice of local error computation strategies. In addition to the default strategy, which uses interpolation to estimate the error, use of the lower order method in the pair can be used to estimate the error with the option. Use of the pair requires approximately 1/2 the run time as the interpolation method, so difficult problems can be computed more quickly.

•

A new option, implicit=true, has been added to dsolve(..., numeric), and indicates that the input ODE system should not be isolated for the high order derivatives until the values of the inpedendent and dependent variables are known. This option is limited to linear systems, and the rkf45, dverk78, gear, lsode, and classical methods.


This option is useful for ODE systems with complicated coefficients, as it avoids the symbolic isolation of the leading derivatives of the system at the cost of a linear system solve at each function evaluation.


One example of where this option is needed is in use of a homotopy method to cast a problem in multivariate polynomial root finding to an ODE problem. In this case, the coefficients of the leading derivatives are large, and symbolic isolation of the leading derivatives results in significant resource use, and often huge illconditioned expressions.



Improvements


•

The core integrators for the default stiff and nonstiff methods have been redesigned to work in different operating modes. These can work storing the solution, not storing the solution, or computing solution values at a predetermined set of points.


The storage of the solution was the default in Maple 7, and is most useful for smaller problems, while earlier versions of Maple did not store the solution, which was more useful for longtime integration problems.


Now both capabilities are available through use of the range option.

•

The odeplot routine in the plots package has been enhanced to take advantage of the changes to the core integrators, and can now produce detailed plots much more efficiently.

•

Additionally, miscellaneous adjustments have been made: dsolve can now work properly with Heaviside and piecewise functions, and complex problem handling, and the stop condition feature have been improved.



