\documentstyle[11pt,reduce,fancyheadings]{article}
\title{
\author{Neil Langmead \\
Konrad-Zuse-Zentrum f\"ur Informationstechnik (ZIB) \\
Takustrasse 7 \\
D- 14195 Berlin Dahlem \\
Berlin
Germany}
\date{January 1997}
\def\foottitle{Rational Integration in \small{REDUCE}}
\pagestyle{fancy}
\lhead[]{{\footnotesize\leftmark}{}}
\rhead[]{\thepage}
\setlength{\headrulewidth}{0.6pt}
\setlength{\footrulewidth}{0.6pt}
\addtolength{\oddsidemargin}{-20 mm}
\addtolength{\textwidth}{25 mm}
\pagestyle{fancy}
\setlength{\headrulewidth}{0.6pt}
\setlength{\footrulewidth}{0.6pt}
\setlength{\topmargin}{1 mm}
\setlength{\footskip}{10 mm}
\setlength{\textheight}{220 mm}
\cfoot{}
\rfoot{\small\foottitle}
\def\exprlist {exp$_{1}$,exp$_{2}$, \ldots ,exp$_{{\tt n}}$}
\def\lineqlist {lin\_eqn$_{1}$,lin\_eqn$_{2}$, \ldots ,lin\_eqn$_{n}$}
\pagestyle{fancy}
\begin{document}
\maketitle{
\begin{center} \Large Package to Integrate Rational Functions using the minimal algebraic extension to the constant field \end{center}}
\pagebreak
\tableofcontents
\pagebreak
\section{Rational Integration}
\normalsize
This package implements the Horowitz/ Rothstein/ Trager algorithms ~\cite{Ged92} for the integration of rational functions in \small{REDUCE}.\normalsize We work within a field $K$ of characteristic $0$ and functions $p,q \in K[x]$. $K$ is normally the field $Q$ of rational numbers, but not always. These procedures return $\int \frac{p}{q} dx.$
The aim is to be able to integrate any function of the form $p/q$ in $x$, where $p$ and $q$ are polynomials in the field $Q$. The algorithms used avoid algebraic number extensions wherever possible, and in general, express the integral using the minimal algebraic extension field. \\
\subsection{Syntax of $ratint$}
This function has the following syntax:
\begin{center} \bf{ratint(p,q,var)} \end{center}
where $ p/q$ is a rational function in $var$. The output of $ratint$ is a list of two elements: the first is the polynomial part of the integral, the second is the logarithmic part. The integral is the sum of these parts.
\subsection{Examples}
consider the following examples in \small{REDUCE}:
\begin{verbatim}
ratint(1,x^2-2,x);
sqrt(2)*x-2 sqrt(2)*x+2
log(-------------) - log(-------------)
sqrt(2) sqrt(2)
{ 0, --------------------------------------- }
2*sqrt(2)
p:=441*x^7+780*x^6-2861*x^5+4085*x^4+7695*x^3+3713*x^2-43253*x
+24500;
q:=9*x^6+6*x^5-65*x^4+20*x^3+135*x^2-154*x+49;
ratint(p,q,x);
49 6 226 5 268 4 1332 3 2809 2 752 256
---*(x + ---*x - ---*x + ----*x - ----*x - ---*x + ---)
2 147 49 49 147 21 9
{----------------------------------------------------------- , 0 }
4 2 3 2 7
x - ---*x - 4*x + 6*x - ---
3 3
k:=36*x^6+126*x^5+183*x^4+(13807/6)*x^3-407*x^2-(3242/5)*x+(3044/15);
l:=(x^2+(7/6)*x+(1/3))^2*(x-(2/5))^3;
ratint(k,l,x);
5271 3 39547 2 31018 7142
------*(x + -------*x - -------*x + -------)
5 52710 26355 26355
{------------------------------------------------,
4 11 3 11 2 2 4
x + ----*x - ----*x - ----*x + ----
30 25 25 75
37451 2 91125 2 128000 1
-------*(log(x - ---) + -------*log(x + ---) - --------*log(x + ---))}
16 5 37451 3 37451 2
ratint(1,x^2+1,x);
2 1
{0,log_sum(beta,beta + ---,0,log(2*beta*x - 1)*beta)}
4
\end{verbatim}
The meaning of the log\_sum function will be explained later.
\pagebreak
\section{The Algorithm}
The following main algorithm is used:
procedure $ratint(p,q,x);$
% p and q are polynomials in $x$, with coefficients in the %constant field Q
solution\_list $\leftarrow HorowitzReduction(p,q,x)$
\\
$c/d \leftarrow$ part(solution\_list,1)\\
$poly\_part \leftarrow$ part(solution\_list,2)
\\
$rat\_part \leftarrow$ part(solution\_list,3)
\\
$rat\_part \leftarrow LogarithmicPartIntegral(rat\_part,x) $
\\
return($rat\_part+c/d +poly\_part$)
\\
end
The algorithm contains two subroutines, $HorowitzReduction$ and $rt$. $HorowitzReduction$ is an implementation of Horowitz' method to reduce a given rational function into a polynomial part and a logarithmic part. The integration of the polynomial part is a trivial task, and is done by the $int$ operator in \small{REDUCE}. The integration of the logarithmic part is done by the routine $rt$, which is an impementation of the Rothstein and Trager method. These two answers are outputed in a list, the complete answer being the sum of these two parts.
\\
These two algorithms are as follows:
procedure $how(p,q,x)$
for a given rational function $p/q$ in $x$, this algorithm calculates the
reduction of $ \int(p/q)$ into a polynomial part and logarithmic part. \\
$ poly\_part \leftarrow quo(p,q); \hspace{3 mm} p \leftarrow rem(p,q)$;
$d \leftarrow GCD(q,q') $; \hspace{3 mm} $b \leftarrow quo(q,d)$; \hspace{3 mm}
$m \leftarrow deg(b)$; \\ $n \leftarrow deg(d)$;
$a \leftarrow \sum_{i=1}^{m-1} a_{i}x^{i}$; \hspace{3 mm}
$ c \leftarrow \sum_{i=1}^{n-1} c_{i}x^{i}$; \\
$r \leftarrow b*c'-quo(b*d',d)+d*a; $\\
\begin{tabbing}
for $i$ from \= $0$ \= to $m+n-1$ do \\
\> \{ \\
\> \> $ eqns(i) \leftarrow coeff(p,i)=coeff(r,i)$; \\
\> \};
\end{tabbing}
$solve(eqns,\{a(0),....,a(m-1),c(0),....,c(n-1)\});$
return($c/d+\int poly\_part + \int a/b$);
\\
end;
\newpage
procedure RothsteinTrager($a,b,x$)
\% Given a rational function $a/b$ in $x$ with $deg(a)<deg(b)$, with b monic and square free, we calculate $ \int(a/b)$ \\
$R(z) \leftarrow residue(a-zb',b)$ \\
$(r_{1}(z)...r_{k}(z)) \leftarrow factors(R(z))$ \\
integral $\leftarrow 0 $\\
\begin{tabbing}
for $i$ from $1$ \= to \= $k$ do \\
\> \{ \\
\> \> $ d \leftarrow degree(r_{i}(z))$ \\
\> \> if $d=1$ then \= \\
\> \> \> \{ \\
\> \> \> c $\leftarrow solve(r_{i}(z)=0,z)$ \\
\> \> \> v $\leftarrow GCD(a-cb',b)$\\
\> \> \> v $\leftarrow v/lcoeff(v) $\\
\> \> \> $integral \leftarrow integral+c*log(v)$ \\
\> \> \> \}\\
\> \> else \= \{ \\
\> \> \> \% we need to do a GCD over algebraic number field\\
\> \> \> v $\leftarrow GCD(a-\alpha*b',b) $ \\
\> \> \> v $\leftarrow v/lcoff(v) $, \hspace{3 mm}
where $\alpha=roof\_of(r_{i}(z)) $\\
\> \> if d=2 then \= \{ \\
\> \> \> \% give answer in terms of radicals \\
\> \> \> c $\leftarrow solve(r_{i}(z)=0,z) $ \\
\> \> \> for j from 1 to 2 do \= \{ \\
\> \> \> $v[j] \leftarrow substitute(\alpha=c[j],v) $ \\
\> \> \> $integral \leftarrow integral+c[j]*log(v[j]) $ \\
\> \> \> \> \> \} \\
\> \> \> \} \\
\> \> \> else \= \{ \\
\> \> \> \% Need answer in terms of root\_of notation \\
\> \> \> for j from 1 to d do \= \{ \\
\> \> \> v[j] $\leftarrow substitute(\alpha=c[j],v) $ \\
\> \> \> integral $ \leftarrow integral+c[j]*log(v[j]) $ \\
\> \> \> \% where $c[j]=root\_of(r_{i}(z))$ \} \\
\> \> \> \> \> \> \} \\
\> \> \> \} \\
\> \> \} \\
\> \} \\
return(integral) \\
end
\end{tabbing}
\pagebreak
\section{The log\_sum operator}
The algorithms above returns a sum of terms of the form
\[ \sum_{\alpha \mid R(\alpha)=0} \log(S(\alpha,x)), \]
where $R \in K[z]$ is square free, and $S \in K[z,x]$. In the cases where the degree of $R(\alpha)$ is less than two, this is merely a sum of logarithms. For cases where the degree is two or more, I have chosen to adopt this notation as the answer to the original problem of integrating the rational function. For example,
consider the integral
\[ \int \frac{a}{b}=\int \frac{2x^5-19x^4+60x^3-159+x^2+50x+11}{x^6-13x^5+58x^4-85x^3-66x^2-17x+1}\, dx \]
Calculating the resultant $R(z)=res_x(a-zb',b)$ and factorising gives
\[ R(z)=-190107645728000(z^3-z^2+z+1)^{2} \]
Making the result monic, we have
\[ R_2(z)=z^3-z^2+z+1 \]
which does not split over the constant field $Q$.
Continuting with the Rothstein Trager algorithm, we now calculate
\[ gcd(a-\alpha\,b',b)=z^2+(2*\alpha-5)*z+\alpha^2, \] where $\alpha$ is a root of $R_2(z)$. \\
Thus we can write
\[ \int \frac{a}{b}= \sum_{\alpha \mid \alpha^3-\alpha^2+\alpha+1=0} \alpha*\log(x^2+2\alpha x-5x+\alpha^2), \]
and this is the answer now returned by \small{REDUCE}, via a function called $log\_sum$. This has the following syntax:
\begin{center}$ log\_sum(\alpha,eqn(\alpha),0,sum\_term,var)$ \end{center}
where $\alpha$ satisfies $eqn=0$, and $sum\_term$ is the term of the summation in the variable $var$. Thus in the above example, we have
\[ \int \frac{a}{b}\,dx= log\_sum(\alpha,\alpha^3-\alpha^2+\alpha+1,0,\alpha*\log(x^2+2\alpha x-5x+\alpha^2),x) \]
Many rational functions that could not be integrated by \small{REDUCE} previously can now be integrated with this package. The above is one example; some more are given on the next page.
\pagebreak
\subsection{More examples}
\begin{eqnarray*}
\int \frac{1}{x^5+1} \, dx & = &\frac{1}{5}\log(x + 1) \\
& & \mbox{} + 5log\_sum(\beta,\beta^4+\frac{1}{5}\beta^3+\frac{1}{25}\beta^2+\frac{1}{125}\beta+\frac{1}{625},0,\log(5*\beta+x)*\beta)
\end{eqnarray*}
which should be read as
\[
\int \frac{1}{x^5+1}\,dx = \frac{1}{5}\log(x+1)+\sum_{\beta \mid \beta^4+\frac{1}{5}\beta^3+\frac{1}{25}\beta^2+\frac{1}{125}\beta+\frac{1}{625}=0} \log(5*\beta+x)\beta \]
\vspace{5 mm}
\begin{eqnarray*}
\lefteqn{\int \frac{7x^{13}+10x^8+4x^7-7x^6-4x^3 -4x^2+3x+3}{x^{14}-2x^8-2x^7-2x^4-4x^3-x^2+2x+1} \, dx =} \\
& & log\_sum(\alpha,\alpha^2 -\alpha -\frac{1}{4},0,log( - 2\alpha x^2 - 2\alpha x + x^7 + x^2 - 1)*\alpha,x) ,
\end{eqnarray*}
\[ \int \frac{1}{x^3+x+1} \, dx = log\_sum(\beta,\beta^3-\frac{3}{31}\beta^2-\frac{1}{31},0,\beta \log(-\frac{62}{9}\beta^2+\frac{31}{9} \beta +x+\frac{4}{9})). \]
\section{Options}
There are several alternative forms that the answer to the integration problem can take. One output is the $log\_sum$ form shown in the examples above. There is an option with this package to convert this to a "normal" sum of logarithms in the case when the degree of $eqn$ in $\alpha$ is two, and $\alpha$ can be expressed in surds. To do this, use the function $convert$, which has the following syntax:
\begin{center} convert(exp) \end{center}
If exp is free of $log\_sum$ terms, then $exp$ itself is returned. If $exp$ contains $log\_sum$ terms, then $\alpha$ is represented as surds, and substituted into the $log\_sum$ expression. For example, using the last example, we have in \small{REDUCE}:
\begin{verbatim}
2: ratint(a,b,x);
{0,
2 1
log_sum(alpha,alpha - alpha - ---,0,
4
2 7 2
log( - 2*alpha*x - 2*alpha*x + x + x - 1)*alpha,x)}
3: convert(ws);
1 2 7
---*(sqrt(2)*log( - sqrt(2)*x - sqrt(2)*x + x - x - 1)
2
2 7
- sqrt(2)*log(sqrt(2)*x + sqrt(2)*x + x - x - 1)
2 7
+ log( - sqrt(2)*x - sqrt(2)*x + x - x - 1)
2 7
+ log(sqrt(2)*x + sqrt(2)*x + x - x - 1))
\end{verbatim}
\subsection{LogtoAtan function}
The user could then combine these to form a more elegant answer, using the switch combinelogs if one so wished. Another option is to convert complex logarithms to real arctangents ~\cite{Bron97}, which is recommended if definite integration is the goal. This is implemented in \small{REDUCE} via a function $convert\_log$, which has the following syntax:
\begin{center} \bf{convert\_log(exp)}, \end{center}
where $exp$ is any log\_sum expression.\\
The procedure to convert complex logarithms to real arctangents is based on an algorithm by Rioboo. Here is what it does: \\
Given a field $K$ of characteristic 0 such that $\sqrt(-1) \not\in K$ and
$A, B \in K[x]$ with $B \not = 0$, return a sum $f$ of arctangents of polynomials in $K[x]$ such that
\[ \frac{df}{dx}=\frac{d}{dx} i \log(\frac{A+ i B}{A- i B}) \]
Example:
\[ \int \frac{x^4-3*x^2+6}{x^6-5*x^4+5*x^2+4} \, dx = \sum_{ \alpha \mid 4\alpha+1=0} \alpha \log(x^3+2\alpha x^2-3 x-4 \alpha) \]
Substituting $\alpha=i/2$ and $\alpha=-i/2$ gives the result
\[ \frac{i}{2} \log(\frac{(x^3-3 x)+i (x^2-2)}{(x^3-3 x)-i (x^2-2)}) \]
Applying logtoAtan now with $A=x^3-3 x$, and $B=x^2-2$ we obtain
\[ \int \frac{x^4-3*x^2+6}{x^6-5*x^4+5*x^2+4} \, dx = \arctan(\frac{x^5-3 x^3+x}{2})+\arctan(x^3)+\arctan(x) , \]
and this is the formula which should be used for definite integration. \\
Another example in \small{REDUCE} is given below:
\begin{verbatim}
1: ratint(1,x^2+1,x);
*** Domain mode rational changed to arnum
2 1
{0,log_sum(beta,beta + ---,0,log(2*beta*x - 1)*beta)}
4
13: part(ws,2);
2 1
log_sum(beta,beta + ---,0,log(2*beta*x - 1)*beta)
4
14: on combinelogs;
15: convertlog(ws);
1 - i*x + 1
---*log(------------)*i
2 i*x + 1
logtoAtan(-x,1,x);
2*atan(x)
\end{verbatim}
\section{Hermite's method}
The package also implements Hermite's method to reduce the integral into its polynomial and logarithmic parts, but occasionally, \small{REDUCE} returns the incorrect answer when this algorithm is used. This is due to the REDUCE operator pf, which performs a complete partial fraction expansion when given a rational function as input. Work is presently being done to give the pf operator a facility which tells it that the input is already factored. This would then enable REDUCE to perform a partial fraction decomposition with respect to a square free denominator, which may not necessarily be fully factored over Q.
\newline
For a complete explanation of this and the other algorithms used in this package, including the theoretical justification and proofs, please consult ~\cite{Ged92}.
\section{Tracing the $ratint$ program}
The package includes a facility to trace in some detail the inner workings of the $ratint$ program. Messages are given at the key stages of the algorithm, together with the results obtained. These messages are displayed when the switch $traceratint$ is on, which is done in \small{REDUCE} \normalsize with the command
\begin{verbatim}
on traceratint;
\end{verbatim}
This switch is off by default. Here is an example of the output obtained with this switch on:
\begin{verbatim}
Loading image file: /silo/tony/red/lisp/psl/solaris/red/reduce.img
REDUCE Development Version, 21-May-97 ...
1: load_package ratint;
2: on traceratint;
3: ratint(1+x,x^2-2*x+1,x);
x + 1
performing Howoritz reduction on --------------
2
x - 2*x + 1
- 2 1
Howoritz gives: {-------,0,-------}
x - 1 x - 1
1
computing Rothstein Trager on -------
x - 1
integral in Rothstein T is log(x - 1)
- 2
{-------,log(x - 1)}
x - 1
\end{verbatim}
\section{Bugs, suggestions and comments}
This package was written when the author was working as a placement student at ZIB Berlin. All comments should therefore be reported to Winfried Neun, ZIB, Takustrasse 7, D 14195 Berlin Dahlem, Germany \\ (email: neun@zib.de).
\pagebreak
\begin{thebibliography}{999999}
\normalsize
\bibitem[Bron97]{Bron97} Bronstein, Manuel,
{\it Symbolic Integration I: Transendental Functions},
Springer-Verlag, Heidelberg, 1997.
\bibitem[Dav88]{Dav88} Davenport, James H. et al,
{\it Computer Algebra- Systems and Algorithms for Algebraic Computation},
Academic Press, 1988.
\bibitem[Ged92]{Ged92} Geddes, K.O. et al,
{\it Algorithms for Computer Algebra}, Klewer Academic \mbox{Publishers}, 1992.
\bibitem[Red36]{Red36} Hearn, Anthony C. and Fitch, John F.
{\it REDUCE User's Manual 3.6}, RAND Corporation, 1995
\end{thebibliography}
\end{document}