Re: How to do quickest
- To: mathgroup at smc.vnet.net
- Subject: [mg123102] Re: How to do quickest
- From: DrMajorBob <btreat1 at austin.rr.com>
- Date: Wed, 23 Nov 2011 07:07:32 -0500 (EST)
- Delivered-to: l-mathgroup@mail-archive0.wolfram.com
- References: <201111210929.EAA14830@smc.vnet.net>
- Reply-to: drmajorbob at yahoo.com
Defining poly as before, here's another small improvement: Clear[x, y, x3, y2] poly2 = Collect[ Eliminate[{poly == 0, x3 == x^3, y2 == y^2}, {y}] /. {Power[x, 3] -> x3, Power[x, 6] -> x3^2} /. Equal -> Subtract, m] 3645 m^12 - 2916 m^10 x + x3^2 + m^6 (270 x3 - 270 y2) - 2 x3 y2 + y2^2 Last@Reap[Do[x3 = x^3; y = Round@Sqrt@x3; y2 = y^2; ! IrreduciblePolynomialQ@poly2 && x3 != y2 && Sow@x, {x, 2, 10^6}]] // Timing {876.934, {{1942, 2878, 3862, 6100, 8380, 11512, 15448, 18694, 31228, 93844, 111382, 117118, 129910, 143950, 186145, 210025, 375376, 445528, 468472, 575800, 844596, 950026}}} The previous timing was 979 seconds. Bobby On Tue, 22 Nov 2011 14:58:11 -0600, Artur <grafix at csl.pl> wrote: > Congratulations Bob! Very good timing 1174s in range 10^6, 2*10^6 on my > computer. > Thank You very much for help! > Best wishes > Artur > > W dniu 2011-11-22 21:17, DrMajorBob pisze: >> The Block construction causes Round[Sqrt[x^3]] to be recomputed each >> time y is mentioned, or so it would seem, and timing bears that out: >> >> Collect[poly = >> Eliminate[{4*m^2 + 6*m*n + n^2 == >> x, (19*m^2 + 9*m*n + n^2)*Sqrt[m^2 + n^2] == y}, {n}] /. >> Equal -> Subtract, m] >> >> 3645 m^12 - 2916 m^10 x + x^6 - 2 x^3 y^2 + y^4 + >> m^6 (270 x^3 - 270 y^2) >> >> Block[{y = Round[Sqrt[x^3]]}, >> Reap[Do[x^3 - y^2 != 0 && >> Or[OddQ[x^6 - 2*x^3*y^2 + y^4], >> Mod[x^6 - 2*x^3*y^2 + y^4, 4] == 0] && ! >> IrreduciblePolynomialQ[poly] && Sow@x, {x, 2, >> 10^4}]][[2]]] // Timing >> >> {10.98, {{1942, 2878, 3862, 6100, 8380}}} >> >> First@Last@ >> Reap[Do[y = Round[Sqrt[x^3]]; >> x^3 - y^2 != 0 && >> Or[OddQ[x^6 - 2*x^3*y^2 + y^4], >> Mod[x^6 - 2*x^3*y^2 + y^4, 4] == 0] && ! >> IrreduciblePolynomialQ[poly] && Sow@x, {x, 2, 10^4}]] // Timing >> >> {8.36221, {1942, 2878, 3862, 6100, 8380}} >> >> This is slightly faster (surprisingly?): >> >> First@Last@ >> Reap[Do[y = Round[Sqrt[x^3]]; ! IrreduciblePolynomialQ[poly] && >> x^3 - y^2 != 0 && >> Or[OddQ[x^6 - 2*x^3*y^2 + y^4], >> Mod[x^6 - 2*x^3*y^2 + y^4, 4] == 0] && Sow@x, {x, 2, >> 10^4}]] // Timing >> >> {8.26819, {1942, 2878, 3862, 6100, 8380}} >> >> And that suggests the Eisenstein criterion is unhelpful after all, so: >> >> First@Last@ >> Reap[Do[y = Round[Sqrt[x^3]]; ! IrreduciblePolynomialQ[poly] && >> x^3 - y^2 != 0 && Sow@x, {x, 2, 10^6}]] // Timing >> >> {978.795, {1942, 2878, 3862, 6100, 8380, 11512, 15448, 18694, 31228, >> 93844, 111382, 117118, 129910, 143950, 186145, 210025, 375376, >> 445528, 468472, 575800, 844596, 950026}} >> >> Bobby >> >> On Tue, 22 Nov 2011 12:40:48 -0600, Andrzej Kozlowski >> <akoz at mimuw.edu.pl> wrote: >> >>> I think the correct use of Eisenstein's criterion is this: >>> >>> >>> Block[{y = Round[Sqrt[x^3]]}, >>> Reap[Table[ >>> If[x^3 - y^2 != 0 && >>> Not[Mod[x^6 - 2*x^3*y^2 + y^4, 2] == 0 && >>> Mod[x^6 - 2*x^3*y^2 + y^4, 4] != 0] && ! >>> IrreduciblePolynomialQ[poly], Sow[{x, y}]], {x, 2, >>> 1000000}]][[2]]] // Timing >>> >>> This is because 3645 is not divisible by 2, but 2916 and 270 both >>> are. So if the term x^6-2 x^3 y^2+y^4 is divisible by 2 but not >>> divisible by 4, the polynomial is irreducible and there is no need to >>> test it further. Only when this isn't the case, we need to use >>> IrreduciblePolynomialQ. >>> >>> I am not sure if adding this will speed up the code. I only have one >>> copy of Mathematica and I need to use it so I can't afford the time to >>> run these programs now, when I am a hurry. >>> >>> Andrzej >>> >>> >>> >>> On 22 Nov 2011, at 19:30, Andrzej Kozlowski wrote: >>> >>>> My memory of Eisenstein's criterion was wrong (also, I was too much >>>> in a hurry to look it up). Rather than correcting it I got rid of it >>>> altogether since I think Mathematica probably uses it anyway. I then >>>> get: >>>> >>>> In[1]:= Collect[ >>>> poly = Eliminate[{4*m^2 + 6*m*n + n^2 == >>>> x, (19*m^2 + 9*m*n + n^2)*Sqrt[m^2 + n^2] == y}, {n}] /. >>>> Equal -> Subtract, m] >>>> >>>> Out[1]= 3645 m^12-2916 m^10 x+m^6 (270 x^3-270 y^2)+x^6-2 x^3 y^2+y^4 >>>> >>>> In[2]:= Block[{y = Round[Sqrt[x^3]]}, >>>> Reap[Table[ >>>> If[x^3 - y^2 != 0 && Not[IrreduciblePolynomialQ[poly]], >>>> Sow[{x, y}]], {x, 2, 1000000}]][[2]]] // Timing >>>> >>>> >>>> Out[2]= {1089.54,({1942,85580} {2878,154396} {3862,240004} >>>> {6100,476425} {8380,767125} {11512,1235168} >>>> {15448,1920032} {18694,2555956} {31228,5518439} >>>> {93844,28748141} {111382,37172564} {117118,40080716} >>>> {129910,46823500} {143950,54615700} {186145,80311375} >>>> {210025,96251275} {375376,229985128} {445528,297380512} >>>> {468472,320645728} {575800,436925600} {844596,776199807} >>>> {950026,925983476} >>>> >>>> )} >>>> >>>> This gets all the numbers but is much slower (I guess it will be >>>> better to add the Eisenstein criterion after all, but of course, in >>>> correct form). >>>> >>>> Andrzej Kozlowski >>>> >>>> >>>> >>>> >>>> On 22 Nov 2011, at 18:02, Artur wrote: >>>> >>>>> Dear Andrzej, >>>>> Your procedure omiited some points >>>>> 6100 >>>>> 8380 Best wishes >>>>> Artur >>>>> >>>>> W dniu 2011-11-22 13:23, Andrzej Kozlowski pisze: >>>>>> On 22 Nov 2011, at 10:07, Andrzej Kozlowski wrote: >>>>>> >>>>>> >>>>>>> On 22 Nov 2011, at 10:06, Andrzej Kozlowski wrote: >>>>>>> >>>>>>> >>>>>>>> On 21 Nov 2011, at 10:29, Artur wrote: >>>>>>>> >>>>>>>> >>>>>>>>> Dear Mathematica Gurus, >>>>>>>>> How to do quickest following procedure (which is very slowly): >>>>>>>>> >>>>>>>>> qq = {}; Do[y = Round[Sqrt[x^3]]; >>>>>>>>> If[(x^3 - y^2) != 0, >>>>>>>>> kk = m /. Solve[{4 m^2 + 6 m n + n^2 == >>>>>>>>> x, (19 m^2 + 9 m n + n^2) Sqrt[m^2 + n^2] == y}, {m, >>>>>>>>> n}][[1]]; >>>>>>>>> ll = CoefficientList[MinimalPolynomial[kk][[1]], #1]; >>>>>>>>> lll = Length[ll]; >>>>>>>>> If[lll < 12, Print[{x/(x^3 - y^2)^2, kk, x, y, x^3 - y^2}]; >>>>>>>>> If[Length[ll] == 3, Print[{kk, x, y}]]]], {x, 2, 1000000}]; >>>>>>>>> qq >>>>>>>>> >>>>>>>>> >>>>>>>>> (*Best wishes Artur*) >>>>>>>>> >>>>>>>>> >>>>>>>> I think it would be better to send not only the code but also the >>>>>>>> mathematical problem, as there may be a way to do it in a >>>>>>>> different way. Unless I am misunderstanding something, what you >>>>>>>> are trying to do is the same as this: >>>>>>>> >>>>>>>> In[31]:= Block[{y = Round[Sqrt[x^3]]}, >>>>>>>> Reap[Table[ >>>>>>>> If[(x^3 - y^2) != 0 && Not[IrreduciblePolynomialQ[poly]], >>>>>>>> Sow[{x, y}]], {x, 2, 1000000}]][[2]]] // Timing >>>>>>>> >>>>>>>> Out[31]= {721.327,{}} >>>>>>>> >>>>>>>> This ought to be a lot faster than your code, but I have not >>>>>>>> tried to run yours to the end. Also, it is possible that using >>>>>>>> the Eisenstein Test explicitly may be somewhat faster: >>>>>>>> >>>>>>>> Block[{y = Round[Sqrt[x^3]]}, >>>>>>>> Reap[Table[ >>>>>>>> If[x^3 - y^2 != 0 && Mod[x^6 - 2*x^3*y^2 + y^4, 4] == 0 && >>>>>>>> ! IrreduciblePolynomialQ[poly], Sow[{x, y}]], {x, 2, >>>>>>>> 1000000}]][[2]]] >>>>>>>> >>>>>>>> {} >>>>>>>> >>>>>>>> but I forgot to use Timing and don't want to wait again, >>>>>>>> particularly that the answer is the empty set. >>>>>>>> >>>>>>>> Andrzej Kozlowski >>>>>>>> >>>>>>> I forgot to include the definition of poly: >>>>>>> >>>>>>> Collect[poly = Eliminate[{4*m^2 + 6*m*n + n^2 == x, >>>>>>> (19*m^2 + 9*m*n + n^2)*Sqrt[m^2 + n^2] == y}, {n}] /. Equal -> >>>>>>> Subtract, m] >>>>>>> >>>>>>> 3645*m^12 - 2916*m^10*x + m^6*(270*x^3 - 270*y^2) + x^6 - >>>>>>> 2*x^3*y^2 + y^4 >>>>>>> >>>>>>> Andrzej Kozlowski >>>>>>> >>>>>> >>>>>> Strange but I run this code with a fresh kernel and got the >>>>>> following answers: >>>>>> >>>>>> In[1]:= >>>>>> Collect[poly=Eliminate[{4*m^2+6*m*n+n^2==x,(19*m^2+9*m*n+n^2)*Sqrt[m^2+n^2]==y},{n}]/.Equal->Subtract,m] >>>>>> Out[1]= 3645 m^12-2916 m^10 x+m^6 (270 x^3-270 y^2)+x^6-2 x^3 >>>>>> y^2+y^4 >>>>>> >>>>>> In[2]:= >>>>>> Block[{y=Round[Sqrt[x^3]]},Reap[Table[If[x^3-y^2!=0&&Mod[x^6-2*x^3*y^2+y^4,4]==0&&!IrreduciblePolynomialQ[poly],Sow[{x,y}]],{x,2,1000000}]][[2]]]//Timing >>>>>> >>>>>> Out[2]= {766.05,({1942,85580} {2878,154396} {3862,240004} >>>>>> {11512,1235168} {15448,1920032} {18694,2555956} >>>>>> {111382,37172564} {117118,40080716} {129910,46823500} >>>>>> {143950,54615700} {186145,80311375} {210025,96251275} >>>>>> {375376,229985128} {445528,297380512} {468472,320645728} >>>>>> {575800,436925600} {950026,925983476} >>>>>> >>>>>> )} >>>>>> >>>>>> >>>>>> I tested the first one and it does seem to be a solution to your >>>>>> problem. >>>>>> >>>>>> {x, y} = {950026, 925983476}; >>>>>> >>>>>> y == Round[Sqrt[x^3]] >>>>>> >>>>>> True >>>>>> >>>>>> x^3 - y^2 != 0 >>>>>> >>>>>> True >>>>>> >>>>>> kk = >>>>>> m /. Solve[{4 m^2 + 6 m n + n^2 == >>>>>> x, (19 m^2 + 9 m n + n^2) Sqrt[m^2 + n^2] == y}, {m, n}][[1]] >>>>>> >>>>>> Out[12]= -Sqrt[-(198/5)-(44 I Sqrt[11])/5] >>>>>> >>>>>> ll = CoefficientList[MinimalPolynomial[kk][[1]], #1]; >>>>>> >>>>>> Length[ll] >>>>>> >>>>>> 5 >>>>>> >>>>>> I don't know why I got no answers the first time round, perhaps one >>>>>> of the variables had values assigned. >>>>>> >>>>>> Andrzej >>>>>> >>>>>> >>>>>> >>>>>> >>>> >>> >> >> -- DrMajorBob at yahoo.com
- References:
- How to do quickest
- From: Artur <grafix@csl.pl>
- How to do quickest