Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

Re: bug ?

• To: mathgroup at smc.vnet.net
• Subject: [mg122420] Re: bug ?
• From: Peter Pein <petsie at dordos.net>
• Date: Fri, 28 Oct 2011 05:33:55 -0400 (EDT)
• Delivered-to: l-mathgroup@mail-archive0.wolfram.com
• References: <j862b6$5ls$1@smc.vnet.net>

Am 25.10.2011 12:17, schrieb swiss:
> Following rational function has a double pole matched by a double root at (4+\sqrt{17})/4. The first result is correct, the other ones are false, but come without warning. Can somebody explain this to me?
> 
> In[50]:= Limit[-(((s (4 s - 5)^2 - (5 + 3 Sqrt[17])/
>       8 (4 s + 1)) ((4 s - 5)^2 - 
>      4 (5 + 3 Sqrt[17])/8 (4 s - 7)))/(16 s^2 - 32 s - 1)^2), 
>  s -> 1/4 (4 + Sqrt[17])]
> 
> Out[50]= 1/544 (331 - 19 Sqrt[17])
> 
> In[51]:= 1/544 (331 - 19 Sqrt[17]) // N
> 
> Out[51]= 0.46445
> 
> In[52]:= -(((s (4 s - 5)^2 - (5 + 3 Sqrt[17])/
>        8 (4 s + 1)) ((4 s - 5)^2 - 
>       4 (5 + 3 Sqrt[17])/8 (4 s - 7)))/(16 s^2 - 32 s - 1)^2) /. 
>   s -> 1/4 (4 + Sqrt[17]) // Simplify
> 
> Out[52]= 1/16 (1 + 2 Sqrt[17])
> 
> In[53]:= 1/16 (1 + 2 Sqrt[17]) // N
> 
> Out[53]= 0.577888
> 
> In[54]:= -(((s (4 s - 5)^2 - (5 + 3 Sqrt[17])/
>       8 (4 s + 1)) ((4 s - 5)^2 - 
>      4 (5 + 3 Sqrt[17])/8 (4 s - 7)))/(16 s^2 - 32 s - 1)^2) /. 
>  s -> 1/4 (4 + Sqrt[17.])
> 
> Out[54]= 0.09375
> 

Well, if Cancel[], Factor[] and/or Apart[] can't do that (and I would
consider this a bug), let's do it by hand:

In[1]:= {n, d} = Through[{Numerator, Denominator}[
    Together[-(((s*(4*s - 5)^2 - (1/8)*(5 + 3*Sqrt[17])*(4*s + 1))*
        ((4*s - 5)^2 - (4/8)*(5 + 3*Sqrt[17])*(4*s - 7)))/
       (16*s^2 - 32*s - 1)^2)]]]

    s0 = (1/4)*(4 + Sqrt[17]);

Out[1]= {(-(85 + 21*Sqrt[17] - 100*s - 12*Sqrt[17]*s + 32*s^2))*
   (-5 - 3*Sqrt[17] + 180*s - 12*Sqrt[17]*s - 320*s^2 + 128*s^3),
  16*(-1 - 32*s + 16*s^2)^2}

In[3]:= gcd = PolynomialGCD[n, d, Extension -> Automatic]
Out[3]= (4 + Sqrt[17] - 4*s)^2

In[4]:= pqr = (PolynomialQuotientRemainder[#1, gcd, s] & ) /@ {n, d}
Out[4]= {{408 - 88*Sqrt[17] + (-928 + 144*Sqrt[17])*s +
         (928 - 32*Sqrt[17])*s^2 - 256*s^3, 0},
        {528 - 128*Sqrt[17]    +(-512 + 128*Sqrt[17])*s + 256*s^2, 0}}

In[5]:= y = Divide @@ pqr[[All,1]];
In[6]:= Limit[y, s -> s0]
        N[%]

Out[6]= (1/544)*(331 - 19*Sqrt[17])
Out[7]= 0.4644503549876185

In[8]:= y /. s -> s0
        N[%]
        Simplify[%%]

Out[8]= (408 - 88*Sqrt[17] +(1/16)*(928 - 32*Sqrt[17])*(4 + Sqrt[17])^2
-  4*(4 + Sqrt[17])^3 + (1/4)*(4 + Sqrt[17])*(-928 + 144*Sqrt[17]))/
  (528 - 128*Sqrt[17] + 16*(4 + Sqrt[17])^2 + (1/4)*(4 + Sqrt[17])*
  (-512 + 128*Sqrt[17]))
Out[9]= 0.4644503549876185
Out[10]= (1/544)*(331 - 19*Sqrt[17])

In[11]:= y /. s -> N[s0]
Out[11]= 0.46445035498761833

In[12]:= RootApproximant[%, 2]
Out[12]= (1/544)*(331 - 19*Sqrt[17])


IMHO these steps should belong to Cancel[].

Peter

• Follow-Ups:
  • Re: Re: bug ?
    • From: Andrzej Kozlowski <akoz@mimuw.edu.pl>