axiom-math
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [Axiom-math] Axiom's integration non-deterministic?


From: William Sit
Subject: Re: [Axiom-math] Axiom's integration non-deterministic?
Date: Sun, 26 Mar 2006 14:11:19 -0500

More intrique. But I think I am getting close to the reason (which is NOT what I
said before -- those were wrong guesses).  But, first, some more "strange"
behavior that led me to the chase. Please bear with me.

The integrate(sqrt(1+x^(-2/3)),x) command is run by selecting integrate from
FSINT(INT, EXPR INT), which when applied to the given function, I believe is
equivalent to:

    integrate(f, x) ==
      not real? f => complexIntegrate(f, x)
      -- [... ] irrelevant for this example
      rec := rischNormalize(realElementary(f, x), x)
      g   := rootSimp(rec.func) -- same as g:=f for this example
      --- [... ] more irrelevant code
      i:IR
      if (comp := goComplex?(rtg, tg, ltg)) then -- false for this example
       --- [... ] irrelevant code
       else i := lfintegrate(g, x) -- this is what should be executed
      --- [... ] 
      (u := rinteg(i, f, x, el and ht, comp)) case F => --- true, u := i
        postSubst(u::F, rec.vals, rec.kers, comp, ltg, x)  --- output u
      --- [... ] irrelevant code

When the system first started, as we know:

(1) -> a:=integrate(sqrt(1+x^(-2/3)),x)
                                +---------+
            3+-+2    3+-+       |3+-+2         3+-+2     3+-+     2
        (4x \|x   + 3\|x  + 7x)\|\|x   + 1  + 6\|x   + 9x\|x  + 4x  + 1
   (1)  ---------------------------------------------------------------
                                  +---------+
                       3+-+2      |3+-+2         3+-+
                     (4\|x   + 1)\|\|x   + 1  + 3\|x  + 4x
                                          Type: Union(Expression Integer,...)

comparing this with (new session):

(1) -> complexIntegrate(sqrt(1+x^(-2/3)),x)

                                +---------+
            3+-+2    3+-+       |3+-+2         3+-+2     3+-+     2
        (4x \|x   + 3\|x  + 7x)\|\|x   + 1  + 6\|x   + 9x\|x  + 4x  + 1
   (1)  ---------------------------------------------------------------
                                  +---------+
                       3+-+2      |3+-+2         3+-+
                     (4\|x   + 1)\|\|x   + 1  + 3\|x  + 4x
                                                     Type: Expression Integer
(2) -> complexIntegrate(sqrt(1+x^(-2/3)),x)

                    +---------+
         3+-+2      |3+-+2
   (2)  (\|x   + 1)\|\|x   + 1
                                                     Type: Expression Integer

On the other hand (new session):
(1) -> real? sqrt(1+x^(-2/3))

   (1)  true
                                                                Type: Boolean

so the first branch in 'integrate' should simply fall through! It is not clear
whether 'true' is correct since sqrt(1+x^(-2/3)) can take on complex values. But
that does not seem to be the reason either.

The simplified integral can be done by (new session):

(1) -> q:=sqrt(1+x^(-2/3));
(2) -> p:=lfintegrate(q,x::Symbol)$ElementaryIntegration(INT, EXPR INT)

                      +---------+
                      |3+-+2
           3+-+2      |\|x   + 1
        (x \|x   + x) |---------
                      |  3+-+2
                     \|  \|x
   (2)  -------------------------
                  3+-+2
                  \|x
                                   Type: IntegrationResult Expression Integer
(3) -> rootSimp p

                    +---------+
         3+-+2      |3+-+2
   (3)  (\|x   + 1)\|\|x   + 1
                                                     Type: Expression Integer
(4) -> a:=integrate(sqrt(1+x^(-2/3)),x)
                    +---------+
         3+-+2      |3+-+2
   (4)  (\|x   + 1)\|\|x   + 1
                                          Type: Union(Expression Integer,...)

Now (4) gave the simplifed answer (which by tracing integrate from FSINT(INT,
EXPR INT), as far as I can tell, follows (1), (2).  So does
complexIntegrate(sqrt(1+x^(-2/3)),x). In fact, even tracing complexIntegrate,
which leads to internalIntegrate and internalIntegrate0 will lead back to
lfintegrate for this integrand. It is no longer clear whether loading CHVAR has
anything to do with the simplification. (I used ')set expose add constructor
CHVAR' at the beginning and it does not seem to help.) More digging.

Further investigation with )trace shows that the intermediate steps involve
different integrands (calls to lfintegrate).

                        AXIOM Computer Algebra System
              Version of Tuesday November 30, 2004 at 21:11:14
-----------------------------------------------------------------------------
   Issue )copyright to view copyright notices.
   Issue )summary for a summary of useful system commands.
   Issue )quit to leave AXIOM and return to shell.
-----------------------------------------------------------------------------

(1) -> )set mess autoload off
(1) -> )trace lfintegrate
lfintegrate is not a function

   Nothing is traced now.

(1) -> )trace lfintegrate$INTEF(INT, EXPR INT)

   AXIOM does not understand the use of $elt(INTEF INT (EXPR INT))
      lfintegrate here.
(1) -> )trace INTEF )math

   Parameterized constructors traced:
      INTEF
(1) -> q:=sqrt(1+x^(-2/3))

         +---------+
         |3+-+2
         |\|x   + 1
   (1)   |---------
         |  3+-+2
        \|  \|x
                                                     Type: Expression Integer
(2) -> complexIntegrate(q,x)
1<enter ElementaryIntegration.lfintegrate,58 :
        +---------+
        |3+-+2
       \|\|x   + 1
 arg1= ------------
           3+-+
           \|x
 arg2= x
 2<enter ElementaryIntegration.lfintegrate,58 :
            4      2
        - %F  + 3%F  + 8x %F - 4
  arg1= ------------------------
           4      2
         %F  - 3%F  - 8x %F + 2
  arg2= x
 2>exit  ElementaryIntegration.lfintegrate,58 :
     4        3
  3%F  + 4x %F  + 1
  -----------------
            3
         4%F
1>exit  ElementaryIntegration.lfintegrate,58 :
                         +---------+
     3+-+2    3+-+       |3+-+2         3+-+2     3+-+     2
 (4x \|x   + 3\|x  + 7x)\|\|x   + 1  + 6\|x   + 9x\|x  + 4x  + 1
 ---------------------------------------------------------------
                           +---------+
                3+-+2      |3+-+2         3+-+
              (4\|x   + 1)\|\|x   + 1  + 3\|x  + 4x

                                +---------+
            3+-+2    3+-+       |3+-+2         3+-+2     3+-+     2
        (4x \|x   + 3\|x  + 7x)\|\|x   + 1  + 6\|x   + 9x\|x  + 4x  + 1
   (2)  ---------------------------------------------------------------
                                  +---------+
                       3+-+2      |3+-+2         3+-+
                     (4\|x   + 1)\|\|x   + 1  + 3\|x  + 4x
                                                     Type: Expression Integer
(3) -> complexIntegrate(q,x)
1<enter ElementaryIntegration.lfintegrate,58 :
        +---------+
        |3+-+2
       \|\|x   + 1
 arg1= ------------
           3+-+
           \|x
 arg2= x
 2<enter ElementaryIntegration.lfintegrate,58 :
               5         4       3          2         2
          - 3%O  - 32x %O  + 63%O  + 840x %O  + (4352x  - 486)%O - 4032x
  arg1= -----------------------------------------------------------------
           5         4        3           2           2
        9%O  + 96x %O  - 189%O  - 2520x %O  + (- 3840x  + 972)%O + 12096x
  arg2= x
 2>exit  ElementaryIntegration.lfintegrate,58 :
              5            4          2           3             2
      - 1455%O  - 15520x %O  + (49664x  + 27936)%O  + 407400x %O
    +
              2                       3           2
      (769792x  - 164997)%O - 1390592x  + 3454976x  - 1882188x - 182196
 /
            2
    1787904x  - 94284
1>exit  ElementaryIntegration.lfintegrate,58 :
                  +---------+
     3+-+2        |3+-+2
 (873\|x   + 873)\|\|x   + 1  + 1687
 -----------------------------------
                 873

                    +---------+
         3+-+2      |3+-+2
   (3)  (\|x   + 1)\|\|x   + 1
                                                     Type: Expression Integer
(4) -> complexIntegrate(q,x)
1<enter ElementaryIntegration.lfintegrate,58 :
        +---------+
        |3+-+2
       \|\|x   + 1
 arg1= ------------
           3+-+
           \|x
 arg2= x
 2<enter ElementaryIntegration.lfintegrate,58 :
              5          4         3            2            2
         224%W  + 441x %W  - 2304%W  - 10360x %W  + (- 18081x  + 9728)%W + 19536
x
  arg1= ------------------------------------------------------------------------
-
             5          4         3            2            2
        448%W  + 882x %W  - 4608%W  - 20720x %W  + (- 17640x  + 11264)%W + 39072
x
  arg2= x
 2>exit  ElementaryIntegration.lfintegrate,58 :
              5             4              2             3               2
      495840%W  + 976185x %W  + (- 1366659x  - 4495616)%W  - 22932600x %W
    +
                2                          3            2
    (- 23623677x  + 14280192)%W + 17766567x  - 95304951x  + 35386448x + 42151936

 /
             2
    19133226x  - 8462336
1>exit  ElementaryIntegration.lfintegrate,58 :
                    +---------+
      3+-+2         |3+-+2
 (2066\|x   + 2066)\|\|x   + 1  - 10291
 --------------------------------------
                  2066

                    +---------+
         3+-+2      |3+-+2
   (4)  (\|x   + 1)\|\|x   + 1
                                                     Type: Expression Integer
(5) -> integrate(q,x)
1<enter ElementaryIntegration.lfintegrate,58 :
        +---------+
        |3+-+2
       \|\|x   + 1
 arg1= ------------
           3+-+
           \|x
 arg2= x
 2<enter ElementaryIntegration.lfintegrate,58 :
                5          4        3           2         2
          - 3%BE  - 32x %BE  + 63%BE  + 840x %BE  + (4352x  - 486)%BE - 4032x
  arg1= ----------------------------------------------------------------------
            5          4         3            2           2
        9%BE  + 96x %BE  - 189%BE  - 2520x %BE  + (- 3840x  + 972)%BE + 12096x
  arg2= x
 2>exit  ElementaryIntegration.lfintegrate,58 :
               5             4          2            3              2
      - 1455%BE  - 15520x %BE  + (49664x  + 27936)%BE  + 407400x %BE
    +
              2                        3           2
      (769792x  - 164997)%BE - 1390592x  + 3454976x  - 1882188x - 182196
 /
            2
    1787904x  - 94284
1>exit  ElementaryIntegration.lfintegrate,58 :
                  +---------+
     3+-+2        |3+-+2
 (873\|x   + 873)\|\|x   + 1  + 1687
 -----------------------------------
                 873

                    +---------+
         3+-+2      |3+-+2
   (5)  (\|x   + 1)\|\|x   + 1
                                          Type: Union(Expression Integer,...)

So it seems that because the integrand is algebraic over Q(x), Axiom tries to
find a primitive element of the algebraic extension which contains the
integrand. But a primitive elment is not unique, and there might be some
randomness in selecting this primitive element (which involving building a tower
of algebraic extensions to successively remove the radicals appearing in the
integrand, and so there could be some choice of order that is not deterministic
or some factoring algorithm involved in getting to the minimal polynomials at
each stage, with a different integral also at each stage). Tracing
primitiveElement:

)trace FSPRIMELT )math 

in addition to 

)trace INTEF )math

in repeated calling 

complexIntegrate(q,x) did show this to be observed (output omitted). I am
running out of time right now to trace the primitiveElement algorithm, but I
believe that is the reason (still could be wrong of course). Here are at least
two such primitive elements, each with different result. It seems the ground
field used here is Q (rational) or maybe C (complex field).

1<enter FunctionSpacePrimitiveElement.primitiveElement,93 :
       3+-+
 arg1= \|x
        +---------+
        |3+-+2
 arg2= \|\|x   + 1
1>exit  FunctionSpacePrimitiveElement.primitiveElement,93 :
             +---------+
             |3+-+2        3+-+
 [primelt= 3\|\|x   + 1  + \|x ,

   pol1 = ... ,

   pol2 = ... ,

         6      4        3       2                2
  prim= ?  - 27?  - 56x ?  + 243?  + 432x ? - 512x  - 729]

and 

 1<enter FunctionSpacePrimitiveElement.primitiveElement,93 :
        3+-+
  arg1= \|x
         +---------+
         |3+-+2
  arg2= \|\|x   + 1
 1>exit  FunctionSpacePrimitiveElement.primitiveElement,93 :
                +---------+
                |3+-+2        3+-+
  [primelt= - 4\|\|x   + 1  + \|x ,

   pol1 = ... ,

   pol2 = ... ,
          6      4        3       2                  2
   prim= ?  - 48?  - 98x ?  + 768?  + 1440x ? - 3375x  - 4096]


The other branch of integrate, using lfintegrate, also leads to computing
primitive elements (but less complicated, and apparently CONSISTENT, because the
ground field is Q(x)):

(26) -> lfintegrate(q,x::Symbol)$ElementaryIntegration(INT, EXPR INT)
1<enter ElementaryIntegration.lfintegrate,58 :
        +---------+
        |3+-+2
        |\|x   + 1
 arg1=  |---------
        |  3+-+2
       \|  \|x
 arg2= x
 1<enter FunctionSpacePrimitiveElement.primitiveElement,93 :
        3+-+
  arg1= \|x
         +---------+
         |3+-+2
         |\|x   + 1
  arg2=  |---------
         |  3+-+2
        \|  \|x
 1>exit  FunctionSpacePrimitiveElement.primitiveElement,93 :
             +---------+
             |3+-+2                                                       2
             |\|x   + 1           2                    6     4     2   - x  - 1
  [primelt=  |--------- ,pol1= x ?  - x,pol2= ?,prim= ?  - 3?  + 3?  + --------]

             |  3+-+2                                                      2
            \|  \|x                                                       x
 2<enter ElementaryIntegration.lfintegrate,58 :
  arg1= %CJ
  arg2= x
 2>exit  ElementaryIntegration.lfintegrate,58 :
       3
  x %CJ
1>exit  ElementaryIntegration.lfintegrate,58 :
               +---------+
               |3+-+2
    3+-+2      |\|x   + 1
 (x \|x   + x) |---------
               |  3+-+2
              \|  \|x
 -------------------------
           3+-+2
           \|x

                       +---------+
                       |3+-+2
            3+-+2      |\|x   + 1
         (x \|x   + x) |---------
                       |  3+-+2
                      \|  \|x
   (26)  -------------------------
                   3+-+2
                   \|x
                                   Type: IntegrationResult Expression Integer

Hope this give you sufficient leads to go on.


William




reply via email to

[Prev in Thread] Current Thread [Next in Thread]