BEGIN {
$Math::Prime::Util::ZetaBigFloat::AUTHORITY
=
'cpan:DANAJ'
;
$Math::Prime::Util::ZetaBigFloat::VERSION
=
'0.14'
;
}
my
@_Riemann_Zeta_Table
= (
'0.64493406684822643647241516664602518921894990'
,
'0.20205690315959428539973816151144999076498629'
,
'0.082323233711138191516003696541167902774750952'
,
'0.036927755143369926331365486457034168057080920'
,
'0.017343061984449139714517929790920527901817490'
,
'0.0083492773819228268397975498497967595998635606'
,
'0.0040773561979443393786852385086524652589607906'
,
'0.0020083928260822144178527692324120604856058514'
,
'0.00099457512781808533714595890031901700601953156'
,
'0.00049418860411946455870228252646993646860643576'
,
'0.00024608655330804829863799804773967096041608846'
,
'0.00012271334757848914675183652635739571427510590'
,
'0.000061248135058704829258545105135333747481696169'
,
'0.000030588236307020493551728510645062587627948707'
,
'0.000015282259408651871732571487636722023237388990'
,
'0.0000076371976378997622736002935630292130882490903'
,
'0.0000038172932649998398564616446219397304546972190'
,
'0.0000019082127165539389256569577951013532585711448'
,
'0.00000095396203387279611315203868344934594379418741'
,
'0.00000047693298678780646311671960437304596644669478'
,
'0.00000023845050272773299000364818675299493504182178'
,
'0.00000011921992596531107306778871888232638725499778'
,
'0.000000059608189051259479612440207935801227503918837'
,
'0.000000029803503514652280186063705069366011844730920'
,
'0.000000014901554828365041234658506630698628864788168'
,
'0.0000000074507117898354294919810041706041194547190319'
,
'0.0000000037253340247884570548192040184024232328930593'
,
'0.0000000018626597235130490064039099454169480616653305'
,
'0.00000000093132743241966818287176473502121981356795514'
,
'0.00000000046566290650337840729892332512200710626918534'
,
'0.00000000023283118336765054920014559759404950248298228'
,
'0.00000000011641550172700519775929738354563095165224717'
,
'0.000000000058207720879027008892436859891063054173122605'
,
'0.000000000029103850444970996869294252278840464106981987'
,
'0.000000000014551921891041984235929632245318420983808894'
,
'0.0000000000072759598350574810145208690123380592648509256'
,
'0.0000000000036379795473786511902372363558732735126460284'
,
'0.0000000000018189896503070659475848321007300850305893096'
,
'0.00000000000090949478402638892825331183869490875386000099'
,
'0.00000000000045474737830421540267991120294885703390452991'
,
'0.00000000000022737368458246525152268215779786912138298220'
,
'0.00000000000011368684076802278493491048380259064374359028'
,
'0.000000000000056843419876275856092771829675240685530571589'
,
'0.000000000000028421709768893018554550737049426620743688265'
,
'0.000000000000014210854828031606769834307141739537678698606'
,
'0.0000000000000071054273952108527128773544799568000227420436'
,
'0.0000000000000035527136913371136732984695340593429921456555'
,
'0.0000000000000017763568435791203274733490144002795701555086'
,
'0.00000000000000088817842109308159030960913863913863256088715'
,
'0.00000000000000044408921031438133641977709402681213364596031'
,
'0.00000000000000022204460507980419839993200942046539642366543'
,
'0.00000000000000011102230251410661337205445699213827024832229'
,
'0.000000000000000055511151248454812437237365905094302816723551'
,
'0.000000000000000027755575621361241725816324538540697689848904'
,
'0.000000000000000013877787809725232762839094906500221907718625'
,
'0.0000000000000000069388939045441536974460853262498092748358742'
,
'0.0000000000000000034694469521659226247442714961093346219504706'
,
'0.0000000000000000017347234760475765720489729699375959074780545'
,
'0.00000000000000000086736173801199337283420550673429514879071415'
,
'0.00000000000000000043368086900206504874970235659062413612547801'
,
'0.00000000000000000021684043449972197850139101683209845761574010'
,
'0.00000000000000000010842021724942414063012711165461382589364744'
,
'0.000000000000000000054210108624566454109187004043886337150634224'
,
'0.000000000000000000027105054312234688319546213119497764318887282'
,
'0.000000000000000000013552527156101164581485233996826928328981877'
,
'0.0000000000000000000067762635780451890979952987415566862059812586'
,
'0.0000000000000000000033881317890207968180857031004508368340311585'
,
'0.0000000000000000000016940658945097991654064927471248619403036418'
,
'0.00000000000000000000084703294725469983482469926091821675222838642'
,
'0.00000000000000000000042351647362728333478622704833579344088109717'
,
'0.00000000000000000000021175823681361947318442094398180025869417612'
,
'0.00000000000000000000010587911840680233852265001539238398470699902'
,
'0.000000000000000000000052939559203398703238139123029185055866375629'
,
'0.000000000000000000000026469779601698529611341166842038715592556134'
,
'0.000000000000000000000013234889800848990803094510250944989684323826'
,
'0.0000000000000000000000066174449004244040673552453323082200147137975'
,
'0.0000000000000000000000033087224502121715889469563843144048092764894'
,
'0.0000000000000000000000016543612251060756462299236771810488297723589'
,
'0.00000000000000000000000082718061255303444036711056167440724040096811'
,
'0.00000000000000000000000041359030627651609260093824555081412852575873'
,
'0.00000000000000000000000020679515313825767043959679193468950443365312'
,
'0.00000000000000000000000010339757656912870993284095591745860911079606'
,
'0.000000000000000000000000051698788284564313204101332166355512893608164'
,
'0.000000000000000000000000025849394142282142681277617708450222269121159'
,
'0.000000000000000000000000012924697071141066700381126118331865309299779'
,
'0.0000000000000000000000000064623485355705318034380021611221670660356864'
,
'0.0000000000000000000000000032311742677852653861348141180266574173608296'
,
'0.0000000000000000000000000016155871338926325212060114057052272720509148'
,
'0.00000000000000000000000000080779356694631620331587381863408997398684847'
,
'0.00000000000000000000000000040389678347315808256222628129858130379479700'
,
'0.00000000000000000000000000020194839173657903491587626465673047518903728'
,
'0.00000000000000000000000000010097419586828951533619250700091044144538432'
,
'0.000000000000000000000000000050487097934144756960847711725486604360898735'
,
'0.000000000000000000000000000025243548967072378244674341937966175648398693'
,
'0.000000000000000000000000000012621774483536189043753999660777148710632765'
,
'0.0000000000000000000000000000063108872417680944956826093943332037500694712'
,
'0.0000000000000000000000000000031554436208840472391098412184847972814371270'
,
'0.0000000000000000000000000000015777218104420236166444327830159601782237092'
,
'0.00000000000000000000000000000078886090522101180735205378276604136878962534'
,
'0.00000000000000000000000000000039443045261050590335263935513575963608141044'
,
'0.00000000000000000000000000000019721522630525295156852383215213909988473843'
,
'0.000000000000000000000000000000098607613152626475748329967604159218377505181'
,
'0.000000000000000000000000000000049303806576313237862187667644776975622245754'
,
'0.000000000000000000000000000000024651903288156618927101395103287812527732549'
,
'0.000000000000000000000000000000012325951644078309462219884645277065145764150'
,
'0.0000000000000000000000000000000061629758220391547306663380205162648609383631'
,
'0.0000000000000000000000000000000030814879110195773651853009095507130250105264'
,
'0.0000000000000000000000000000000015407439555097886825433610878728841686496904'
,
'0.00000000000000000000000000000000077037197775489434125525075496895150086398231'
,
'0.00000000000000000000000000000000038518598887744717062214878116197893873445220'
,
'0.00000000000000000000000000000000019259299443872358530924885847349054449873362'
,
'0.000000000000000000000000000000000096296497219361792654015918534245633717541108'
,
'0.000000000000000000000000000000000048148248609680896326805122366289604787579935'
,
'0.000000000000000000000000000000000024074124304840448163334948882867065229914248'
,
'0.000000000000000000000000000000000012037062152420224081644937008007620275295506'
,
'0.0000000000000000000000000000000000060185310762101120408149560261951727031681191'
,
'0.0000000000000000000000000000000000030092655381050560204049738538280405431094080'
,
'0.0000000000000000000000000000000000015046327690525280102016522071575050028177934'
,
'0.00000000000000000000000000000000000075231638452626400510054786365991407868525313'
,
'0.00000000000000000000000000000000000037615819226313200255018118519034423181524371'
,
'0.00000000000000000000000000000000000018807909613156600127505967704863451341028548'
,
'0.000000000000000000000000000000000000094039548065783000637519533342138055875645097'
,
'0.000000000000000000000000000000000000047019774032891500318756331610342627662060287'
,
'0.000000000000000000000000000000000000023509887016445750159377020784929180405960294'
,
'0.000000000000000000000000000000000000011754943508222875079688128719050545728002924'
,
'0.0000000000000000000000000000000000000058774717541114375398439371350539247056872356'
,
'0.0000000000000000000000000000000000000029387358770557187699219261593698463000750878'
,
'0.0000000000000000000000000000000000000014693679385278593849609489436325511324487536'
,
'0.00000000000000000000000000000000000000073468396926392969248046975979881822702829326'
,
'0.00000000000000000000000000000000000000036734198463196484624023330922692333378216377'
,
'0.00000000000000000000000000000000000000018367099231598242312011613105596640698043218'
,
'0.000000000000000000000000000000000000000091835496157991211560057891008818116853335663'
,
'0.000000000000000000000000000000000000000045917748078995605780028887331354029547708393'
,
'0.000000000000000000000000000000000000000022958874039497802890014424274658671814201226'
,
'0.000000000000000000000000000000000000000011479437019748901445007205673656554920549667'
,
'0.0000000000000000000000000000000000000000057397185098744507225036006822706837980911955'
,
'0.0000000000000000000000000000000000000000028698592549372253612517996229494773449843879'
,
'0.0000000000000000000000000000000000000000014349296274686126806258995720794504878051247'
,
'0.00000000000000000000000000000000000000000071746481373430634031294970624129584900687276'
,
'0.00000000000000000000000000000000000000000035873240686715317015647482652117145953820656'
,
'0.00000000000000000000000000000000000000000017936620343357658507823740439409357478069335'
,
'0.000000000000000000000000000000000000000000089683101716788292539118699241549402394210037'
,
'0.000000000000000000000000000000000000000000044841550858394146269559348635608906198392806'
,
'0.000000000000000000000000000000000000000000022420775429197073134779673989415854766292332'
,
'0.000000000000000000000000000000000000000000011210387714598536567389836885245061272178142'
,
'0.0000000000000000000000000000000000000000000056051938572992682836949184061349085990997301'
,
'0.0000000000000000000000000000000000000000000028025969286496341418474591909049136205534180'
,
'0.0000000000000000000000000000000000000000000014012984643248170709237295913982765839445600'
,
'0.00000000000000000000000000000000000000000000070064923216240853546186479434774488319489698'
,
'0.00000000000000000000000000000000000000000000035032461608120426773093239672340797200498749'
,
'0.00000000000000000000000000000000000000000000017516230804060213386546619821154916280500674'
,
'0.000000000000000000000000000000000000000000000087581154020301066932733099055722973670007705'
,
'0.000000000000000000000000000000000000000000000043790577010150533466366549511177617590838630'
,
'0.000000000000000000000000000000000000000000000021895288505075266733183274750027519047364241'
,
'0.000000000000000000000000000000000000000000000010947644252537633366591637373159996274330429'
,
'0.0000000000000000000000000000000000000000000000054738221262688166832958186859620770540479841'
,
'0.0000000000000000000000000000000000000000000000027369110631344083416479093427750648326515819'
,
'0.0000000000000000000000000000000000000000000000013684555315672041708239546713188745182016542'
,
'0.00000000000000000000000000000000000000000000000068422776578360208541197733563655129305944821'
,
'0.00000000000000000000000000000000000000000000000034211388289180104270598866781064699118259780'
,
'0.00000000000000000000000000000000000000000000000017105694144590052135299433390278061047559013'
,
'0.000000000000000000000000000000000000000000000000085528470722950260676497166950542676865892145'
,
'0.000000000000000000000000000000000000000000000000042764235361475130338248583474988795642311765'
,
'0.000000000000000000000000000000000000000000000000021382117680737565169124291737400216890944447'
,
'0.000000000000000000000000000000000000000000000000010691058840368782584562145868668714802068411'
,
'0.0000000000000000000000000000000000000000000000000053455294201843912922810729343238928532329351'
,
'0.0000000000000000000000000000000000000000000000000026727647100921956461405364671584582440160440'
,
'0.0000000000000000000000000000000000000000000000000013363823550460978230702682335780663944745475'
,
'0.00000000000000000000000000000000000000000000000000066819117752304891153513411678864562139278223'
,
'0.00000000000000000000000000000000000000000000000000033409558876152445576756705839419361874822728'
,
'0.00000000000000000000000000000000000000000000000000016704779438076222788378352919705374539139236'
,
'0.000000000000000000000000000000000000000000000000000083523897190381113941891764598512518034789088'
,
'0.000000000000000000000000000000000000000000000000000041761948595190556970945882299251474130425513'
,
'0.000000000000000000000000000000000000000000000000000020880974297595278485472941149624142102889746'
,
'0.000000000000000000000000000000000000000000000000000010440487148797639242736470574811539397337203'
,
'0.0000000000000000000000000000000000000000000000000000052202435743988196213682352874055924806327115'
,
'0.0000000000000000000000000000000000000000000000000000026101217871994098106841176437027371676377257'
,
'0.0000000000000000000000000000000000000000000000000000013050608935997049053420588218513488929259862'
,
'0.00000000000000000000000000000000000000000000000000000065253044679985245267102941092566788283203421'
,
);
@_Riemann_Zeta_Table
=
map
{ Math::BigFloat->new(
$_
) }
@_Riemann_Zeta_Table
;
my
@_Riemann_Zeta_Premult
;
my
$_Riemann_Zeta_Premult_accuracy
;
my
$_Borwein_n
= 55;
my
@_Borwein_dk
= (
'1'
,
'6051'
,
'6104451'
,
'2462539971'
,
'531648934851'
,
'71301509476803'
,
'6504925195108803'
,
'429144511928164803'
,
'21392068013887742403'
,
'832780518854440804803'
,
'25977281563850106233283'
,
'662753606729324750201283'
,
'14062742362385399866745283'
,
'251634235316509414702211523'
,
'3841603462178827861104812483'
,
'50535961819850087101900022211'
,
'577730330374203014014104003011'
,
'5782012706584553297863989289411'
,
'50984922488525881477588707205571'
,
'398333597655022403279683908035011'
,
'2770992240330783259897072664469955'
,
'17238422988353715312442126057365955'
,
'96274027751337344115352100618133955'
,
'484350301573059857715727453968687555'
,
'2201794236784087151947175826243477955'
,
'9068765987529892610841571032285864387'
,
'33926582279822401059328069515697217987'
,
'115535262182820447663793177744255246787'
,
'358877507711760077538925500462137369027'
,
'1018683886695854101193095537014797787587'
,
'2646951832121008166346437186541363159491'
,
'6306464665572570713623910486640730071491'
,
'13799752848354341643763498672558481367491'
,
'27780237373991939435100856211039992177091'
,
'51543378762608611361377523633779417047491'
,
'88324588911945720951614452340280439890371'
,
'140129110249040241501243929391690331218371'
,
'206452706984942815385219764876242498642371'
,
'283527707823296964404071683165658912154051'
,
'364683602811933600833512164561308162744771'
,
'441935796522635816776473230396154031661507'
,
'508231717051242054487234759342047053767107'
,
'559351463001010719709990637083458540691907'
,
'594624787018881191308291683229515933311427'
,
'616297424973434835299724300924272199623107'
,
'628083443816135918099559567176252011864515'
,
'633714604276098212796088600263676671320515'
,
'636056734158553360761837806887547188568515'
,
'636894970116484676875895417679248215794115'
,
'637149280289288581322870186196318041432515'
,
'637213397278310656625865036925470191411651'
,
'637226467136294189739463288384528579584451'
,
'637228536449134002301138291602841035366851'
,
'637228775173095037281299181461988671775171'
,
'637228793021615488494769154535569803469251'
,
'637228793670652595811622608101881844621763'
,
);
sub
_Recompute_Dk {
my
$nterms
=
shift
;
$_Borwein_n
=
$nterms
;
@_Borwein_dk
= ();
foreach
my
$k
(0 ..
$nterms
) {
my
$dsum
= Math::BigFloat->bzero;
$dsum
->accuracy(2
*$_Borwein_n
);
my
$n
= Math::BigInt->new(
$nterms
-1)->bfac;
my
$d
= Math::BigInt->new(
$nterms
)->bfac;
foreach
my
$i
(0 ..
$k
) {
my
$term
= Math::BigFloat->bone;
$term
->accuracy(2
*$_Borwein_n
);
$term
->bmul(
$n
)->bdiv(
$d
);
$dsum
+=
$term
;
$n
->bmul(
$nterms
+
$i
)->bmul(4);
$d
->bdiv(
$nterms
-
$i
)->bmul(2
*$i
+1)->bmul(2
*$i
+2);
}
my
$dk
= (
$nterms
*
$dsum
+ 1e-20)->as_int;
$_Borwein_dk
[
$k
] =
$dk
;
}
}
sub
RiemannZeta {
my
(
$x
) =
@_
;
$x
= new Math::BigFloat
"$x"
if
ref
(
$x
) ne
'Math::BigFloat'
;
my
$xdigits
=
$x
->accuracy || Math::BigFloat->accuracy() || Math::BigFloat->div_scale();
return
$_Riemann_Zeta_Table
[
int
(
$x
)-2]
if
$x
==
int
(
$x
)
&&
defined
$_Riemann_Zeta_Table
[
int
(
$x
)-2]
&&
$xdigits
<= 44;
my
$tol
= 0.0 +
"1e-$xdigits"
;
my
$superx
= 1;
my
$subx
=
$x
->copy;
while
(
$subx
> 8) {
$superx
*= 2;
$subx
/= 2;
}
if
(
$x
> 30) {
my
$sum
= Math::BigFloat->bzero;
$sum
->accuracy(
$xdigits
);
for
my
$k
(4 .. 1000) {
my
$term
= (
$k
** -
$subx
) **
$superx
;
$sum
+=
$term
;
last
if
$term
< (
$sum
*$tol
);
}
$sum
+= ( 3 ** -
$subx
) **
$superx
;
$sum
+= ( 2 ** -
$subx
) **
$superx
;
return
$sum
;
}
{
my
$dig
=
int
(
$_Borwein_n
/ 1.3)+1;
_Recompute_Dk(
int
(
$xdigits
* 1.3) + 4 )
if
$dig
<
$xdigits
;
}
if
(
ref
$_Borwein_dk
[0] ne
'Math::BigInt'
) {
@_Borwein_dk
=
map
{ Math::BigInt->new(
"$_"
) }
@_Borwein_dk
;
}
my
$n
=
$_Borwein_n
;
my
$intermediate_accuracy
=
undef
;
my
$one
= Math::BigFloat->bone;
$one
->accuracy(
$intermediate_accuracy
)
if
defined
$intermediate_accuracy
;
my
$d1
= Math::BigFloat->new(2);
$d1
->accuracy(
$intermediate_accuracy
)
if
defined
$intermediate_accuracy
;
$d1
=
$d1
** (
$one
-
$x
);
my
$divisor
=
$one
->copy->bsub(
$d1
)->bmul(-
$_Borwein_dk
[
$n
]);
$tol
=
$divisor
->copy->bmul(
$tol
)->babs();
my
$sum
= Math::BigFloat->bzero;
$sum
->accuracy(
$intermediate_accuracy
)
if
defined
$intermediate_accuracy
;
foreach
my
$k
(1 ..
$n
-1) {
my
$term
= Math::BigFloat->new(
$_Borwein_dk
[
$k
] -
$_Borwein_dk
[
$n
] );
$term
*= -1
if
$k
% 2;
$term
->accuracy(
$intermediate_accuracy
)
if
defined
$intermediate_accuracy
;
my
$den
= Math::BigFloat->new(
$k
+1);
$den
->accuracy(
$intermediate_accuracy
)
if
defined
$intermediate_accuracy
;
$den
= (
$den
**
$subx
) **
$superx
;
$term
/=
$den
;
$sum
+=
$term
;
last
if
$term
->copy->babs() <
$tol
;
}
$sum
+= Math::BigFloat->new(
$one
-
$_Borwein_dk
[
$n
] );
$sum
->bdiv(
$divisor
);
$sum
->bsub(1);
return
$sum
;
}
sub
RiemannR {
my
(
$x
) =
@_
;
$x
= new Math::BigFloat
"$x"
if
ref
(
$x
) ne
'Math::BigFloat'
;
my
$xdigits
=
$x
->accuracy || Math::BigFloat->accuracy() || Math::BigFloat->div_scale();
my
$tol
= 0.0 +
"1e-$xdigits"
;
if
( (
scalar
@_Riemann_Zeta_Premult
== 0) || (
$_Riemann_Zeta_Premult_accuracy
<
$xdigits
) ) {
$_Riemann_Zeta_Premult_accuracy
=
$xdigits
;
@_Riemann_Zeta_Premult
=
map
{
my
$v
= Math::BigFloat->bone;
$v
->accuracy(
$xdigits
);
$v
/ (
$_Riemann_Zeta_Table
[
$_
-1] *
$_
+
$_
) }
(1 ..
@_Riemann_Zeta_Table
);
}
my
$sum
= Math::BigFloat->bone;
my
$flogx
=
log
(
$x
);
my
$part_term
= Math::BigFloat->bone;
for
my
$k
(1 .. 10000) {
my
$zeta_term
=
$_Riemann_Zeta_Premult
[
$k
-1];
if
(!
defined
$zeta_term
) {
my
$zeta
=
$_Riemann_Zeta_Table
[
$k
-1];
if
(!
defined
$zeta
) {
my
$kz
= Math::BigFloat->new(
$k
+1);
$kz
->accuracy(
$xdigits
);
if
(
$kz
>= 100 &&
$xdigits
<= 40) {
$zeta
= Math::BigFloat->new(3)->bpow(-
$kz
)
+ Math::BigFloat->new(2)->bpow(-
$kz
);
}
else
{
$zeta
= Math::Prime::Util::ZetaBigFloat::RiemannZeta(
$kz
);
}
}
$zeta_term
= Math::BigFloat->bone / (
$zeta
*
$k
+
$k
);
}
$part_term
*=
$flogx
/
$k
;
my
$term
=
$part_term
*
$zeta_term
;
$sum
+=
$term
;
last
if
$term
< (
$tol
*$sum
);
}
return
$sum
;
}
1;