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] >> `