John Cremona on Fri, 06 Nov 2015 15:46:11 +0100


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

Re: real polynomials


I am so embarrassed at having such trivially elementary questions, so
let me first say that I do not ask on the list until spending at least
an hour experimenting and reading the library reference manual.

If I have a rational polynomial f created by

GEN f = mkpoln(ncoeffs,g0,g1,g2,g3,g4);

where the g0, ..., g4 are rationals, and then I do

GEN fdash = derivpol(f);
GEN g = ggcd(f,fdash);

I get a run-time error:

  ***   incorrect type in Q_divi_to_int (t_FRAC).

which I do not understand.  Obstacles to debugging include (1) not
knowing how to output a GEN;  (2) asking for the type (or "typ") of a
GEN gives a number (both the above have type 10) when I was expecting
t_POL or similar.

While I am here, in the code above I am clearly defining a polynomial
of degree 4.  I really want to have polynomials of arbitrary degree
but can I do that with the mkpoln() function?  i.e. can I replace have
n separate GEN arguments for the coefficients with a single GEN*
array?

Clearly I should not have given up so easily the first time I tried
programming with the pari library which was some time in the
mid-1990s....

John

On 6 November 2015 at 13:43, John Cremona <john.cremona@gmail.com> wrote:
> Thanks for the helpful comments.  I don't need multiplicities, and in
> fact all I need is to know whether there are any roots at all in
> (-oo,oo) or (-oo,0] or [0,oo).
> At present the coeffs are converted from C doubles.  My pari
> programming skills are rather basic so almost all the program uses
> plain C types and I only convert to GEN for this one test.  Perhaps I
> should be more brave...
>
> John
>
> On 6 November 2015 at 12:36, Karim Belabas
> <Karim.Belabas@math.u-bordeaux.fr> wrote:
>> * John Cremona [2015-11-06 12:13]:
>>> I am using the library functions sturm() and sturmpart() to get the
>>> number of real roots of a real polynomial in the whole real line or
>>> just the positive or negative half-lines.  But occasionally it fails
>>> with
>>>
>>>  ***   domain error in polsturm: issquarefree(pol) = 0
>>>
>>> which I expect is only happening if the discriminant is exactly 0.
>>
>> No, it can occur due to round-off errors whenever the input is inexact
>> and two roots are "close enough".
>>
>> The documentation (??polsturm) specifies that the input must be
>> squarefree. (And, in 2.8.*, that it better be exact...)
>>
>>> Is there a reliable way of taking the squarefreepart, or
>>> the radical, or a real polynomial?
>>
>> Not if the input is inexact.
>>
>>> It is in fact true that in my program the polynomial coefficients are
>>> always rational with denominator a power of 2, so will be known
>>> exactly.
>>
>> In that case, the above problems don't occur. Do you need the number of
>> real roots including multiplicity ? If not, P /= gcd(P,P') is enough.
>>
>>> I would like to avoid having to keep track of the
>>> denominators and using (scaled) integer coefficients.
>>
>> Just use rational coefficients, sturm/sturmpart will cope.
>>
>> If you're using 2.8.* and library mode, you should try
>> ZX_sturmpart / ZX_sturm  (sturmpart itself is now deprecated).
>>
>> Cheers,
>>
>>     K.B.
>> --
>> Karim Belabas, IMB (UMR 5251)  Tel: (+33) (0)5 40 00 26 17
>> Universite de Bordeaux         Fax: (+33) (0)5 40 00 69 50
>> 351, cours de la Liberation    http://www.math.u-bordeaux.fr/~kbelabas/
>> F-33405 Talence (France)       http://pari.math.u-bordeaux.fr/  [PARI/GP]
>> `