       Re: Re: How to get an answer as a Root object?

• To: mathgroup at smc.vnet.net
• Subject: [mg57189] Re: [mg57156] Re: [mg57137] How to get an answer as a Root object?
• From: DrBob <drbob at bigfoot.com>
• Date: Fri, 20 May 2005 04:43:23 -0400 (EDT)
• References: <200505170520.BAA25934@smc.vnet.net> <200505190708.DAA13114@smc.vnet.net>
• Reply-to: drbob at bigfoot.com
• Sender: owner-wri-mathgroup at wolfram.com

```Two questions, Daniel:

1) Is this the missing definition of sroot?

sroot = FindRoot[expr == 0, {s, 2.7}]

2) Doesn't Drop do nothing at all in this statement?

fax2 = Drop[Map[First,FactorList[poly, Extension->Automatic]]];

Bobby

On Thu, 19 May 2005 03:08:33 -0400 (EDT), Daniel Lichtblau <danl at wolfram.com> wrote:

> David W. Cantrell wrote:
>> This question must have an obvious answer!
>>
>> A certain expression, given below my signature, involves nothing worse than
>> some nested square roots. It equals 0 for a value of the variable s near
>> 2.707, easily approximated using FindRoot. However, I would like to get
>> that value expressed precisely as a Root object. Trying to use Solve or
>> Reduce to do that, Mathematica seems to take "forever".
>>
>> Of course, I could do this by hand; in the days before CASs -- before I got
>> spoiled by them! -- I wouldn't have hesitated to do so. But surely there
>> must be an easy way to get Mathematica to give that value as a Root object.
>> How?
>>
>> David Cantrell
>>
>>
>> (-1 + Sqrt - Sqrt[(1 + 2*Sqrt - s)*(-1 + s)] + s)*(-1 - Sqrt +
>> Sqrt*s - Sqrt[-1 - 2*Sqrt + 2*(2 + Sqrt - s)*s]) -
>> (-1 - Sqrt + Sqrt*s + Sqrt[-1 - 2*Sqrt + 2*(2 + Sqrt - s)*s])*
>> (-2 - 2*Sqrt + 2*s + 2*Sqrt[1 - (1/4)*(-1 - Sqrt + Sqrt[(1 +
>> 2*Sqrt - s)*(-1 + s)] + s)^2])
>
> There are a few reasonable methods. Here are some.
>
> expr = (-1 + Sqrt - Sqrt[(1 + 2*Sqrt - s)*(-1 + s)] + s)*(-1 -
> Sqrt +
>    Sqrt*s - Sqrt[-1 - 2*Sqrt + 2*(2 + Sqrt - s)*s]) -
>    (-1 - Sqrt + Sqrt*s + Sqrt[-1 - 2*Sqrt + 2*(2 + Sqrt -
> s)*s])*
>    (-2 - 2*Sqrt + 2*s + 2*Sqrt[1 - (1/4)*(-1 - Sqrt + Sqrt[(1 +
>    2*Sqrt - s)*(-1 + s)] + s)^2]);
>
> Method 1.
>
> First we use numerical approximation to recognize a root polynomial.
> Due, I suspect, to nested radicals, we'll get a small spurious imaginary
> part which I remove using Re.
>
> In:= InputForm[srootval = Re[s/.sroot]]
> Out//InputForm=
> 2.706945837045705738269073202812169790793196895164013983001575249106155514352\
> 266542032260316431664549796229819466076004645196713561849037`120.
>
> To find a candidate polynomial we'll use Recognize (uses a method based
> on lattice reduction), then factor, and find the factor with the right root.
>
> <<NumberTheory`Recognize`
> fax = Drop[Map[First,FactorList[Recognize[srootval, 30, s]]],1];
> InputForm[minpoly = First[Select[fax, Chop[#/.sroot,10^(-acc)]===0&]]]
>
> Why is the factorization needed? Because the larger product polynomial
> has smaller vector norm in the lattice than does the actual minimal
> polynomial.
>
> Method 2.
>
> We'll compute a Groebner basis with 's' the main variable. The hope,
> born out in this case, is that internal handling of algebraics will
> force to the fore a "polynomial" in s.
>
> gb = GroebnerBasis[expr, s];
>
> I do not show it here, but the polynomial we want is readily seen to be
> the first one. Again we'll factor and find the one with the correct
> root. As there are algebraics in the coefficients of this polynomial we
> use Extension->Automatic in the factoring.
>
> poly = First[gb];
> fax2 = Drop[Map[First,FactorList[poly, Extension->Automatic]]];
>
> In:= InputForm[minpoly2 = First[Select[fax2,
>    Chop[#/.sroot,10^(-acc)]===0&]]]
> Out//InputForm=
> -136 - 96*Sqrt + (512 + 364*Sqrt)*s + (-804 - 578*Sqrt)*s^2 +
>   (676 + 486*Sqrt)*s^3 + (-329 - 220*Sqrt)*s^4 + (96 +
> 44*Sqrt)*s^5 -
>   16*s^6
>
> How to represent as a polynomial over the integers? Just multiply by an
> algebraic conjugate.
>
> sqrtConjugate[ee_] := ee /. aa_^(1/2)->-aa^(1/2)
>
> In:= InputForm[minpoly2 = Expand[minpoly2*sqrtConjugate[minpoly2]]]
> Out//InputForm=
> 64 + 512*s - 6112*s^2 + 21024*s^3 - 32136*s^4 + 10832*s^5 + 43568*s^6 -
>   86152*s^7 + 81425*s^8 - 46080*s^9 + 15872*s^10 - 3072*s^11 + 256*s^12
>
>
> Method 3.
>
> Replace radicals by hand, use defining relations explicitly, formulate
> Groebner basis, find polynomial in s alone.
>
> rep1 = {Sqrt->rad0};
> rep2 = {Sqrt[(1+2*Sqrt-s)*(-1+s)] -> rad1} /. rep1;
> rep3 = ({Sqrt[-1-2*Sqrt+2*(2+Sqrt-s)*s] -> rad2,
>    Sqrt[1-(1/4)*(-1-Sqrt+Sqrt[(1+2*Sqrt-s)*(-1+s)]+s)^2] ->
>      rad3} /. rep1 /. rep2)
>
> repssq = Map[#[]^2-#[]^2&, Flatten[{rep1,rep2,rep3}]]
> expr2 = (expr /. rep1) /. reps;
>
> In:= InputForm[Timing[poly4 =
> First[GroebnerBasis[Join[repssq,{expr2}], s,
>
> Out//InputForm=
> {0.30000000000000115*Second, 200704 + 3956736*s - 16568320*s^2 -
>    259776512*s^3 + 2455161856*s^4 - 9217882112*s^5 + 16736476672*s^6 -
>    4378572288*s^7 - 50113953728*s^8 + 119132812800*s^9 - 109419876224*s^10 -
>    38011997376*s^11 + 225422464224*s^12 - 259623250704*s^13 +
>    80225006088*s^14 + 142924054128*s^15 - 196411664327*s^16 +
>    53713807776*s^17 + 132853129136*s^18 - 220218741248*s^19 +
>    191796620928*s^20 - 115209149440*s^21 + 51080841728*s^22 -
>    17021165568*s^23 + 4236362752*s^24 - 766771200*s^25 + 95551488*s^26 -
>    7340032*s^27 + 262144*s^28}
>
> Again one would factor and find the one with the correct root.
>
> I'm sure there is a resultant formulation as well that will get a
> polynomial in s for the root in question.
>
>
> Daniel Lichtblau
> Wolfram Research
>
>
>
>
>

--
DrBob at bigfoot.com

```

• Prev by Date: Re: Modifying displayed form of an expression?
• Next by Date: Re: Modifying displayed form of an expression?
• Previous by thread: Re: How to get an answer as a Root object?
• Next by thread: Re: Re: How to get an answer as a Root object?