Wolfram：New Methods for Computing Algebraic Integrals

 

In 2004, I became obsessed with computing integrals, using both elementary techniques known to calculus students and the Risch algorithm and its extensions by Davenport, Trager and Bronstein [1, 2, 3, 4, 5]. Taking inspiration from George Beck’s step-by-step differentiation program, WalkD, I decided to write a step-by-step program for computing indefinite integrals.

For example, it’s easy to write a simple set of replacement rules in Mathematica for integrating (expanded) polynomials:

 

 

 

 

 

Using these rules, we can now integrate polynomials, for example：

 

 

 

 

 

 

 

 

 

After spending far too many late nights entering integration techniques as pattern-matching rules into Mathematica, I had the code at a reasonable state and I sent it to Wolfram Research for possible inclusion in Mathematica. Soon after, I was offered a job and ended up working for Wolfram for several years, predominantly on Wolfram|Alpha. My step-by-step integrator is still computing many integrals, some of which I have most likely forgotten how to do myself.

 

As an aside, the idea of using rule-based programming to compute indefinite integrals dates back to 1961, with the Symbolic Automatic Integrator (SAINT) by James Slagle at MIT. SAINT could solve “symbolic integration problems approximately at the level of a good college freshman and, in fact, uses many of the same methods (including heuristics) used by a freshman” [6]. The step-by-step integrator I wrote used around 350 rules and could integrate more than 99% of integrals in calculus textbooks. Since then, Albert Rich used Mathematica to create the Rule-based Integrator (Rubi), which uses over 6,700 rules.

 

Fast forward to 2020 and I hadn’t looked at integrals for a decade. I decided to see what advances had been made in the last 10 or so years. I was particularly interested in a Gröbner basis–based algorithm developed by Manuel Kauers and its extensions by Brian Miller that could seemingly outperform the algebraic case of the Risch algorithm in the AXIOM computer algebra system on many integrals [7, 8]. For example:

 

 

 

It’s trivial to check the result:

 

 

Once again, I quickly became hooked on integrals, or more specifically, algorithmic solutions to indefinite integrals.

When I last looked at symbolic integration, I was interested in the transcendental case of the Risch algorithm, of which Mathematica has a near-complete implementation. For example, the following is a simple integral for the transcendental case of the Risch algorithm:

 

 

 

Once again, I quickly became hooked on integrals, or more specifically, algorithmic solutions to indefinite integrals.

When I last looked at symbolic integration, I was interested in the transcendental case of the Risch algorithm, of which Mathematica has a near-complete implementation. For example, the following is a simple integral for the transcendental case of the Risch algorithm:

 

If we’re happy to play it fast and loose with branch cuts, then we can write this integral as:

 

For this integral, we can substitute ，and we get:

Substituting back, we get:

 

This answer has branch cut issues; however, we can fix this by writingasThen we have the correct antiderivative:

 

 

 

I wondered if this method of using a Laurent polynomial substitution to simplify an algebraic integral was just a trick that worked for this integral or a hint to a more general method. It turns out this trick works for many integrals; for example, the integral we tried previously on Kauers’s algorithm

 

 

 

can be reduced to

 

 

with the substitution u = x⁴ + x^–4. Once corrected for branch cut issues, the solution is given by:

 

 

 

A general method would seek a substitution of the form u =such that

 

 

 

where R1(u), R2(u) are rational functions of u and ，are undetermined coefficients.

We start by using SolveAlways to compute the undetermined coefficients in the u substitution:

 

 

 

So we have a candidate substitution that fits the radicand part of the integrand. Does this substitution fit the rest of the integrand? We can compute this as follows:

 

We have made the same simplification to the integral that we made by hand previously, namely

 

 

 

where u = .

 

This method is implemented with IntegrateAlgebraic at the Wolfram Function Repository. (In 2020, I further investigated the computation of pseudo-elliptic integrals in terms of elementary functions [9].) Given the simplicity of this method, it can integrate a wide range of integrals.

Here are some examples:

 

 

 

 

 

 

 

 

 

 

 

 

 

Unlike the algebraic case of the Risch algorithm, this technique can quickly solve many integrals involving parameters:

 

 

 

 

 

 

 

 

 

 

 

 

What about the following integral?

 

 

Unlike the previous examples, this integral is not solvable with a Laurent polynomial substitution.

In 1882, Günther developed a method for computing some otherwise difficult algebraic integrals [10]. Given the integral

 

 

 

where p(x) and q(x) are polynomials in x, Günther made the substitution

 

 

 

such that the integral becomes

 

 

 

where s(x) = v_Q–P+R–1 x^Q–P+R–1 + v_Q–P+R–2 x^Q–P+R–2 + … + v1 x + v0, where P = degx(p(x)), Q = degx(q(x)), R = degx(r(x)) and v0, v1, …, v_Q–P+R–1, c0, c1,c2, c3 are undetermined coefficients and R(un) is a rational function in uⁿ with undetermined coefficients.

We can use Günther’s method to solve this integral in Mathematica as follows. The substitution is of the form:

 

 

 

And we assume the integrand in u is of the form:

 

 

 

Then the integrand in x is given by:

 

 

 

Now we need to solve for the undetermined coefficients in the substitution (v0, v1, v2) and in the rational integrand (w0, w1):

 

 

We can substitute this solution into our integrand in u and substitution:

 

 

 

Now we can use Integrate to compute the resulting integral:

 

 

 

Then substitute back to solve our original integral:

 

 

 

A quick check that our solution is correct:

 

 

 

A generalized version of Günther’s method is implemented in IntegrateAlgebraic. This method can solve many otherwise difficult integrals. Here are some examples:

 

 

 

 

 

 

 

 

 

 

 

 

 

This method also handles integrals containing parameters:

 

 

 

 

 

If we combine this method with integrating term by term after expanding the integrand into a sum of terms, we can handle more exotic algebraic integrals:

 

 

Combining the Laurent polynomial substitution method with the generalized Günther method and integrating term by term allows us to compute even more complex integrals:

 

 

 

In this case, we wrote the integral as:

 

 

 

Then the integralThen the integralwith the substitution u =  = 1 – x ³，while the integralwas reduced towith the substitution s  = 。

Integrating expressions containing nested radicals has always been a tricky business. A well-known example is

 

 

 

 

which can be computed using the substitution We can make this substitution with GroebnerBasis as follows:

 

 

 

We need to express this relationship in terms of Dt[y]:

 

 

 

We can now integrate the rational function of u and substitute back for x:

 

 

This method can solve much more difficult integrals involving nested radicals. For example:

 

 

 

 

 

A generalization of this approach is used within IntegrateAlgebraic. Here are some challenging examples:

 

 

 

 

 

 

Like the other methods in IntegrateAlgebraic, we readily handle integrals involving parameters:

 

 

 

All of the integrals in this post contain polynomial radicands; however, these methods generalize to rational radicands. For example:

 

 

 

 

 

There are still many algebraic integrals that these methods will not compute. For example, the following integral possesses an elementary (albeit enormous) solution: 

 

 

However, compared to the algebraic case of the Risch–Trager–Bronstein algorithm, which is not completely implemented in any computer algebra system, these methods are fast, simple and complement the existing integration capabilities of Mathematica’s Integrate function. We are currently considering including IntegrateAlgebraic within Integrate in an upcoming release.