(ProQuest: ... denotes non-US-ASCII text omitted.)
Recommended by B. V. Rathish Kumar
Departamento de Matemática Aplicada, Universidad de Alicante, Carretera San Vicente del Raspeig s/n, 03690 San Vicente del Raspeig, Alicante, Spain
Received 26 October 2011; Revised 25 March 2012; Accepted 26 March 2012
1. Introduction
The origin of symbolic manipulation derives from the sheer magnitude of the work involved in the building of perturbation theories, which made inevitable that scientific community became interested in the possibility of constructing those theories with the help of computers.
Perturbation theories for differential equations containing a small parameter ... are quite old. The small perturbation theory originated by Sir Isaac Newton has been highly developed by many others, and an extension of this theory to the asymptotic expansion, consisting of a power series expansion in the small parameter, was devised by Poincaré (1892) [1]. The main point is that for the most of the differential equations, it is not possible to obtain an exact solution. In cases where equations contain a small parameter, we can consider it as a perturbation parameter to obtain an asymptotic expansion of the solution. In practice, the work involved in the application of this approach to compute the solution to a differential equation cannot be performed by hand, and algebraic processors seem to be a very useful tool.
As explained in [2], the first symbolic processors were developed to work with Poisson series, that is, multivariate Fourier series whose coefficients are multivariate Laurent series, [figure omitted; refer to PDF] where Ci1 ,...,inj1 ,...,jm ∈... , i1 ,...,in , j1 ,...,jm ∈... , and x1 ,...,xn and [varphi]1 ,...,[varphi]m are called polynomial and angular variables, respectively. These processors were applied to problems in nonlinear mechanics or nonlinear differential equations problems, in the field of Celestial Mechanics.
One of the first applications of these processors was concerned with the theory of the Moon. Delaunay invented his perturbation method to treat it and spent 20 years doing algebraic manipulations by hand to apply it to the problem. Deprit et al. [3, 4] prolongated the solution of Delaunay's work with the help of a special purpose symbolic processor, and Henrard [5] pushed it to order 25. This solution was improved by iteration by Chapront-Touzé [6], and planetary perturbations were also introduced by Chapront-Touzé [7]. At present, the most complete solution, ELP (Ephemeride Lunaire Parisien) contains more than 50000 periodic terms.
But the motion of the Moon is not the only application of algebraic processors. There are many problems where the facilities provided by Poisson series processors can lead rather quickly to very accurate results. As examples, we would like to mention planetary theories, the theory of the rotation of the Earth (e.g., [8]), and artificial satellite theories (AST). Abad et al. [9] have analysed the convenience of developing specific computer algebra systems in order to deal with AST. As it is explained in detail in their work, the series involved in the computation of the solution through the application of a Lie transformation have a total amount of almost 2000000 terms.
Nowadays, many general purpose computer algebra packages--as, for example, Maple, Mathematica, and Matlab--contain tools for the calculation of the solution of certain classes of ODEs. All these packages have the advantage of being very general, so they can deal with a lot of problems of different nature. However, if one is interested in higher order solutions, the most common perturbation methods tend to produce expressions containing thousands of terms, and their treatment with those general processors becomes inefficient.
To achieve better accuracies in the applications of analytical theories, high orders of the approximate solution must be computed, making a continuous maintenance and revision of the existing symbolic manipulation systems necessary, as well as the development of new packages adapted to the peculiarities of the problem to be treated.
In order to contribute to the solution of this problem, we have developed a symbolic computation package (called MQSP) based on the object-oriented philosophy which manipulates objects of the form [figure omitted; refer to PDF] where nν ∈..., αν ,ων ,λν ,μν ∈... , σν ∈... , and τ1 ,...,τQ are real constants with unknown value. We will refer to those elements as modified quasipolynomials [10]. The kernel of this symbolic processor has been developed in C++. The operations on series implemented in the manipulator are the usual operations of the Algebra of quasipolynomials: addition, subtraction, multiplication, multiplication by a scalar, differentiation with respect to t , and substitution of a quasipolynomial into an undetermined coefficient.
We have also constructed a set of subroutines to deal with the solution to the perturbed differential equation [figure omitted; refer to PDF]
where a0 ,a1 ,t0 ,x0 ,x...0 ∈... , u(t) is given by (1.2), and [figure omitted; refer to PDF]
In a previous contribution, the author employed the kernel of this processor to compute periodic solutions in equations of type (1.3) via the Poincaré-Lindstedt method [11, 12]. If the unperturbed equation ( ...=0 ) has periodic solutions and ... is a measure of the size of the perturbing terms, then the trajectories for the full system will remain pretty close to those of the nonperturbed system, for any finite period of time t0 <t<t0 +α ( α>0 ) with an error not larger than O(α) . In general, even a small perturbation is enough to destroy periodicity, that is, nonlinearity will end with most of the periodic orbits of the unperturbed system, but some of them may persist. To calculate those periodic orbits, the solution and the modified frequency are expanded with respect to the small parameter, allowing to kill secular terms which appear in the recursive scheme.
The aim of this paper is to construct an algorithm to implement the asymptotic expansion method [13]. This new implementation is general and does not depend on the function f(x,x...) as given in (1.4), that is, the user does not need to programme the algorithm described in this paper, and only has to introduce the adequate parameters when calling the corresponding routine of the package. The code of this specific symbolic system is not available on internet but it can be provided by contacting its author.
2. Data Structure
The algorithms that can be implemented to perform the basic manipulation on a series and their efficiency depend on the way a series is coded. An overcoded structure that makes good use of memory generally requires complex algorithms, which increase the computational cost in terms of time. On the other hand, an undercoded computational representation of the terms generates simple algorithms, because the location of all the coefficients can be obtained directly. However, this scheme presents the inconvenience of being very wasteful in the memory resources required for the storage of the series [2].
In this section, we follow San-Juan and Abad [14] to introduce the representation of a mathematical object in a computer. To that purpose, let us introduce the concepts of normal and canonical functions. Let E be a set of symbolic objects, and let ~ be an equivalence relation in E , defined as follows: a~b if a=b , with a,b∈E . Here, the operator = is considered as the equality on the mathematical object level. Moreover, a...1;b if a and b are identical as symbolic objects. A function f:E[arrow right]E is said to be normal in (E,~) if f(a)~a for all a∈E , and f is said to be canonical in (E,~) if it is normal and a~b[implies]f(a)=f(b) for all a,b∈E . Thus, a canonical function provides identical objects when objects are equivalent, that is, when they represent the same mathematical object.
For the sake of simplicity, let us focus our attention on the set of quasipolynomials in the independent variable Q[t] . A quasipolynomial is a map u:...[arrow right]... , defined by [figure omitted; refer to PDF] where nν ∈... , αν , ων , λν , and μν ∈... . Let us now consider the set of quasipolynomials in the independent variable t , Q[t] .
2.1. QuasiPolynomials as Symbolic Objects
We look for a canonical representation for each equivalence class defined in Q[t] . For that purpose, the following operations must be performed over each quasipolynomial.
(1) : Let us consider a quasipolynomial u(t)=∑ν...5;0tnνeαν t (λν cos (ων t)+μν sin (ων t)) . If ων <0 , the following rules must be applied:
[figure omitted; refer to PDF]
(2) : The terms of a quasipolynomial will be ordered as follows: let us consider two term of a quasipolynomial: τ1 =tn1eα1 t (λ1 cos ω1 t+μ2 sin ω1 t) and τ2 =tn2eα2 t (λ2 cos ω2 t+μ2 cos ω2 t) . We say that τ1 <τ2 if ( n1 <n2 ) or ( n1 =n2 and α1 <α2 ) or ( n1 =n2 , α1 =α2 , and ω1 <ω2 ) or ( n1 =n2 , α1 =α2 , ω1 =ω2 , and λ1 <λ2 ) or ( n1 =n2 , α1 =α2 , ω1 =ω2 , λ1 <λ2 , and μ1 <μ2 ).
(3) : The terms of a quasipolynomial
[figure omitted; refer to PDF]
with identical values of nν , αν , and ων must be grouped together.
2.2. QuasiPolynomials as Computational Objects
In this section, we will consider the basic information which characterizes a quasipolynomial, as well as the data structure to store it in the computer. This must be done preserving the canonical representation we have chosen.
Each quasipolynomial is collected in a sorted dynamical list: a sorted list is one in which the order of the items is defined by some collating sequence. The codification of each term of the list contained in a quasipolynomial is statical, and given by the following elements.
(1) : λ,μ∈... are the coefficients of the term.
(2) : n∈... is the degree of the monomial t .
(3) : α∈... is the exponent of the exponential part of the term.
(4) : ω...5;0 is the frequency of the periodic part.
In Scheme 1 we represent the way a quasipolynomial is stored in memory by using a simple list.
Scheme 1: [figure omitted; refer to PDF]
As pointed out in [14], most of the operations involving a series are based on navigating and searching through the structure that represents the series. For example, the addition of two quasipolynomials is equivalent to inserting each term of one series into the other one. So, a good choice of the data structure results in simple and efficient algorithms. The binary tree resulting seems to be a very useful data structure for rapidly storing sorted data and rapidly retrieving saved data. A binary tree is composed of parent nodes, or leaves, each of which stores data and also links to up to two other child nodes (leaves), which can be visualized spatially as below the first node with one placed to the left and the other to the right. In this structure, the relationship between the linked leaves and the linking leaf makes the binary tree an efficient data structure: the leaf on the left has a lesser key value, and the leaf on the right has an equal or greater key value.
A special type of tree is the red-black tree. In a red-black tree, each node has a color attribute, the value of which is either red or black. In addition to the ordinary requirements imposed on binary search trees, the following additional requirements of any valid red-black tree apply: A node is either red or black. The root is black. All leaves are black, even when the parent is black. Both children of every red node are black. Every simple path from a node to a descendant leaf contains the same number of black nodes. A critical property of red-black trees is enforced by these constraints: the longest path from the root to a leaf is no more than twice as long as the shortest path from the root to a leaf in that tree. The result is that the tree is roughly balanced. Since operations such as inserting, deleting, and finding values requires worst-case time proportional to the height of the tree, this fact makes the red-black tree efficient, for instance, the search-time results to be O(log n) .
With the use of this structure, the complexity of the algorithms for addition, multiplication, derivation, and integration of quasipolynomials is significantly reduced. Unfortunately, this structure, which results to be ideal for Poisson series, cannot be applied directly in our case due to the fact that the numbers which identify each term of a quasipolynomial are not indexed arrays. However, an alternative aggrupation of terms in a quasipolynomial can be performed in order to introduce this balanced structure. To do that, let us express a quasipolynomial as follows: [figure omitted; refer to PDF]
where Cp,ν (t) and Sq,ν (t) are polynomials in the variable t with constant coefficients, of degree p and q respectively, being p,q∈... : [figure omitted; refer to PDF] with λ0 ,λ1 ,...,λp , μ0 ,μ1 ,...,μq ∈... . If we aggrupate terms of a quasipolynomial in such a manner, we can use a tree structure to store it, saving only significant terms. In Scheme 2 we show the tree structure in which the quasipolynomial is stored.
Scheme 2: [figure omitted; refer to PDF]
In Scheme 2 we show how we store in memory the quasipolynomial [figure omitted; refer to PDF] where [figure omitted; refer to PDF]
with ...C,j and ...S,j being sets of indices for any j=1,2 , ν1,i ...5;0 and λ1,i ...0;0 for any i∈...C,1 , ν2,i ...5;0 and λ2,i ...0;0 for any i∈...C,2 , σ1,i ...5;0 and μ1,i ...0;0 for any i∈...S,1 , and σ2,i ...5;0 and μ2,i ...0;0 for any i∈...S,2 . Moreover, νj,i-1 <νj,i <νj,i+1 and σj,i-1 <σj,i <σj,i+1 for any j=1 , 2 and i belonging to the corresponding set of indices.
To illustrate the way a quasipolynomial is stored by means of the tree structure described above, we show in Scheme 3 the structure associated to the quasipolynomial given by [figure omitted; refer to PDF]
Scheme 3: [figure omitted; refer to PDF]
2.3. Modified QuasiPolynomials as Computational Objects
In our symbolic system, we represent a modified quasipolynomial by an ordered and dynamical list of terms, keeping in memory only significant terms. The codification of a term is statical, and given by the following elements.
(1) λ,μ∈... are the coefficients of the term (for cos and sin , resp.).
(2) (σ1 ,...,σQ )∈...Q . For each 0...4;ν...4;Q , σν represents the exponent of the undetermined coefficient τν .
(3) n∈... is the degree of the monomial t .
(4) α∈... is the exponent of the exponential part of the term.
(5) ω...5;0 is the frequency of the periodic part.
We have included the two following additional vectors, which are common to all the terms of a modified quasipolynomial:
(1) : (m1 ,...,mQ )∈{0,1}Q . For each 0...4;ν...4;Q , mν =0 indicates that τν is a coefficient with unknown value, while mν =1 implies that τν has a real value assigned,
(2) : τ=(τ1 ,...,τQ )∈...Q . If the coefficient τν has a real value assigned, its content is given by τν , for each 0...4;ν...4;Q .
It is also absolutely essential to store the number of undetermined coefficients that are currently in use, to generate correctly new undetermined constants if needed. In Table 1 we illustrate the way a term is coded with a few examples, where it has been assumed that Q=4 .
Table 1: Some examples of representation of terms.
Term | λ | μ | σ1 | σ2 | σ3 | σ4 | n | α | ω |
3τ12τ4 cos (5t) | 3 | 0 | 2 | 0 | 0 | 1 | 0 | 0 | 5 |
-τ1τ2τ32t7e-t sin (3t) | 0 | -1 | 1 | 1 | 2 | 0 | 7 | -1 | 3 |
2τ2τ32 te2t cos (3t) | 2 | 0 | 0 | 1 | 2 | 0 | 1 | 2 | 3 |
3. Design of the Symbolic System in C++
As pointed out in [15], object-oriented programming using C++ provides many advantages in the design of computer algebra systems, as this programming technique combines both the data and the functions that operate on that data into a single unit (called class). The main reasons given by Hardy et al. to use C++ to implement a symbolic system are as follows.
(1) : C++ allows the introduction of abstract data types. Thus, we can define a modified quasipolynomial as an abstract data type.
(2) : The language C++ supports encapsulation, inheritance, polymorphism, and operator overloading. Consequently, we can overload the operators +, -, and * for modified quasipolynomials, as well as * and / for multiplication and division of modified quasipolynomials by real numbers.
Some symbolic computation systems have been constructed in C++. MuPAD is a computer algebra system developed by the MuPAD research group at the University of Paderborn, in Germany. This symbolic system manipulates formulas symbolically and provides packages for linear algebra, differential equations, number theory, statistics, and functional programming, as well as an interactive graphic system that supports animations and transparent areas in 3D. MuPAD also offers a programming language that supports object-oriented programming and functional programming [16].
Symbolic C++ also uses C++ to develop a computer algebra system. This package introduces the Symbolic class which is used for all symbolic computation, and provides almost all of the features required for symbolic computation including symbolic terms, substitution, noncommutative multiplication, and vectors and matrices.
Both symbolic systems could have been used to implement the algorithm over a general symbolic class. However, as the goal of the symbolic system we are developing is to handle modified quasipolynomial to apply perturbation methods to solve some types of differential equations, we have constructed it directly over C++, instead of using some other symbolic processor.
The specific symbolic processor designed is written in clean C++ code, is very portable, it can compile stand-alone, and is easily embeddable. It implements a new data type, called MQ, which represents a series of the form given by (1.2). The class MQ is defined as an ordered and dynamical list of terms. To that end, we have also defined a class associated to the structure of a term of a modified quasipolynomial, called Term. The definition of these two classes in C++ code is shown in Algorithms 1 and 2, respectively.
Algorithm 1: C++ code for the definition of a modified quasipolynomial.
typedef Term *serie;
class MQ
{
public:
[definition of functions]
private:
serie first;
};
Algorithm 2: C++ code for the definition of a term of a modified quasipolynomial.
class Term
{
public:
Term ();
Term (double L, double M,
int *sigm, int m, double alph, double omeg,
Term *nxt);
~Term ();
private:
double lambda, mu;
int *sigma, n;
double alpha, omega;
Term *next;
friend class MQ;
};
The set of routines developed includes: addition, subtraction and product of series, multiplication of a series by a-real number, differentiation with respect to t , and computation of the solution to a linear second-order differential equation of type [figure omitted; refer to PDF] where u(t) is a modified quasipolynomial presenting undetermined coefficients, and computation of the solution to (1.3) via the asymptotic expansion method. These two algorithms are described in detail in Sections 4 and 5, respectively.
The basic algebra associated to quasipolynomials is easily implemented because of the undercoded data scheme chosen for their computational representation. Thus, for example, the addition of two quasipolynomials is performed by directly juxtaposing both quasipolynomials, arranging the resulting series, and joining terms with equal elements. In a similar way, the rest of algebraic operations are simply implemented.
4. Solution of a Linear Second-Order ODE
The general solution to a nonhomogeneous differential equation can be expressed as the sum of general solutions to the corresponding homogenous, linear differential equation and any solution to the complete equation [17]. The solution to the homogeneous ODE is expressed in terms of the roots of the characteristic equation, α2 +a1 α+a0 =0 , and it is well-known. We will resume now the formulae that are required to construct a particular solution to a complete ordinary differential equation of second-order (3.1). Without loss of generality, we will assume that u(t) is written as follows: [figure omitted; refer to PDF] where α,ω∈... , and pν,n (t,τ) , and qν,m (t,τ) are n th and m th degree polynomials in t with undetermined coefficients respectively, of the form [figure omitted; refer to PDF] being [figure omitted; refer to PDF]
with uν,ρ* ,vν,ρ* ∈... , σν,ρ,i ,σν,ρ,i* ∈... , 1...4;i...4;Q , 0...4;ρ...4;n for uν,ρ , and 0...4;ρ...4;m for vν,ρ .
The principle of superposition is applied to calculate the particular solution, so we can focus our attention on the computation of a particular solution to the equation [figure omitted; refer to PDF] where [figure omitted; refer to PDF] being [figure omitted; refer to PDF] with uρ* ,vρ* ∈... , σρ,i ,σρ,i* ∈... , 1...4;i...4;Q , 0...4;ρ...4;n for uρ , and 0...4;ρ...4;m for vρ .
At this point, we will distinguish two cases depending on if ω=0 or ω...0;0 .
Case 1.
Let us consider the second-order ODE (4.4), with ω=0 . Then, the equation is written as [figure omitted; refer to PDF] where pn (t,τ) is given by (4.5).
Subcase 1.1.
If α2 +a1 α+a0 ...0;0 , the particular solution to the complete differential equation is expressed as
[figure omitted; refer to PDF] Thus, substituting x(t) and its derivatives in (4.7), we get that [figure omitted; refer to PDF] and, for any p<n-1 , [figure omitted; refer to PDF]
Subcase 1.2.
If α2 +a1 α+a0 =0 and a12 ...0;4a0 (i.e., α...0;-a1 /2 ), the particular solution can be written as
[figure omitted; refer to PDF] Now, the substitution of x(t) and its derivatives into (4.7) leads to the following formula for αn (τ) , [figure omitted; refer to PDF] and, for any p...4;n-1 , [figure omitted; refer to PDF]
Subcase 1.3.
If α2 +a1 α+a0 =0 and a12 =4a0 (i.e., α=-a1 /2 ), the particular solution can be written as
[figure omitted; refer to PDF] Now, the substitution of x(t) and its derivatives into (4.7) leads to the following formula for αp (τ) , [figure omitted; refer to PDF]
Case 2.
Let us consider the second-order differential equation (4.4), where now ω...0;0 . There are two possible situations:
(1) : α±iω is not a root of the characteristic equation, that is, α ...0;-a1 /2 or ω...0;4a0 -a12 /2 ,
(2) : α ± iω is a root of the characteristic equation, that is, α=-a1 /2 and ω=4a0 -a12 /2 .
Subcase 2.1.
Let us call N=max {n,m} . Then,
[figure omitted; refer to PDF] Now, if we assume pN (t,τ) and qN (t,τ) to be [figure omitted; refer to PDF] and substitute x(t) and its derivatives into (4.4), we get that [figure omitted; refer to PDF] where Δ=(a0 +a1 α+α2 -ω2)2 +ω2 (a1 +2α)2 . Note that Δ...0;0 , because a0 +a1 α+α2 -ω2 ...0;0 or ω(2α+a1 )...0;0 . Now, we can compute αN-1 (τ) and ρN-1 (τ) by solving the system [figure omitted; refer to PDF] As before, this system can be solved by applying the Cramer's rule. Finally, for any p<N-1 , we have to solve the system in αp (τ) and ρp (τ) , as follows: [figure omitted; refer to PDF]
Subcase 2.2.
Now we consider the case where α=-a1 /2 and ω=4a0 -a12 /2 . This implies that a1 +2α=0 and α2 -ω2 +a1 α+a0 =0 . In this case, the particular solution adopts the form
[figure omitted; refer to PDF] and substituting x(t) , x...(t) , and x...(t) into (4.4), we obtain [figure omitted; refer to PDF] for any p<N .
From (4.8)-(4.22), it is a straightforward task; the derivation of an algorithm for the computation of the solution to (4.4).
5. Computation of the Solution to the Perturbed Problem
The symbolic package is thought to compute the solution to (1.3). The standard approach [13] is to try a power series solution of the form [figure omitted; refer to PDF] This series is inserted into the governing equation and initial conditions, and coefficients of same powers of ... are then grouped to obtain a collection of equations for the coefficient functions xi (t) , which are then solved in a sequential manner. The resulting series need not converge for any value of ... ; nevertheless, the solution x(t,...) can be useful in approximating the function when ... is small.
Considering the zero-order term in ... yields [figure omitted; refer to PDF] the so-called nonperturbed problem. The symbolic manipulation system calculates the solution to a differential equation of the form (5.2), and arranges it as a quasipolynomial, as it has been described in detail in the previous section.
The coefficient xq (t) of the solution to the order q...5;1 is computed by solving the equation [figure omitted; refer to PDF] where the notation (xνx...κ-ν)q refers to the q th order term of the series xνx...κ-ν .
At each order of the solution, the series (xν)q , (x...ν)q , and (xνx...κ-ν)q must be computed once the order q has been solved, for each q...5;0 , following the formulae [figure omitted; refer to PDF]
According to this, the algorithm to apply the asymptotic expansion method to solve the initial value problem given by (1.3), consists of the following steps.
(1) : Define a three-dimensional array of quasipolynomials X(ρ1 , ρ2 , q) , where ρ1 ,ρ2 ,q∈... , 0...4;ρ1 , ρ2 ...4;M , and 0...4;q...4;Q , Q being the order of the asymptotic expansion.
(2) : Define an array of quasipolynomials x(ρ) , where 0...4;ρ...4;Q .
(3) : Initialize X(0,0,0)=1 , and proceed as follows:
(3.1) : Compute X(1,0,0) as the solution to (5.2).
(3.2) : Compute X(0,1,0)=d/dt(X(1,0,0)) .
(3.3) : Calculate, for each ρ such that 2...4;ρ...4;M , [figure omitted; refer to PDF]
(3.4) : For each ρ1 , ρ2 such that 1...4;ρ1 , ρ2 ...4;M , determine the quasipolynomial [figure omitted; refer to PDF]
(4) : For each q such that 1...4;q...4;Q , do the following.
(4.1) : Compute the quasipolynomial [figure omitted; refer to PDF]
(4.2) : Calculate X(1,0,q) as the solution to (5.3), [figure omitted; refer to PDF]
(4.3) : Compute X(0,1,q)=d/dt(X(1,0,q)) .
(4.4) : Calculate, for each ρ such that 2...4;ρ...4;M , [figure omitted; refer to PDF]
(4.5) : For each ρ1 ,ρ2 such that 1...4;ρ1 , ρ2 ...4;M , compute [figure omitted; refer to PDF]
(5) : For each ρ such that 0...4;ρ...4;Q , [figure omitted; refer to PDF]
The input arguments of the algorithm consist of the order Q of the asymptotic approximation, the coefficients a1 , a0 of the differential equation, the real values t0 , x0 , and x...0 which define the initial conditions of the problem, the quasipolynomial u(t) , and the perturbed part of the equation f(x,x...) . As f can be written as given in (1.4), [figure omitted; refer to PDF]
it can be specified by a real (M+1)×(M+1) matrix, [figure omitted; refer to PDF]
The parameter M is also required as input.
The output argument of the algorithm is the array x(ρ) containing the coefficients of the asymptotic expansion. The asymptotic approximation to the order Q is given by [figure omitted; refer to PDF]
6. On the Implementation of the Algorithm
In Algorithms 3, 4, 5, 6, and 7 we show the C++ code for the implementation of the algorithm as described in Section 5, omitting error control sentences for the sake of simplicity. The head of the definition of the function that implements the algorithm is shown below.
Algorithm 3: Implementation of the algorithm. Steps (1) and (2).
MQ xdx [M+1] [M+1] [ORDER+1] ; Define a three-dimensional array of quasi
polynomials X(ρ1 ,ρ2 ,q) , where ρ1 ,ρ2 ,q∈
N , 0...4;ρ1 ,ρ2 ...4;M and 0...4;q...4;Q , being
Q the order of the asymptotic expansion.
MQ *s; Define an array of quasipolynomials x(ρ) ,
s=new MQ [ORDER+1] ; where 0...4;ρ...4;Q .
int i, j, k, p, q, nu ; Definition of auxiliary variables.
MQ U;
Algorithm 4: Implementation of the algorithm. Steps (3.1) and (3.2).
xdx [0] [0] [0] .add Initialize X(0,0,0)=1 .
(1., 0., 0, 0., 0.) ;
xdx [1] [0] [0] = u.solvePVI Compute X(1,0,0) as the solution to (5.2),
(a1, a0, t0, x0, dx0); and X(0,1,0)=d/dt(X(1,0,0)) .
xdx [0] [1] [0] = xdx [1] [0] [0] .der ();
Algorithm 5: Implementation of the algorithm. Steps (3.3) and (3.4).
for (i=2; i<=m; i++) { Calculate, for each ρ such that 2...4;ρ...4;M ,
xdx [i] [0] [0] =
xdx [1] [0] [0] *xdx [i-1] [0] [0] ; X(ρ,0,0)=X(1,0,0)×X(ρ-1,0,0) ,
xdx [0] [i] [0] =
xdx [0] [1] [0] *xdx [0] [i-1] [0] ; X(0,ρ,0)=X(0,1,0)×X(0,ρ-1,0) .
}
for (i=1; i<=m; i++) For each ρ1 ,ρ2 such that 1...4;ρ1 ,ρ2 ...4;M ,
for (j=1; j<=m; j++) determine the quasipolynomial
xdx [i] [j] [0] =
xdx [i] [0] [0] *xdx [0] [j] [0] ; X(ρ1 ,ρ2 ,0)=X(ρ1 ,0,0)×X(0,ρ2 ,0) .
Algorithm 6: Implementation of the algorithm. Step (4).
for (q=1; q<=order; q++) For each q such that 1...4;q...4;Q , do:
{
U=U*0.; Compute the quasipolynomial
for (k=0; k<=m; k++)
for (nu=0; nu<=k; nu++) U=∑κ=0M ∑0...4;ρ...4;κfρ,κ-ρ ×X(ρ,κ-ρ,q-1).
U = U + xdx [nu] [k-nu] [q-1]
*f [nu] [k-nu];
xdx [1] [0] [q] = U.solvePVI Calculate X(1,0,q) as the solution to (5.3),
(a1, a0, t0, 0., 0.);
xdx [0] [1] [q] = xdx [1] [0] [q] .der (); x...q +a1x...q +a0xq =U ,
xq (t0 )=0, x...q (t0 )=0 ,
and X(0,1,q)=d/dt(X(1,0,q)) .
for (i=2; i<=m; i++) Calculate, for each ρ such that 2...4;ρ...4;M ,
for (nu=0; nu<=q; nu++) {
xdx [i] [0] [q] = xdx [i] [0] [q] X(ρ,0,q)=∑0...4;p...4;q X(ρ-1,0,p)
+ xdx [i-1] [0] [nu]
*xdx [1] [0] [q-nu] ; ×X(1,0,q-p) ,
xdx [0] [i] [q] = xdx [0] [i] [q] X(0,ρ,q)=∑0...4;p...4;q X(0,ρ-1,p)
+ xdx [0] [i-1] [nu]
*xdx [0] [1] [q-nu] ; ×X(0,1,q-p) .
}
for (i=1; i<=m; i++) For each ρ1 ,ρ2 such that 1...4;ρ1 ,ρ2 ...4;M ,
for (j=1; j<=m; j++) compute
for (nu=0; nu<=q; nu++)
xdx [i] [j] [q] = X(ρ1 ,ρ2 ,q)
xdx [i] [j] [q]
+ xdx [i] [0] [nu] =∑0...4;p...4;q X(ρ1 ,0,p)×X(0,ρ2 ,q-p) .
*xdx [0] [j] [q-nu] ;
}
Algorithm 7: Implementation of the algorithm. Step (5).
for (i=0; i<=order; i++){ For each ρ such that 0...4;ρ...4;Q ,
s [i] = xdx [1] [0] [i] ;
s [i].normalice (); x(ρ)=X(1,0,ρ), and x(t)=∑0...4;ρ...4;Q...ρ x(ρ).
s [i].order ();
s [i].join ();
s [i].neglect (PREC);
}
MQ *solveAM (double a1, double a0, double t0, double x0, double dx0, MQ u, double **f, int m, int order)
7. Example
Let us consider the initial value problem given by [figure omitted; refer to PDF] where ...=10-1 . Here, a1 =0 , a0 =1 , t0 =0 , x0 =1 , x...0 =0 , M=3 , and f is defined by the following 4 × 4 matrix: [figure omitted; refer to PDF]
In Table 2 we show the output generated from the symbolic manipulator, just to order nine for the sake of simplicity, as the number of terms of the asymptotic expansion increases very quickly with the order of the expansion.
Table 2: Some coefficients of the asymptotic expansion of the solution to (7.1).
x0 = | +1.000000cos (t) |
|
|
| |||
x1 = | -0.031250cos (t) | +0.031250cos (3t) | +0.500000sin (t) |
| -0.500000tcos (t) | -0.375000tsin (t) |
|
| |||
x2 = | +0.022461cos (t) | -0.023438cos (3t) | +0.000977cos (5t) |
| +0.226563sin (t) | +0.070313sin (3t) | -0.390625tcos (t) |
| -0.046875tcos (3t) | -0.031250tsin (t) | -0.035156tsin (3t) |
| +0.054688t2 cos (t) | +0.375000t2 sin (t) |
|
| |||
x3 = | +0.047760cos (t) | -0.046326cos (3t) | -0.001465cos (5t) |
| +0.000031cos (7t) | +0.122640sin (t) | -0.005859sin (3t) |
| +0.003988sin (5t) | -0.201660tcos (t) | +0.079102tcos (3t) |
| -0.002441tcos (5t) | -0.133057tsin (t) | -0.059692tsin (3t) |
| -0.001831tsin (5t) | +0.103760t2 cos (t) | +0.015381t2 cos (3t) |
| +0.140625t2 sin (t) | +0.070313t2 sin (3t) | +0.084635t3 cos (t) |
| -0.194336t3 sin (t) |
|
|
| |||
x4 = | +0.031557cos (t) | -0.025802cos (3t) | -0.005688cos (5t) |
| -0.000069cos (7t) | +0.000001cos (9t) | +0.314784sin (t) |
| +0.001236sin (3t) | -0.003563sin (5t) | +0.000181sin (7t) |
| -0.329605tcos (t) | +0.018219tcos (3t) | +0.009552tcos (5t) |
| -0.000107tcos (7t) | -0.441519tsin (t) | +0.023666tsin (3t) |
| -0.006212tsin (5t) | -0.000080tsin (7t) | +0.340805t2 cos (t) |
| -0.113068t2 cos (3t) | +0.001335t2 cos (5t) | +0.323486t2 sin (t) |
| -0.007141t2 sin (3t) | +0.005493t2 sin (5t) | -0.113566t3 cos (t) |
| +0.031860t3 cos (3t) | -0.168498t3 sin (t) | -0.064362t3 sin (3t) |
| -0.081533t4 cos (t) | +0.060547t4 sin (t) |
|
| |||
x5 = | +0.101090cos (t) | -0.103018cos (3t) | +0.002327cos (5t) |
| -0.000396cos (7t) | -0.000003cos (9t) | +0.000000cos (11t) |
| +0.567049sin (t) | +0.011443sin (3t) | -0.003304sin (5t) |
| -0.000300sin (7t) | +0.000007sin (9t) | -0.630124tcos (t) |
| +0.042347tcos (3t) | +0.004305tcos (5t) | +0.000649tcos (7t) |
| -0.000004tcos (9t) | -0.861960tsin (t) | -0.136819tsin (3t) |
| +0.014968tsin (5t) | -0.000407tsin (7t) | -0.000003tsin (9t) |
| +0.769419t2 cos (t) | +0.056821t2 cos (3t) | -0.019625t2 cos (5t) |
| +0.000082t2 cos (7t) | +0.653250t2 sin (t) | -0.009047t2 sin (3t) |
| -0.002311t2 sin (5t) | +0.000320t2 sin (7t) | -0.271820t3 cos (t) |
| +0.077614t3 cos (3t) | +0.003465t3 cos (5t) | -0.512211t3 sin (t) |
| +0.066607t3 sin (3t) | -0.007243t3 sin (5t) | +0.084625t4 cos (t) |
| -0.054769t4 cos (3t) | +0.110530t4 sin (t) | +0.027557t4 sin (3t) |
| +0.046019t5 cos (t) | -0.005018t5 sin (t) |
|
| |||
x6 = | +0.041699cos (t) | -0.035635cos (3t) | -0.006540cos (5t) |
| +0.000498cos (7t) | -0.000022cos (9t) | -0.000000cos (11t) |
| +0.000000cos (13t) | +1.753506sin (t) | -0.124680sin (3t) |
| +0.000293sin (5t) | -0.000409sin (7t) | -0.000018sin (9t) |
| +0.000000sin (11t) | -1.732879tcos (t) | +0.361830tcos (3t) |
| -0.007317tcos (5t) | +0.000416tcos (7t) | +0.000035tcos (9t) |
| -0.000000tcos (11t) | -1.950362tsin (t) | +0.090323tsin (3t) |
| -0.009679tsin (5t) | +0.001762tsin (7t) | -0.000022tsin (9t) |
| -0.000000tsin (11t) | +1.753296t2 cos (t) | -0.260570t2 cos (3t) |
| +0.014763t2 cos (5t) | -0.001785t2 cos (7t) | +0.000004t2 cos (9t) |
| +1.911049t2 sin (t) | +0.270222t2 sin (3t) | -0.022385t2 sin (5t) |
| -0.000220t2 sin (7t) | +0.000016t2 sin (9t) | -1.003710t3 cos (t) |
| -0.074659t3 cos (3t) | +0.018743t3 cos (5t) | +0.000255t3 cos (7t) |
| -1.212436t3 sin (t) | -0.124235t3 sin (3t) | +0.019212t3 sin (5t) |
| -0.000552t3 sin (7t) | +0.305366t4 cos (t) | -0.009847t4 cos (3t) |
| -0.008565t4 cos (5t) | +0.511480t4 sin (t) | -0.077788t4 sin (3t) |
| +0.004177t4 sin (5t) | -0.061540t5 cos (t) | +0.046341t5 cos (3t) |
| -0.067702t5 sin (t) | +0.007443t5 sin (3t) | -0.018889t6 cos (t) |
| -0.007806t6 sin (t) |
|
|
| |||
x7 = | -0.030141cos (t) | +0.026721cos (3t) | +0.003535cos (5t) |
| -0.000160cos (7t) | +0.000046cos (9t) | -0.000001cos (11t) |
| -0.000000cos (13t) | +0.000000cos (15t) | +5.481497sin (t) |
| +0.136811sin (3t) | -0.022934sin (5t) | +0.000391sin (7t) |
| -0.000032sin (9t) | -0.000001sin (11t) | +0.000000sin (13t) |
| -5.726072tcos (t) | -0.099462tcos (3t) | +0.047957tcos (5t) |
| -0.002159tcos (7t) | +0.000029tcos (9t) | +0.000002tcos (11t) |
| -0.000000tcos (13t) | -5.240757tsin (t) | +0.347041tsin (3t) |
| +0.020306tsin (5t) | -0.000669tsin (7t) | +0.000136tsin (9t) |
| -0.000001tsin (11t) | -0.000000tsin (13t) | +4.898568t2 cos (t) |
| -0.633927t2 cos (3t) | -0.017751t2 cos (5t) | +0.002044t2 cos (7t) |
| -0.000122t2 cos (9t) | +0.000000t2 cos (11t) | +5.910751t2 sin (t) |
| -0.423093t2 sin (3t) | +0.042183t2 sin (5t) | -0.003593t2 sin (7t) |
| -0.000015t2 sin (9t) | +0.000001t2 sin (11t) | -3.289085t3 cos (t) |
| +0.621827t3 cos (3t) | -0.037683t3 cos (5t) | +0.002208t3 cos (7t) |
| +0.000015t3 cos (9t) | -3.616263t3 sin (t) | -0.288847t3 sin (3t) |
| +0.002661t3 sin (5t) | +0.002305t3 sin (7t) | -0.000034t3 sin (9t) |
| +1.238058t4 cos (t) | -0.002983t4 cos (3t) | -0.002782t4 cos (5t) |
| -0.000822t4 cos (7t) | +1.664171t4 sin (t) | +0.252574t4 sin (3t) |
| -0.030580t4 sin (5t) | +0.000408t4 sin (7t) | -0.311598t5 cos (t) |
| -0.037610t5 cos (3t) | +0.009709t5 cos (5t) | -0.443419t5 sin (t) |
| +0.045286t5 sin (3t) | +0.002005t5 sin (5t) | +0.041874t6 cos (t) |
| -0.023862t6 cos (3t) | +0.045263t6 sin (t) | -0.022910t6 sin (3t) |
| +0.005561t7 cos (t) | +0.006555t7 sin (t) |
|
| |||
x8 = | +0.707697cos (t) | -0.737652cos (3t) | +0.029647cos (5t) |
| +0.000299cos (7t) | +0.000006cos (9t) | +0.000003cos (11t) |
| -0.000000cos (13t) | -0.000000cos (15t) | +17.859287sin (t) |
| +0.537553sin (3t) | +0.020622sin (5t) | -0.001769sin (7t) |
| +0.000063sin (9t) | -0.000002sin (11t) | -0.000000sin (13t) |
| +0.000000sin (15t) | -19.388600tcos (t) | -0.105739tcos (3t) |
| -0.072738tcos (5t) | +0.004113tcos (7t) | -0.000258tcos (9t) |
| +0.000002tcos (11t) | +0.000000tcos (13t) | -0.000000tcos (15t) |
| -16.767762tsin (t) | -0.906868tsin (3t) | +0.089594tsin (5t) |
| +0.000388tsin (7t) | -0.000052tsin (9t) | +0.000008tsin (11t) |
| -0.000000tsin (13t) | -0.000000tsin (15t) | +16.120089t2 cos (t) |
| +0.421931t2 cos (3t) | -0.092803t2 cos (5t) | +0.001418t2 cos (7t) |
| +0.000192t2 cos (9t) | -0.000007t2 cos (11t) | +0.000000t2 cos (13t) |
| +19.608050t2 sin (t) | -0.420458t2 sin (3t) | -0.128984t2 sin (5t) |
| +0.005901t2 sin (7t) | -0.000351t2 sin (9t) | -0.000001t2 sin (11t) |
| +0.000000t2 sin (13t) | -11.042356t3 cos (t) | +0.875871t3 cos (3t) |
| +0.101426t3 cos (5t) | -0.006986t3 cos (7t) | +0.000187t3 cos (9t) |
| +0.000001t3 cos (11t) | -12.224352t3 sin (t) | +0.801292t3 sin (3t) |
| -0.047084t3 sin (5t) | +0.002053t3 sin (7t) | +0.000194t3 sin (9t) |
| -0.000002t3 sin (11t) | +4.715957t4 cos (t) | -0.906645t4 cos (3t) |
| +0.033974t4 cos (5t) | -0.000385t4 cos (7t) | -0.000061t4 cos (9t) |
| +5.809016t4 sin (t) | +0.172731t4 sin (3t) | +0.039748t4 sin (5t) |
| -0.004671t4 sin (7t) | +0.000031t4 sin (9t) | -1.468495t5 cos (t) |
| +0.132440t5 cos (3t) | -0.017574t5 cos (5t) | +0.001171t5 cos (7t) |
| -1.830821t5 sin (t) | -0.280024t5 sin (3t) | +0.025814t5 sin (5t) |
| +0.000249t5 sin (7t) | +0.278969t6 cos (t) | +0.046180t6 cos (3t) |
| -0.006166t6 cos (5t) | +0.377937t6 sin (t) | -0.001896t6 sin (3t) |
| -0.006983t6 sin (5t) | -0.028889t7 cos (t) | +0.004692t7 cos (3t) |
| -0.029867t7 sin (t) | +0.021274t7 sin (3t) | -0.000816t8 cos (t) |
| -0.003373t8 sin (t) |
|
|
| |||
x9 = | +2.175721cos (t) | -2.094928cos (3t) | -0.084502cos (5t) |
| +0.003726cos (7t) | -0.000018cos (9t) | +0.000001cos (11t) |
| +0.000000cos (13t) | -0.000000cos (15t) | -0.000000cos (17t) |
| +64.332742sin (t) | +0.333348sin (3t) | +0.044040sin (5t) |
| +0.002384sin (7t) | -0.000104sin (9t) | +0.000006sin (11t) |
| -0.000000sin (13t) | -0.000000sin (15t) | +0.000000sin (17t) |
| -68.359632tcos (t) | +2.878461tcos (3t) | -0.079121tcos (5t) |
| -0.008828tcos (7t) | +0.000335tcos (9t) | -0.000022tcos (11t) |
| +0.000000tcos (13t) | +0.000000tcos (15t) | -0.000000tcos (17t) |
| -64.047243tsin (t) | -1.871954tsin (3t) | -0.188447tsin (5t) |
| +0.009868tsin (7t) | -0.000169tsin (9t) | -0.000004tsin (11t) |
| +0.000000tsin (13t) | -0.000000tsin (15t) | -0.000000tsin (17t) |
| +61.716104t2 cos (t) | -0.768222t2 cos (3t) | +0.295936t2 cos (5t) |
| -0.011360t2 cos (7t) | +0.000383t2 cos (9t) | +0.000014t2 cos (11t) |
| -0.000000t2 cos (13t) | +0.000000t2 cos (15t) | +69.633323t2 sin (t) |
| +2.498056t2 sin (3t) | -0.113525t2 sin (5t) | -0.012762t2 sin (7t) |
| +0.000674t2 sin (9t) | -0.000026t2 sin (11t) | -0.000000t2 sin (13t) |
| +0.000000t2 sin (15t) | -40.085044t3 cos (t) | -0.028980t3 cos (3t) |
| +0.047373t3 cos (5t) | +0.009129t3 cos (7t) | -0.000811t3 cos (9t) |
| +0.000013t3 cos (11t) | +0.000000t3 cos (13t) | -46.345112t3 sin (t) |
| +0.014214t3 sin (3t) | +0.328854t3 sin (5t) | -0.010719t3 sin (7t) |
| +0.000291t3 sin (9t) | +0.000013t3 sin (11t) | -0.000000t3 sin (13t) |
| +18.758768t4 cos (t) | -1.045383t4 cos (3t) | -0.199353t4 cos (5t) |
| +0.009593t4 cos (7t) | -0.000040t4 cos (9t) | -0.000004t4 cos (11t) |
| +21.493321t4 sin (t) | -0.827863t4 sin (3t) | -0.028074t4 sin (5t) |
| +0.005286t4 sin (7t) | -0.000478t4 sin (9t) | +0.000002t4 sin (11t) |
| -6.218301t5 cos (t) | +0.995799t5 cos (3t) | +0.007797t5 cos (5t) |
| -0.003375t5 cos (7t) | +0.000105t5 cos (9t) | -7.652998t5 sin (t) |
| +0.017785t5 sin (3t) | -0.070424t5 sin (5t) | +0.004954t5 sin (7t) |
| +0.000022t5 sin (9t) | +1.549574t6 cos (t) | -0.229644t6 cos (3t) |
| +0.028038t6 cos (5t) | -0.000898t6 cos (7t) | +1.852341t6 sin (t) |
| +0.213697t6 sin (3t) | -0.008050t6 sin (5t) | -0.001065t6 sin (7t) |
| -0.230798t7 cos (t) | -0.028635t7 cos (3t) | +0.000743t7 cos (5t) |
| -0.309074t7 sin (t) | -0.025324t7 sin (3t) | +0.008028t7 sin (5t) |
| +0.019906t8 cos (t) | +0.004745t8 cos (3t) | +0.017987t8 sin (t) |
| -0.012705t8 sin (3t) | -0.000299t9 cos (t) | +0.001329t9 sin (t) |
In Figure 1, we show a comparison of the solution computed through the symbolic method (to the fourth-order of the solution) presented in this paper, and the numeric solution to the problem calculated by a Runge-Kutta fourth-order method with a step of h=0.001 . Let us stress that the difference between both solutions at any time is smaller than 10-7 . The numerical solution has been calculated with the only purpose of validating the symbolic solution, as it is not the goal of this paper to develop a numerical tool.
Figure 1: Comparison between a fourth order symbolic and a numerical solution.
[figure omitted; refer to PDF]
In Table 3, we compare the CPU time of the algorithm implemented with MQSP and Maple running on an iMac 2.8 GHz Intel Core 2 Duo. For the test problem described by (7.1), we have calculated the solution up to order 14. In Table 3, we show the CPU time for the different orders of the solution from x0 (t) (nonperturbed problem) to x14 (t) . The computation time of x4 (t) with the use of the dsolve command of Maple in the implementation of the algorithm is larger than 100 000 seconds. This is due to the fact that dsolve is a general solver which handles different types of ordinary differential equations. We have also computed the symbolic solution taking into account the linearity of (7.1) and implementing a routine following the algorithm detailed in Section 4. These CPU times are given in the fourth column of Table 3. For orders of the solution larger than 11, Maple does not give any response. We can see that the difference between MQSP and Maple in terms of CPU time is significant.
Table 3: CPU time (in seconds) for the implementation with the MQSP, dsolve, and Maple.
xn | MQSP | Maple dsolve (using dsolve) | Maple |
x0 | 0.000187 | 0.059 | 0.042 |
x1 | 0.000847 | 0.131 | 0.078 |
x2 | 0.004637 | 0.196 | 0.147 |
x3 | 0.020543 | 0.300 | 0.299 |
x4 | 0.070672 | -- | 0.622 |
x5 | 0.209268 | -- | 1.330 |
x6 | 0.538049 | -- | 2.694 |
x7 | 1.255629 | -- | 5.426 |
x8 | 2.711981 | -- | 10.715 |
x9 | 5.421418 | -- | 21.052 |
x10 | 10.367602 | -- | 39.304 |
x11 | 18.494896 | -- | 73.471 |
x12 | 31.738276 | -- | -- |
x13 | 52.468991 | -- | -- |
x14 | 83.748298 | -- | -- |
8. Conclusions
In this paper, we have described an algorithm for the computation of the solution to a perturbed second-order differential equation through the asymptotic expansion technique. This algorithm has been implemented via a symbolic computation system which handles quasipolynomials.
[1] H. Poincaré Les Méthodes nouvelles de la Mécanique Céleste, Tome I , Gauthier-Villars, Paris, France, 1892.
[2] J. Henrard, "A survey of Poisson series processors," Celestial Mechanics , vol. 45, no. 1-3, pp. 245-253, 1988.
[3] A. Deprit, J. Henrard, A. Rom, "La théorie de la lune de Delaunay et son prolongement," Comptes Rendus de l'Academie des Sciences Paris A , vol. 271, pp. 519-520, 1970.
[4] A. Deprit, J. Henrard, A. Rom, "Analytical lunar ephemeris: Delaunay's theory," The Astronomical Journal , vol. 76, pp. 269-272, 1971.
[5] J. Henrard, "A new solution to the main problem of Lunar theory," Celestial Mechanics , vol. 19, no. 4, pp. 337-355, 1979.
[6] M. Chapront-Touzé, "La solution ELP du problème central de la Lune," Astronomy & Astrophysics , vol. 83, pp. 86, 1980.
[7] M. Chapront-Touzé, "Progress in the analytical theories for the orbital motion of the Moon," Celestial Mechanics , vol. 26, no. 1, pp. 53-62, 1982.
[8] J. F. Navarro, J. M. Ferrándiz, "A new symbolic processor for the Earth rotation theory," Celestial Mechanics & Dynamical Astronomy , vol. 82, no. 3, pp. 243-263, 2002.
[9] A. Abad, A. Elipe, J. F. San-Juan, S. Serrano, "Is symbolic integration better than numerical integration in satellite dynamics?," Applied Mathematics Letters , vol. 17, no. 1, pp. 59-63, 2004.
[10] V. I. Arnold Ordinary Differential Equations , The Massachusetts Institute of Technology, Cambridge, Mass, USA, 1973.
[11] J. F. Navarro, "On the implementation of the Poincaré-Lindstedt technique," Applied Mathematics and Computation , vol. 195, no. 1, pp. 183-189, 2008.
[12] J. F. Navarro, "Computation of periodic solutions in perturbed second-order ODEs," Applied Mathematics and Computation , vol. 202, no. 1, pp. 171-177, 2008.
[13] N. Minorsky Nonlinear Oscillations , Krieger, Boca Raton, Fla, USA, 1987.
[14] F. San-Juan, A. Abad, "Algebraic and symbolic manipulation of Poisson series," Journal of Symbolic Computation , vol. 32, no. 5, pp. 565-572, 2001.
[15] Y. Hardy, K. S. Tan, W. H. Steeb Computer Algebra with SymbolicC++ , World Scientific, Hackensack, NJ, USA, 2008.
[16] B. C. Ikoki, M. J. Richard, M. Bouazara, S. Datoussaïd, "Symbolic treatment for the equations of motion for rigid multibody systems," Transactions of the Canadian Society for Mechanical Engineering , vol. 34, no. 1, pp. 37-55, 2010.
[17] P. Hartman Ordinary Differential Equations , John Wiley & Sons Inc., New York, 1964.
You have requested "on-the-fly" machine translation of selected content from our databases. This functionality is provided solely for your convenience and is in no way intended to replace human translation. Show full disclaimer
Neither ProQuest nor its licensors make any representations or warranties with respect to the translations. The translations are automatically generated "AS IS" and "AS AVAILABLE" and are not retained in our systems. PROQUEST AND ITS LICENSORS SPECIFICALLY DISCLAIM ANY AND ALL EXPRESS OR IMPLIED WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES FOR AVAILABILITY, ACCURACY, TIMELINESS, COMPLETENESS, NON-INFRINGMENT, MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. Your use of the translations is subject to all use restrictions contained in your Electronic Products License Agreement and by using the translation functionality you agree to forgo any and all claims against ProQuest or its licensors for your use of the translation functionality and any output derived there from. Hide full disclaimer
Copyright © 2012 Juan F. Navarro. Juan F. Navarro et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Abstract
This paper describes an algorithm for implementing a perturbation method based on an asymptotic expansion of the solution to a second-order differential equation. We also introduce a new symbolic computation system which works with the so-called modified quasipolynomials, as well as an implementation of the algorithm on it.
You have requested "on-the-fly" machine translation of selected content from our databases. This functionality is provided solely for your convenience and is in no way intended to replace human translation. Show full disclaimer
Neither ProQuest nor its licensors make any representations or warranties with respect to the translations. The translations are automatically generated "AS IS" and "AS AVAILABLE" and are not retained in our systems. PROQUEST AND ITS LICENSORS SPECIFICALLY DISCLAIM ANY AND ALL EXPRESS OR IMPLIED WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES FOR AVAILABILITY, ACCURACY, TIMELINESS, COMPLETENESS, NON-INFRINGMENT, MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. Your use of the translations is subject to all use restrictions contained in your Electronic Products License Agreement and by using the translation functionality you agree to forgo any and all claims against ProQuest or its licensors for your use of the translation functionality and any output derived there from. Hide full disclaimer