From: rusin@vesuvius.math.niu.edu (Dave Rusin) Newsgroups: sci.math.symbolic Subject: Re: how can i find information on symbolic integration Date: 4 Mar 1999 23:02:31 GMT Richard J. Fateman wrote: >In article <36DE14C2.E31947B7@lynx.neu.edu>, >Mark Logan wrote: >>i'm wondering how hard it is to write a program to integrate elementary >>functions symbolically. If you can point me to a place where i can find [...] >There is no particular reason for you to write a program, of course. >You can get one free, or buy one. I know that algorithms exist and I understand the basic framework, but I've never looked in detail at how they operate or how concrete they are. So I am curious about the extent to which the theoretical procedures can be, and have been, implemented. For example, can any program give a symbolic antiderivative in the following example? I ask because the underlying question here is not really integration but algebraic geometry, and I don't even know if software to do _that_ problem is available. Here goes: in Maple format, let eq:=-x^2-6*x+3+(-9-3*x^2)*y+(x^2+3)*y^3; 2 2 2 3 -x - 6 x + 3 + (-9 - 3 x ) y + (x + 3) y This is a cubic in y, so we can "explicitly" solve for y in terms of x; I've rigged it up so the discriminant is a perfect square (of 9*(x+1)(x-3)(x^2+3), say), so y has the "simple" form R1^1/3 + R2^1/3 for two rational functions R1, R2 in x. I want to know the integral of this function y with respect to x. Now on the other hand, I know this integral is do-able. More precisely, the given polynomial eq relating x and y describes a curve in the x-y plane which happens to have genus 0, [this is what happens when you hand an algebraic geometer a multivariable polynomial...] which means it is possible to parameterize the curve with rational functions in t. This is clear if you try to solve the original equation for x in terms of y; a quadratic in y must be square, so you parameterize y with a rational function of t, and then get x also in terms of t. Thus \integral y dx = \integral y(t)x'(t) dt is the integral of a rational function of t, and so may be antidifferentiated with the help of partial fractions and so on. [Can you tell what I'm teaching this semester?] The parameter t may conversely be expressed as a rational function of x and y, and so by substituting y = R1^1/3+R2^1/3, we can get an explicit antiderivative of x as a complicated function of x. There's no magic or deep thought here, beyond the ability to recognize and then parameterize a curve of genus zero, but for curves with singular points it's a challenge even to recognize the genus, much less to provide a birational equivalence with a simpler curve. Now, I don't _know_ that it's necessary to be able to do so in order to compute integrals, but I would bet that some simple combination of examples would (perhaps by the Riemann-Roch theorem?) allow us to obtain a parameterization of the curve if we could antidifferentiate enough things. Thus it seems to me that an effective procedure to perform integration would have to be at least as powerful as a parameterization tool. So how 'bout it: does any existing code perform this integral, even assuming plenty of machine resources? What about the general problem of applying algebraic geometry to integration problems? dave ============================================================================== From: Helmut Kahovec Newsgroups: sci.math.symbolic Subject: Re: how can i find information on symbolic integration Date: Fri, 05 Mar 1999 00:39:02 GMT Dave Rusin wrote: > > I know that algorithms exist and I understand the basic framework, but I've > never looked in detail at how they operate or how concrete they are. So I am > curious about the extent to which the theoretical procedures can be, and > have been, implemented. For example, can any program give a symbolic > antiderivative in the following example? I ask because the underlying > question here is not really integration but algebraic geometry, and I > don't even know if software to do _that_ problem is available. > > Here goes: in Maple format, let eq:=-x^2-6*x+3+(-9-3*x^2)*y+(x^2+3)*y^3; > > 2 2 2 3 > -x - 6 x + 3 + (-9 - 3 x ) y + (x + 3) y > > This is a cubic in y, so we can "explicitly" solve for y in terms of > x; I've rigged it up so the discriminant is a perfect square (of > 9*(x+1)(x-3)(x^2+3), say), so y has the "simple" form R1^1/3 + R2^1/3 > for two rational functions R1, R2 in x. > > I want to know the integral of this function y with respect to x. In either Release 4 or Release 5 you may compute the antiderivative of y w.r.t. x as follows. If you set the information level for int() to a higher value than the default then you can also see what is going on: > restart; Algebraic numbers are represented as RootOf() expressions in Maple: > alias(alpha=RootOf(-x^2-6*x+3+(-9-3*x^2)*y+(x^2+3)*y^3,y)); I, alpha You may also try a smaller value than 5 if this gives you too much information: > infolevel[int]:=5; infolevel[int] := 5 > int(alpha,x); int/indef1: first-stage indefinite integration int/indef2: second-stage indefinite integration int/rischnorm: enter Risch-Norman integrator int/rischnorm: exit Risch-Norman integrator int/algrisch/int: Risch/Trager's algorithm for algebraic function int/algrisch/int: entering integrator at time 1.432 int/algrisch/int: function field has degree 3 int/algrisch/int: pol min is (1+3*_X^2)*y^3+(-9*_X^2-3)*y-1-6*_X+3*_X^2 integrand is -y/_X^2 int/algrisch/int: computation of an integral basis: start time 1.633 int/algrisch/int: integral basis is array(1 .. 3, 1 .. 3,[(2, 2)=1/3+_X^2,(1, 2)=0,(1, 1)=1,(2, 1)=0,(2, 3)=-1/2*_X*(5+15*_X^2+3*_X+9*_X^3)/(_X+1)/(3*_X-1),(3, 2)=0,(1, 3)=4/3*(7*_X-5+6*_X^2)/(_X+1)/(3*_X-1),(3, 1)=0,(3, 3)=(1+3*_X^2)/(_X+1)/(3*_X-1)]) int/algrisch/int: computation of an integral basis: end time 3.165 int/algrisch/int: normalization at infinity: start time 3.165 int/algrisch/int: normal integral basis is array(1 .. 3, 1 .. 3,[(2, 2)=-1/2*(3*_X^2+3*_X+1+9*_X^3)/(_X+1)/(3*_X-1),(1, 2)=4/3*(7*_X-5+6*_X^2)/(_X+1)/(3*_X-1),(1, 1)=1,(2, 1)=0,(2, 3)=-2*_X*(1+3*_X^2)/(_X+1)/(3*_X-1),(3, 2)=(1+3*_X^2)/(_X+1)/(3*_X-1),(1, 3)=-4/3*(_X^2-12*_X+6*_X^3+5)/(_X+1)/(3*_X-1),(3, 1)=0,(3, 3)=-(_X+3*_X^3-1-3*_X^2)/(_X+1)/(3*_X-1)]) int/algrisch/int: discriminant of the basis is (1+3*_X^2)^2 int/algrisch/int: normalization at infinity: end time 3.605 int/algrisch/int: genus of the function field 0 int/algrisch/int: derivatives are array(1 .. 3, 1 .. 3,[(2, 2)=9*_X-3,(1, 2)=0,(1, 1)=0,(2, 1)=-28*_X,(2, 3)=3,(3, 2)=-4,(1, 3)=0,(3, 1)=-28*_X-28/3,(3, 3)=9*_X+3]) int/algrisch/int: denominator of the derivatives is 3+9*_X^2 int/algrisch/int: denominator of the integrand is _X^2+3*_X^4 int/algrisch/int: computation of the algebraic part: start time 4.065 int/algrisch/algpart: reducing factor of multplicity 2 int/algrisch/solvemodv1: solving a linear system modulo _X int/algrisch/solvemodv1: solutions are [0, 3, -3] int/algrisch/int: computation of the algebraic part: end time 4.125 int/algrisch/int: coordinates of the algebraic part are [0, 1/_X, -1/_X] int/algrisch/int: coordinates of the transcendental part are [-14, 3, 0] int/algrisch/int: computation of the transcendental part: start time 4.145 int/algrisch/transcpar: computing a basis for the residues at time 4.336 int/algrisch/residues: computing a splitting field at time 4.356 int/algrisch/residues: polynomial has degree 3 int/algrisch/residues: splitting field computed at time 4.887 int/algrisch/residues: number of distinct residues 3 int/algrisch/transcpar: basis for the residues computed at time 5.107 int/algrisch/transcpar: dimension is 2 int/algrisch/transcpar: Z-basis for the residues array(1 .. 2,[(1)=1+RootOf(-6*_Z+24*_Z^3-1)-6*RootOf(-6*_Z+24*_Z^3-1)^2,(2)=-1-2*RootOf(-6*_Z+24*_Z^3-1)+6*RootOf(-6*_Z+24*_Z^3-1)^2]) int/algrisch/transcpar: building divisors at time 7.961 int/algrisch/transcpar: exponents are [[1, 0, -1], [0, 1, -1]] int/algrisch/transcpar: testing divisors for principality at time 11.857 int/algrisch/areprinc: the order of the divisor must be 1 int/algrisch/areprinc: the divisor is principal: time 33.818 int/algrisch/areprinc: the order of the divisor must be 1 int/algrisch/areprinc: the divisor is principal: time 56.912 int/algrisch/transcpar: divisors proven principal at time 56.912 int/algrisch/transcpar: orders are [1, 1] int/algrisch/int: computation of the transcendental part: end time 56.912 int/algrisch/int: exiting integrator for algebraic functions at time 56.992 2 x x - 56/3 ---------------- + 40/3 ---------------- (x + 1) (-3 + x) (x + 1) (-3 + x) 2 1 alpha x - 16 ---------------- - ---------------- (x + 1) (-3 + x) (x + 1) (-3 + x) 3 alpha alpha x alpha x - 3 ---------------- + ---------------- + 3 ---------------- (x + 1) (-3 + x) (x + 1) (-3 + x) (x + 1) (-3 + x) 2 2 2 alpha x alpha - 2 ---------------- - 6 ---------------- + 4 %6 + 4 %6 %1 (x + 1) (-3 + x) (x + 1) (-3 + x) 2 2 - 24 %6 %1 - 4 %5 - 8 %5 %1 + 24 %5 %1 3 %1 := RootOf(-6 _Z + 24 _Z - 1) 2 2 3 %2 := alpha %1 x 2 3 3 %3 := alpha %1 x 2 2 2 %4 := alpha %1 x 2 2 %5 := ln((44405964 %1 + 9425669 x - 1584847 x - 660174 alpha 3 2 + 250239 x - 16410072 %1 - 4754541 + 1692159 alpha 4 3 + 105423552 alpha %1 x - 79888464 alpha %1 x 2 3 2 2 + 64711872 alpha %1 x - 5392656 %4 + 1464216 alpha %1 x 2 4 4 3 - 12893400 alpha %1 x + 486128448 %1 x - 37394496 %1 x 2 3 3 - 14106360 alpha %1 x - 26629488 alpha %1 x 4 3 3 3 2 + 35141184 alpha %1 x - 11110656 %1 x - 84806556 %1 x 2 2 2 + 48547544 %1 x - 5470024 %1 x + 14801988 %1 x 2 3 3 2 2 + 10879340 %1 x + 3532904 %1 x + 564053 alpha x 2 3 - 1692159 alpha x + 1980522 alpha x + 660174 alpha x 2 2 2 3 - 8321088 alpha %1 x - 4297800 alpha %1 x 3 + 3624934 alpha %1 x + 21570624 %3 - 2773696 %2 2 3 2 2 2 - 564053 alpha x - 220058 alpha x + 2041588 alpha %1 x 2 - 621490 alpha %1 x + 10874802 alpha %1 x 2 3 2 2 - 4702120 alpha %1 x - 16177968 alpha %1 2 2 + 4392648 alpha %1 - 1864470 alpha %1 + 6124764 alpha %1 3 - 329225088 %1 x)/((x + 1) (-3 + x))) 2 2 3 %6 := ln((144 %1 + 11 x - 33 + 33 x + 12 alpha - 11 x - 24 %1 2 4 3 + 12 alpha - 1440 alpha %1 x + 768 alpha %1 x 2 3 2 2 + 192 alpha %1 x - 16 %4 - 20 alpha %1 x 2 4 4 3 2 - 54 alpha %1 x + 2304 %1 x + 768 %1 x + 372 alpha %1 x 3 3 4 3 3 3 + 256 alpha %1 x - 480 alpha %1 x + 256 %1 x 2 2 2 2 2 3 - 640 %1 x - 80 %1 x + 200 %1 x - 16 %1 x - 192 %1 x 3 2 3 2 2 2 - 768 %1 x - 96 %1 x + 4 alpha x - 12 alpha x 3 2 2 - 36 alpha x - 12 alpha x + 24 alpha %1 x 2 3 3 - 18 alpha %1 x - 44 alpha %1 x + 64 %3 + 8 %2 2 3 2 2 2 2 - 4 alpha x + 4 alpha x + 2 alpha %1 x + 4 alpha %1 x 2 3 2 2 - 132 alpha %1 x + 124 alpha %1 x - 48 alpha %1 2 2 - 60 alpha %1 + 12 alpha %1 + 6 alpha %1)/( (x + 1) (-3 + x))) Helmut -------------------------------------------------------------- The mail server of my ISP cannot send or receive an email that is attached to a newsgroup article. Please, send such an email separately. Sorry for any inconvenience you might have! -------------------------------------------------------------- ============================================================================== From: rusin@vesuvius.math.niu.edu (Dave Rusin) Newsgroups: sci.math.symbolic Subject: Re: how can i find information on symbolic integration Date: 5 Mar 1999 06:48:12 GMT Helmut Kahovec gave an excellent response to my query. I had written, > > I know that algorithms exist and I understand the basic framework, but I've > never looked in detail at how they operate or how concrete they are. So I am > curious about the extent to which the theoretical procedures can be, and > have been, implemented. For example, can any program give a symbolic > antiderivative in the following example? I ask because the underlying > question here is not really integration but algebraic geometry, and I > don't even know if software to do _that_ problem is available. > > Here goes: in Maple format, let eq:=-x^2-6*x+3+(-9-3*x^2)*y+(x^2+3)*y^3; > > I want to know the integral of this function y with respect to x. I observe that Maple cannot integrate if I try to be too helpful! "Solving" for y as a function of x gives an expression with two cube roots, which is actually ambiguous (of the 3x3=9 combinations of roots, only 3 are actually solutions to eq ); Maple, evidently trying to be cautious, does not assume I necessarily mean the roots to eq when I try to integrate. Kahovec suggests instead integrating with ( infolevel[int]:=5 and ) alias(alpha=RootOf(-x^2-6*x+3+(-9-3*x^2)*y+(x^2+3)*y^3,y)); which works well. In particular, I see in the output the step >int/algrisch/int: genus of the function field 0 which suggests that Maple can actually do more algebraic geometry than I had thought! I would appreciate hearing of its abilities from those who have used it for this. It is interesting to see what else Maple can do this way. I tried the same game with the curve y^2 + x(x-1)(x^2+1) = 0, of genus 2. Since I guess Maple does not know any functions associated with curves of genus greater than 1, it returns an unevaluated integral. Note that the _definite_ integral \integral_{x=0}^{x=1} y(x) dx is also of interest here (part of the period matrix or something -- I don't remember much Riemann surface material). I tried to define alpha as before and compute int(alpha,x=0..1), but it seems to recognize this as a contour integral and then quit. On the other hand, asking for int(sqrt(-x*(x-1)*(x^2+1)),x=0..1) did send it off through some Risch-like contortions, to produce a lovely, complicated answer. I don't really understand the theory well enough to explain why this definite integral can be expected to have a "closed-form" value while the indefinite integral is checked similarly but found wanting, nor do I really see why int(y(x),x=0..1) works but not int(alpha,x=0..1) while in the previous problem int(alpha,x) worked but not int(y(x),x). Well, looks like it's time to really read that differential field theory! dave ============================================================================== From: Helmut Kahovec Newsgroups: sci.math.symbolic Subject: Re: how can i find information on symbolic integration Date: Fri, 05 Mar 1999 12:06:23 GMT Dave Rusin wrote: > > I don't really understand the theory well enough to explain why > this definite integral can be expected to have a "closed-form" value while > the indefinite integral is checked similarly but found wanting, nor > do I really see why int(y(x),x=0..1) works but not int(alpha,x=0..1) while > in the previous problem int(alpha,x) worked but not int(y(x),x). > I also do not understand the theory behind those algorithms. However, if you again set the information level w.r.t. int() to 5 and successfully compute int(y(x),x); and int(y(x),x=1..0); you will see that Maple follows a completely different line of computation than in your previous example. It discards algebraic integration and uses elliptic integration instead. Thus introducing algebraic numbers via RootOf() expressions forces Maple to use Risch/Trager's algorithm for algebraic functions and prevents it from using elliptic integration. As a rule of thumb: if you do not know the type of result beforehand then try both, algebraic and elliptic integration, i.e. in your second example, try RootOf() and sqrt() expressions. Helmut -------------------------------------------------------------- The mail server of my ISP cannot send or receive an email that is attached to a newsgroup article. Please, send such an email separately. Sorry for any inconvenience you might have! --------------------------------------------------------------