rifsimp[overview] - description of algorithm concepts and use

Example 1 : rifsimp form of a linear algebraic system

> with(Rif):
rifsimp([x - 2*y = 1, x + 2*z = 4]);

This is the same result that one would expect from the Gaussian elimination algorithm. Note that rifsimp has chosen to solve for x and y in terms of z . For this problem the order of solving, called a ranking , is x > y > z (the default ranking for rifsimp ). Given the inequality, rifsimp looks at the unknowns in an equation, and chooses the greatest one to solve for. This is a purely illustrative example, because Maple has specialized functions that solve such systems.

Example 2 : rifsimp form of a linear ODE system

> sys:= [3*x(t)-diff(x(t),t)=0, diff(x(t),t,t)-2*x(t)=0];

> rifsimp(sys);

As described in the prior example, rifsimp requires a ranking to determine the unknown for which to solve in an equation. For this ODE system, the ranking (using prime notation for derivatives) defaults to the following:

x  <  x' <  x'' <  x'''  <  ...

So the first equation is solved for x' giving x' = 3*x . The second equation is then differentially simplified with respect to the first, and gives 3*x'-2*x = 0 after the first simplification, and finally 7*x = 0 after the second. This equation is then solved for x , then back substitution into the first equation gives 0 = 0 .

For more general systems, rifsimp accounts for all of the derivative consequences of a system in a systematic and rigorous way.

Example 3 : rifsimp form of a nonlinear algebraic system

> sys := [x^10-2*x^5*y^2+y^4-y^3+1, 2*x^10-4*x^5*y^2+2*y^4+3*y^3-3,
x^5-y^2+y^3-1];

> rifsimp(sys);

The output was returned in the Constraint entries which means algebraic constraints as opposed to differential constraints. For this system, rifsimp gives the Groebner basis form of the algebraic system as output, as can be seen by running the Groebner package on the same system (with the appropriate ranking).

> Groebner[gbasis](sys,tdeg(x,y));

In general, rifsimp on algebraic systems is closer to a factoring case splitting Groebner basis (see rifsimp[nonlinear] , and Maple's Groebner package). rifsimp does not use Maple's Groebner package, but instead its own, optimized for differential elimination.

Example 4 : rifsimp form of a nonlinear ODE system

> sys:= [x(t)^3 + x(t) + 1 = 0,diff(x(t),t,t)-3*diff(x(t),t)=0,
diff(x(t),t,t,t)-2*diff(x(t),t)=0];

> rifsimp(sys);

Now the output is present in both the Solved and Constraint lists. The equation is present in the constraint list because the leading indeterminate of the equation occurs nonlinearly (for information on the leading indeterminate, please see rifsimp[ranking] ). Again, to obtain the answer, a ranking of the derivatives was required. The ranking used ranks "higher derivatives in x(t) before lower derivatives". This is a fairly standard ranking called the total degree ranking .

It is clear that the solution of the constraint equation gives x(t) = constant , so from that point of view, the equation in the derivative of x(t) is redundant, but it needs to be retained for the proper functioning of other algorithms, such as initialdata , and rtaylor . The additional equation in the Solved is a differential consequence of the Constraint , as can be verified through differentiation. This redundant equation is called a spawn equation, as it is a differential constraint spawned (through differentiation by t) from the nonlinear equation. Redundant equations can be removed with the clean or fullclean options as follows (see rifsimp[options] ):

> rifsimp(sys,clean=[spawn]);

> rifsimp(sys,fullclean);

Example 5 : rifsimp form of a linear PDE system

> sys:= [diff(u(x,y),x,x) - y*diff(u(x,y),x) = 0,
diff(u(x,y),y) = 0];

> rifsimp(sys);

In this example, rifsimp used the default ranking (where indices are used to denote differentiation)

u < u[x] < u[y] < u[xx] < u[xy] < u[yy] < ...

This ranking satisfies the following properties:

1. It is preserved under differentiation ( u < u[x] implies u[x] < u[xx] ).

2. u is less than any derivative of u .

These two properties are crucial for the termination of the algorithm. Fortunately, the default rankings chosen by rifsimp obey these properties, so these do not have to be considered. In general, many such orderings are possible, and methods for specifying them are discussed in rifsimp[ranking] and checkrank .

As an illustration, we will manually simplify the above system to rifsimp form. As the first step, we solve each equation for its leading derivative (given the ranking described above). This gives u[xx] = y u[x] and u[y] = 0 . Alas, we do not have our final answer, but no differential simplifications remain. How is this possible? Well, for PDE systems we have one further consideration called integrability conditions.

After looking at our solved system for a time, we may realize that differentiation of the first equation by y and differentiation of the second equation by x twice will give the derivative u[xyy] on the left hand side of both equations. If the solution to the system is differentiable to third order (an assumption implicitly made by rifsimp ), we differentiate the equations to the described order and subtract, giving the new, nontrivial equation

0 = D[y](u[xx])-D[xx](u[y]) = D[y](y u[x])+D[xx](0) = y u[xy]+u[x]

Differential simplification of this equation with respect to the current system yields u[x] = 0 , which is then used to simplify the system and give us the answer returned by rifsimp .

To summarize, for systems of linear PDEs, the rifsimp form of the system can be obtained by the following heavily simplified algorithm, for a given input ranking <:

Care must be taken in how the diff-gauss-elim routine operates on the solved rifsys.

Example 6 : rifsimp form of a nonlinear PDE system

This system is considered in the paper Reid et al. , 1996 (see Rif ) and is used as an example in Rust, 1998. For the ranking u < u[x] < u[y] < u[xx] < u[xy] < u[yy] < ... , the initial system consists of one leading linear PDE, u[xx]-u[xy] = 0 (with leading derivative u[xx] ), and one leading nonlinear, PDE u[y]^2+u[y]-u = 0 (with leading derivative u[y] ). No linear elimination can be done in the single leading linear PDE, so the leading nonlinear PDE is differentiated with respect to both independent variables to spawn the two leading linear PDEs: ( 2 u[y]+1) u[xy]-u[x] = 0 and (2 u[y]+1) u[yy]-u[y] = 0 (for more information on spawning, see rifsimp[nonlinear] ). These are solved for their leading derivatives, yielding u[xy] = u[x]/(2 u[y]+1) and u[yy] = u[y]/(2 u[y]+1) , subject to the pivot or inequation condition 2 u[y]+1 <> 0 . Integrability conditions are then taken across the leading linear PDE and simplified subject to the system, yielding one more PDE: u[x] u[y]-u[x]^2 = 0 . When this equation is spawned, the resulting equations vanish to zero upon reduction by the current system, so the algorithm then terminates. In summary, constraints are treated algebraically, and the leading linear equations are treated differentially. The system is then viewed modulo the constraint equations, using the concept of a relative basis , rigorously described in the thesis of Rust (see references on the Rif help page).

> sys:= [diff(u(x,y),x,x) - diff(u(x,y),x,y) = 0,
diff(u(x,y),y)^2 + diff(u(x,y),y) - u(x,y) = 0];

> rifsimp(sys);

Example 7 : rifsimp form of a nonlinear PDE system with cases

This problem arises from the determination of Lie symmetry groups for different values of a , b for the ODE:

> ODE := diff(y(x),x,x)-2*a^2*y(x)^3+2*a*b*x*y(x)-b;

The Lie symmetry determining equation can be obtained through use of the DEtools[odepde] function as follows:

> sys := {coeffs(numer(DEtools[odepde](ODE, [xi,eta](x,y), y(x))),_y1)};

When viewed as a system in xi and eta , it is a linear system, but if we are also considering a , b to be unknowns, the system becomes effectively nonlinear. Now what we want to accomplish is to find the equations for the symmetry generators for all values of a , b . We can do this by allowing rifsimp to perform case splitting on the system using the casesplit option (see rifsimp[cases] for more detail).

> ans:=rifsimp(sys,casesplit);

The structure of the output of rifsimp has changed form into a recursive table for multiple cases. The returned table has a casecount entry, and individual rifsimp answers for each case indexed by the numbers 1,2,3 . For the three cases, the form of the Solved equations for xi , eta is different. This indicates that for different values of a , b the structure of the Lie symmetry group is different, and provides the systems that must be satisfied for each group. This case splitting also covers all possible values of a , b , so it fully describes the dependence of the determining system on these values, and no other systems for other values of a , b exist.

The structure of the case tree can be observed in a plot using caseplot .

> caseplot(ans);

caseplot , checkrank , Groebner , Rif , rifsimp , rifsimp[cases] , rifsimp[nonlinear] , rifsimp[rankings]