On Tue, Jun 10, 2014 at 09:22:26PM +0200, Julien Puydt wrote: > Hi, > > >As far as I know, the cancellation was an issue before, just not in this > >particular example. > > I find the fact that it was working on 32bits with pari 2.5 and it's > now failing on 32bits with pari 2.7 disturbing...
The loop for(y=1,ncoeff,Gvec[d+y]=1.*sum(j=1,d,subst(recG[j],n,y)*Gvec[d+y-j])); is summing 40 values which all suffer from catastrophic cancellation. On example of such value is a = 1/2*s-89/4 b = 1.862645149230957031E-9*s^29-8.102506399154663086E-7*s^28+0.0001682620495557785034*s^27-0.02220491878688335419*s^26+2.091165667399764061*s^25-149.6241444256156683*s^24+8455.149309849366546*s^23-387268.1340125855058*s^22+14641142.25303431042*s^21-462916819.1110430751*s^20+12357645255.57782636*s^19-280452403262.2972860*s^18+5436948501972.644890*s^17-90315237062854.12058*s^16+1287561913303249.816*s^15-15757185306542214.28*s^14+165358894598946191.8*s^13-1.484481047791345401E18*s^12+1.135705563672453560E19*s^11-7.364757121963575366E19*s^10+4.018851784405922897E20*s^9-1.827902353611132295E21*s^8+6.843733896893160794E21*s^7-2.074859461455565033E22*s^6+4.982814333549402435E22*s^5-9.194289489019875941E22*s^4+1.247059559562901079E23*s^3-1.159723714777000500E23*s^2+6.524473424974792687E22*s-1.646906508829389132E22 The product a*b suffer from catastrophic cancellation, thus the result is meaningless The product of b by 89 give slightly different result due to the rounding fix. > Upstream bugs are relevant... Sure, but it is a bug in an unrelated software that is not packaged in Debian, so I cannot reassign it. > and I'm not convinced yet there isn't > one hiding in there. So are you convinced now ? Cheers, -- Bill. <ballo...@debian.org> Imagine a large red swirl here. -- To UNSUBSCRIBE, email to debian-bugs-dist-requ...@lists.debian.org with a subject of "unsubscribe". Trouble? Contact listmas...@lists.debian.org
