Efficient computation of the Euler-Kronecker constants of prime cyclotomic fields

TL;DR

This paper presents a new, more efficient algorithm for computing Euler-Kronecker constants of prime cyclotomic fields, enabling the evaluation of these constants for very large primes and expanding the known data range.

Contribution

The authors introduce a faster algorithm for calculating Euler-Kronecker constants, allowing computations for larger primes and providing new numerical data.

Findings

Computed Euler-Kronecker constants for large primes, including two new negative values.

Extended the known range of Euler-Kronecker constants up to primes ≤ 10^6.

Discovered a simplified method to compute the difference between Euler-Kronecker constants of a field and its maximal real subfield.

Abstract

We introduce a new algorithm, which is faster and requires less computing resources than the ones previously known, to compute the Euler-Kronecker constants $\mathfrak{G}_q$ for the prime cyclotomic fields $\mathbb{Q}(\zeta_q)$, where $q$ is an odd prime and $\zeta_q$ is a primitive $q$-root of unity. With such a new algorithm we evaluated $\mathfrak{G}_q$ and $\mathfrak{G}_q^+$, where $\mathfrak{G}_q^+$ is the Euler-Kronecker constant of the maximal real subfield of $\mathbb{Q}(\zeta_q)$, for some very large primes $q$ thus obtaining two new negative values of $\mathfrak{G}_q$: $\mathfrak{G}_{9109334831}= -0.248739\dotsc$ and $\mathfrak{G}_{9854964401}= -0.096465\dotsc$ We also evaluated $\mathfrak{G}_q$ and $\mathfrak{G}^+_q$ for every odd prime $q\le 10^6$, thus enlarging the size of the previously known range for $\mathfrak{G}_q$ and $\mathfrak{G}^+_q$. Our method also reveals that difference $\mathfrak{G}_q - \mathfrak{G}^+_q$ can be computed in a much simpler way than both its summands, see Section 3.4. Moreover, as a by-product, we also computed $M_q=\max_{\chi\ne \chi_0} \vert L^\prime/L(1,\chi) \vert $ for every odd prime $q\le 10^6$, where $L(s,\chi)$ are the Dirichlet $L$-functions, $\chi$ run over the non trivial Dirichlet characters mod $q$ and $\chi_0$ is the trivial Dirichlet character mod $q$. As another by-product of our computations, we will also provide more data on the generalised Euler constants in arithmetic progressions. The programs used to performed the computations here described and the numerical results obtained are available at the following web address: \url{http://www.math.unipd.it/~languasc/EK-comput.html}.

Tables

Table 1. Table 1. Comparison of the main characteristics of procedures a ) and b ) to compute 𝔊 q subscript 𝔊 𝑞 \mathfrak{G}_{q} and 𝔊 q + superscript subscript 𝔊 𝑞 \mathfrak{G}_{q}^{+} .

Comparison	Approach a)	Approach a)	Approach b)
	($\psi$ comp. with GSL)	($\psi$ comp. with PARI/GP)
Magnitude of the functions for $x\to 0^{+}$:	$\psi \left(x\right)\sim -1/x$	$\psi \left(x\right)\sim -1/x$	$\log \left(\Gamma \left(x\right)\right)\sim \log \left(1/x\right)$
	$T\left(x\right)\sim \frac{\log \left(1/x\right)}{x}$	$T\left(x\right)\sim \frac{\log \left(1/x\right)}{x}$	$S\left(x\right)\sim {\left(\log x\right)}^2$
Precomputations ($T$ and $S$ with PARI/GP):
needed space for storing precomputed values	$q-1$ values of	$2\left(q-1\right)$ values of	$\left(q-1\right)/2$ values of
($⟨g⟩={\mathbb{Z}}_q^{\ast }$, $a_k:=g^{k}modq$):	$T\left(a_k/q\right)$	$T\left(a_k/q\right)$ and $\psi \left(a_k/q\right)$	$S\left(a_k/q\right)+S\left(1-a_k/q\right)$
number of write operations on hard disks:	$q-1$	$2\left(q-1\right)$	$\left(q-1\right)/2$
number of sumnum or intnum calls:	$q-1$	$q-1$	$\left(q-1\right)/2$
FFT-step (with fftw):
number of read operations on hard disks:	$q-1$	$2\left(q-1\right)$	$\left(q-1\right)/2$
number of FFTs:	$2$	2	$3$
length of FFTs:	both $q-1$	both $q-1$	one of length $q-1$;
			the others of length $\left(q-1\right)/2$
total RAM usage (in number of long
double positions; in-place FFTs):	$2q+2$	$2q+2$	$2q$222But the computation for ${𝔊}_q^{+}$ requires only $\left(3q+5\right)/2$ long double positions; so, reusing a portion of the RAM after the computation of ${𝔊}_q-{𝔊}_q^{+}$, in the second part of the program we essentially have a gain of about $\left(q-1\right)/2$ long double positions for the RAM usage.

Table 2. Table 2. Values of 𝔊 q subscript 𝔊 𝑞 \mathfrak{G}_{q} , 𝔊 q + superscript subscript 𝔊 𝑞 \mathfrak{G}_{q}^{+} and M q subscript 𝑀 𝑞 M_{q} for every odd prime up to 300 300 300 with a precision of 30 30 30 digits; computed with PARI/GP, v. 2.11.4 with trivial summing over a 𝑎 a . Total computation time: for 𝔊 q subscript 𝔊 𝑞 \mathfrak{G}_{q} , 𝔊 q + superscript subscript 𝔊 𝑞 \mathfrak{G}_{q}^{+} : 18 sec. 852 millisec., for M q subscript 𝑀 𝑞 M_{q} : 19 sec., 171 millisec. on the Dell Optiplex machine mentioned before.

$q$	${𝔊}_q$	${𝔊}_q^{+}$	$M_q$
$3$	$0.94549728087168070323974999415\dots$	$0.57721566490153286060651209008\dots$	$0.36828161597014784263323790407\dots$
$5$	$1.72062421251340476169572878865\dots$	$1.40489514161703774859755907976\dots$	$0.82767947671550488799104698967\dots$
$7$	$2.08759407471733013281542471957\dots$	$1.95715645444971475271382186143\dots$	$0.69374325299917902224231637393\dots$
$11$	$2.41542590428326783034287963583\dots$	$2.66207409890433174906654072453\dots$	$0.64960999942397995363690453077\dots$
$13$	$2.61075773741765019699776108857\dots$	$2.89959572414790509559591203013\dots$	$0.69630986299203715584089218352\dots$
$17$	$3.58197604409757765927178812919\dots$	$3.23179164885108167689200470642\dots$	$1.36293176857311326439833395890\dots$
$19$	$4.79040941571428332590703936458\dots$	$3.36702810226943360422911738361\dots$	$1.56821936415476775304938942269\dots$
$23$	$2.61128917618820092550739164964\dots$	$3.56605274186303485506490005633\dots$	$1.07370241439895666993863022504\dots$
$29$	$3.09373170599426872316275179819\dots$	$3.77451272291818155837540505527\dots$	$1.37173438584080190328583030799\dots$
$31$	$4.31444292526747509770757441042\dots$	$3.74063417131631765163927862231\dots$	$1.41315141911004437078399808370\dots$
$37$	$4.30493818995760201798557926417\dots$	$3.88346103237113739135523493388\dots$	$1.29518958101078356915278401821\dots$
$41$	$3.97152162792133216028257040014\dots$	$3.90067243331576039538420460289\dots$	$1.29673609198958173353796568380\dots$
$43$	$4.37862750574695049413775062336\dots$	$4.37462848511375110150884874389\dots$	$1.41176882240051173489451389181\dots$
$47$	$4.79939425890741613452758429988\dots$	$4.78330592374031492736088514964\dots$	$1.39567565425273602292102717603\dots$
$53$	$4.33773685859709231869696082307\dots$	$4.06734814093911422415451881781\dots$	$1.30627572903790815149667975264\dots$
$59$	$5.43351634538500398077634438193\dots$	$5.74977495098717868985714511291\dots$	$1.81899383678937843989348366929\dots$
$61$	$5.07108519057651619595805098113\dots$	$4.71919160448137601223479232791\dots$	$1.41809980889441627035459190983\dots$
$67$	$5.29213930662896260873428461831\dots$	$5.49478574409231087894450914285\dots$	$1.67019193303154369921782607634\dots$
$71$	$5.25525819281894616772013128637\dots$	$5.02459221437013823603453457463\dots$	$1.47455511100236771011015896767\dots$
$73$	$4.06694909044749529201648815625\dots$	$5.56638018904420607773144876527\dots$	$1.78248970799598673447282517891\dots$
$79$	$4.99827631817068010789431392945\dots$	$4.31392816983842153234814442952\dots$	$1.34616837027813468918588610688\dots$
$83$	$3.03313611343607418716403819105\dots$	$4.06119890648015486954960478374\dots$	$1.34527786237910789501875868023\dots$
$89$	$4.16409079888983276880841110372\dots$	$5.44834851555434719261902953243\dots$	$1.61654649274126300156782088673\dots$
$97$	$4.89124074040389666830751468857\dots$	$4.44563411256346738186380452664\dots$	$1.60286118570076458480362218799\dots$
$101$	$5.29701289150966971887860032739\dots$	$5.93364557387726998305789899164\dots$	$1.51871979857079618912367283335\dots$
$103$	$5.14433955125208822113330503220\dots$	$5.53312508630999898815400644939\dots$	$1.56072764165486011343921965820\dots$
$107$	$5.45827420997024503421680245453\dots$	$5.35744691959596839332603590620\dots$	$1.55529418086936504978552066530\dots$
$109$	$6.90663814626423653219469837704\dots$	$6.28639312060842026587282318484\dots$	$1.65357828827908326582841136643\dots$
$113$	$4.02173038257803067578318006617\dots$	$4.71308052553071355344451609738\dots$	$1.51486982889352164427060492878\dots$
$127$	$5.08859912415333449423215636240\dots$	$5.28427526641642291108714895825\dots$	$1.55590143040596443193792941854\dots$
$131$	$2.83682634158837909860285797321\dots$	$4.29182422162389365669036230041\dots$	$1.43797882292531602089564238879\dots$
$137$	$4.93700022614368468691962999711\dots$	$5.17281966401368126952267004684\dots$	$1.53929870904867707257469538680\dots$
$139$	$5.88916863399867186726383730369\dots$	$5.15673467267785693456200640445\dots$	$1.58828875478913218915240825692\dots$
$149$	$5.98342477769515981450242785739\dots$	$6.35744273145487616682151978517\dots$	$1.55933423387754689170927007457\dots$
$151$	$5.04201611352872179914519461022\dots$	$5.66732269410388218441768644382\dots$	$1.48171078244888795642226012230\dots$
$157$	$7.40802206572222729350845201390\dots$	$5.67766459100970078752076942990\dots$	$1.52915091159611605159149879696\dots$
$163$	$5.92966482288720678755499913844\dots$	$5.54289611872522541669860167904\dots$	$2.16832712928352380386400324642\dots$
$167$	$8.03300175268872470467583357802\dots$	$6.80394798958259907108839110755\dots$	$1.56607236656750344030293511154\dots$
$173$	$3.38434753653206190344297798897\dots$	$4.74313680866654143318864467269\dots$	$1.54242401828716131644723995819\dots$
$179$	$3.86236132549903008112126130282\dots$	$5.59074764196693719810304550344\dots$	$1.60085064594072009293300914735\dots$
$181$	$5.14111848776848135810136664257\dots$	$5.52401113238735460988935254057\dots$	$1.65656567095010010041093792977\dots$
$191$	$4.69286990201422664003552434812\dots$	$6.21621633683078754687889560801\dots$	$1.69400806335478035992195123369\dots$
$193$	$5.16342219673915483320078262720\dots$	$6.33516880970302226248749231989\dots$	$1.72106839151430000218016220949\dots$
$197$	$7.55148715896640647886485129372\dots$	$6.72431280547758930911931614898\dots$	$1.58425224704856913591906318269\dots$
$199$	$6.47366513609320738699497459778\dots$	$4.97867314026834059118807347477\dots$	$1.52055512030192431037107983792\dots$
$211$	$7.73613578424586162532810587585\dots$	$5.43928767077706865027592727891\dots$	$1.58887689723521687477342354947\dots$
$223$	$7.81777971785991367471336734851\dots$	$6.97640718267880419790301145060\dots$	$1.57809439787964273689310796956\dots$
$227$	$8.08053156951296218697071193757\dots$	$6.16478105833535800088839052312\dots$	$1.61440476278289514090073256762\dots$
$229$	$7.16298632058099546745778115058\dots$	$5.19368182825228459062582716349\dots$	$1.64391627222705529854073112016\dots$
$233$	$3.11948354485127541303115295258\dots$	$5.48268694035180653761326391137\dots$	$1.56534808865669695863593307680\dots$
$239$	$3.99911017207833249512632297919\dots$	$4.89826038220509731091188200357\dots$	$1.83593237895342242137799671838\dots$
$241$	$6.03752521401034215065709250935\dots$	$6.91099570349028181262249488655\dots$	$1.74483502309356231328685290592\dots$
$251$	$5.04313708502347351042811119022\dots$	$5.85522475367262429906377535883\dots$	$1.60634233356394595761434310531\dots$
$257$	$8.16991391232741391670225155227\dots$	$7.41413126491779482941571986652\dots$	$1.52986363395322517571321794433\dots$
$263$	$7.30343624736815435414348077406\dots$	$6.88761891078185993452639437420\dots$	$1.61873689910065712561008039262\dots$
$269$	$6.26034831666577102735252755712\dots$	$6.33572466741282346876839833227\dots$	$1.58662353583078976012953348699\dots$
$271$	$5.97717804854803304223773905976\dots$	$4.91607375378349595312704873315\dots$	$1.51145118046000075647340279932\dots$
$277$	$4.59280817714077895164777081661\dots$	$6.07306330239530923314413596279\dots$	$1.72974155675277125427451583060\dots$
$281$	$4.66496432366211457505220852623\dots$	$4.99043740542558229612252801406\dots$	$1.60536366070704717918242357661\dots$
$283$	$7.15028579741068251409225231188\dots$	$7.04969230270522888347459792033\dots$	$1.55609186296142373233316514603\dots$
$293$	$3.38438152121953978658468259238\dots$	$5.38438152121953978658468259238\dots$	$1.58515317244284064528356780036\dots$

Table 3. Table 3. A few other values of 𝔊 q subscript 𝔊 𝑞 \mathfrak{G}_{q} and 𝔊 q + superscript subscript 𝔊 𝑞 \mathfrak{G}_{q}^{+} with a precision of 30 30 30 digits; computed with PARI/GP, v. 2.11.4 with trivial summing over a 𝑎 a [m = minutes, s = seconds]. Computation performed on the Dell Optiplex machine mentioned before.

$q$	${𝔊}_q$	${𝔊}_q^{+}$	time	time	time
			$T$	$S$	direct
$1009$	$8.4421351518492992758606946727\dots$	$6.2733540844322103172186250111\dots$	5s.	3s.	14s.
$2003$	$5.7934213690793633280384982162\dots$	$6.9935258611413978746616842142\dots$	10s.	7s.	39s.
$3001$	$8.6474651369683869388023453509\dots$	$8.6459700672984138998934976577\dots$	17s.	11s.	1m. 11s.
$4001$	$7.0034355462031439943568517684\dots$	$8.7805380094230735872867993849\dots$	24s.	17s.	1m. 49s.
$5003$	$5.5492930045816142277368795404\dots$	$7.2440224742791062634412330617\dots$	32s.	23s.	2m. 36s.
$6007$	$8.3116101219984838165629034403\dots$	$9.8742666472425769486896123420\dots$	41s.	30s.	3m. 22s.
$7001$	$8.5052778761008771393168780384\dots$	$9.6833327734910786447084880544\dots$	52s.	38s.	4m. 07s.
$8009$	$11.6868463915493575353450869960\dots$	$11.4431421556247084876087109206\dots$	1m. 03s.	47s.	5m. 00s.
$9001$	$10.1094784318383409358225035802\dots$	$9.4868388831454962767492760006\dots$	1m. 15s.	57s.	5m. 56s.
$10007$	$12.6646120045606923275389356783\dots$	$11.0601624759024741933308283063\dots$	1m. 27s.	1m. 07s.	7m. 12s.
$20011$	$10.7996803112999205186430402899\dots$	$10.5489807692170969459672226221\dots$	4m. 30s.	3m. 43s.	20m. 01s.
$30011$	$10.3330799721240242255136062255\dots$	$11.0127039500540893278498877674\dots$	9m. 19s.	8m. 11s.	37m. 28s.

Table 4. Table 4. A few other values of 𝔊 q subscript 𝔊 𝑞 \mathfrak{G}_{q} and 𝔊 q + superscript subscript 𝔊 𝑞 \mathfrak{G}_{q}^{+} ; computed with PARI/GP, v. 2.11.4 and fftw , v. 3.3.8, with long double precision. The sum over a 𝑎 a was performed using the FFT algorithm on the Dell Optiplex machine mentioned before [s = seconds, ms = milliseconds; precomputations of decimated in frequency S 𝑆 S -values performed on the Optiplex; their computation time is excluded from this table]. ∗ : on the Intel Xeon machine due to a runtime memory error on the Dell Optiplex.

$q$	${𝔊}_q$	${𝔊}_q^{+}$	time
$10007$	$12.664612\dots$	$11.060162\dots$	10ms.
$20011$	$10.799680\dots$	$10.548981\dots$	23ms.
$30011$	$10.333080\dots$	$11.012704\dots$	15ms.
$40009$	$13.146885\dots$	$13.469520\dots$	25ms.
$42611$	$2.499688\dots$	$8.367404\dots$	41ms.
$50021$	$9.910507\dots$	$11.063741\dots$	98ms.∗
$60013$	$12.810360\dots$	$12.671109\dots$	36ms.
$70001$	$12.572765\dots$	$13.428551\dots$	25ms.
$80021$	$14.185633\dots$	$11.617216\dots$	100ms.∗
$90001$	$11.819424\dots$	$9.601757\dots$	33ms.
$100003$	$15.166074\dots$	$14.765926\dots$	69ms.
$305741$	$1.650523\dots$	$8.839799\dots$	198ms.
$1000003$	$17.379970\dots$	$15.298449\dots$	876ms.
$4178771$	$0.922855\dots$	$8.909168\dots$	2s. 613ms.
$6766811$	$1.604045\dots$	$10.961044\dots$	4s. 584ms.
$10000019$	$17.087945\dots$	$15.974742\dots$	6s. 361ms.
$28227761$	$2.361562\dots$	$10.153369\dots$	17s. 996ms.
$75743411$	$2.469939\dots$	$12.234097\dots$	2m. 24s. 217ms.

Table 5. Table 5. A few other values of 𝔊 q subscript 𝔊 𝑞 \mathfrak{G}_{q} and 𝔊 q + superscript subscript 𝔊 𝑞 \mathfrak{G}_{q}^{+} ; computed with PARI/GP, v. 2.11.4 and fftw , v. 3.3.8, with long double precision. Boldfaced results are the ones corresponding to known instances of 𝔊 q < 0 subscript 𝔊 𝑞 0 \mathfrak{G}_{q}<0 . The sum over a 𝑎 a was performed using the FFT algorithm on the Intel Xeon machine or, for q = 251160191 𝑞 251160191 q=251160191 , 212634221 212634221 212634221 , 1139803271 1139803271 1139803271 , 7079770931 7079770931 7079770931 , 9109334831 9109334831 9109334831 , 9854964401 9854964401 9854964401 on the CAPRI infrastructure mentioned before. [m = minutes, s = seconds; precomputations of decimated in frequency S 𝑆 S -values performed on the cluster; their computation time is excluded from this table].

$q$	${𝔊}_q$	${𝔊}_q^{+}$	time
$193894451$	$0.662110\dots$	$9.607705\dots$	4m. 29s.
$212634221$	$1.435141\dots$	$11.883540\dots$	4m. 28s.
$251160191$	$1.912681\dots$	$11.785574\dots$	2m. 53s.
$538906601$	$1.474911\dots$	$12.957235\dots$	11m. 56s.
$\mathrm{𝟗𝟔𝟒𝟒𝟕𝟕𝟗𝟎𝟏}$	$-0.182374\dots$	$10.402224\dots$	23m. 13s.
$1139803271$	$0.768538\dots$	$8.313111\dots$	27m. 56s.
$1217434451$	$0.877596\dots$	$12.946690\dots$	29m. 16s.
$1806830951$	$0.880396\dots$	$11.973128\dots$	47m. 48s.
$2488788101$	$0.424880\dots$	$12.248837\dots$	103m. 08s.
$2830676081$	$1.254528\dots$	$12.438044\dots$	89m. 59s.
$2918643191$	$0.302793\dots$	$12.573983\dots$	87m. 49s.
$7079770931$	$1.544698\dots$	$14.301772\dots$	742m. 09s.
$\mathrm{𝟗𝟏𝟎𝟗𝟑𝟑𝟒𝟖𝟑𝟏}$	$-0.248739\dots$	$12.128187\dots$	311m. 28s.
$\mathrm{𝟗𝟖𝟓𝟒𝟗𝟔𝟒𝟒𝟎𝟏}$	$-0.096465\dots$	$12.807752\dots$	326m. 03s.

Table 6. Table 6. Computation of the generalised Euler constants γ n subscript 𝛾 𝑛 \gamma_{n} , 0 ≤ n ≤ 30 0 𝑛 30 0\leq n\leq 30 , with a precision of at least 40 digits; computed with PARI/GP, v. 2.11.4.

$n$	${\gamma }_n$
$0$	$0.5772156649015328606065120900824024310\dots$
$1$	$-0.0728158454836767248605863758749013191\dots$
$2$	$-0.0096903631928723184845303860352125293\dots$
$3$	$0.0020538344203033458661600465427533842\dots$
$4$	$0.0023253700654673000574681701775260680\dots$
$5$	$0.0007933238173010627017533348774444448\dots$
$6$	$-0.0002387693454301996098724218419080042\dots$
$7$	$-0.0005272895670577510460740975054788582\dots$
$8$	$-0.0003521233538030395096020521650012087\dots$
$9$	$-0.0000343947744180880481779146237982273\dots$
$10$	$0.0002053328149090647946837222892370653\dots$
$11$	$0.0002701844395439035266729020820679556\dots$
$12$	$0.0001672729121051401933535015433411834\dots$
$13$	$-0.0000274638066037601588600076036933551\dots$
$14$	$-0.0002092092620592999458371396973445849\dots$
$15$	$-0.0002834686553202414466429344749971269\dots$
$16$	$-0.0001996968583089697747077845632032403\dots$
$17$	$0.0000262770371099183366994665976305101\dots$
$18$	$0.0003073684081492528265927547519486256\dots$
$19$	$0.0005036054530473556290555964377171600\dots$
$20$	$0.0004663435615115594494005948244335505\dots$
$21$	$0.0001044377697560001158107956743677204\dots$
$22$	$-0.0005415995822039977016551961731741055\dots$
$23$	$-0.0012439620904082457792997415995371658\dots$
$24$	$-0.0015885112789035615619061966115211158\dots$
$25$	$-0.0010745919527384888247242919873531730\dots$
$26$	$0.0006568035186371544315047730033562152\dots$
$27$	$0.0034778369136185382090073595742588115\dots$
$28$	$0.0064000685317006294581072282219458636\dots$
$29$	$0.0073711517704722391344124024235594021\dots$
$30$	$0.0035577288555731609479135377489084026\dots$

Equations

$${\zeta }_K(s)=\frac{c_{-1}}{s-1}+c_0+\mathcal{O}\left(s-1\right),$$

$$\lim_{s\to 1}\Bigl{(}\frac{\zeta_{K}(s)}{c_{1}}-\frac{1}{s-1}\Bigr{)}=\frac{c_{0}}{c_{-1}}.$$

$${\mathfrak{G}}_q:=\gamma +\sum _{\chi \ne {\chi }_0}\frac{L^{\prime }}{L}(1,\chi ),$$

$${\mathfrak{G}}_q^{+}:=\gamma +\sum _{\begin{array}{c}\chi \ne {\chi }_0\\ \chi \text{even}\end{array}}\frac{L^{\prime }}{L}(1,\chi ).$$

$$M_{q}:=\max_{\chi\neq\chi_{0}}\Bigl{|}\frac{L^{\prime}}{L}(1,\chi)\Bigr{|},$$

$${\gamma }_k(a,q)$$

$$\displaystyle=-\frac{1}{q}\Bigl{(}\frac{(\log q)^{k+1}}{k+1}+\sum_{n=0}^{k}\binom{k}{n}(\log q)^{k-n}\psi_{n}\bigl{(}\frac{a}{q}\bigr{)}\Bigr{)},$$

$$\psi_{n}(z):=-\gamma_{n}-\frac{(\log z)^{n}}{z}-\sum_{m=1}^{+\infty}\Bigl{(}\frac{(\log(m+z))^{n}}{m+z}-\frac{(\log m)^{n}}{m}\Bigr{)}$$

$$\gamma_{n}:=\lim_{N\to+\infty}\Bigl{(}\sum_{j=1}^{N}\frac{(\log j)^{n}}{j}-\frac{(\log N)^{n+1}}{n+1}\Bigr{)}=\sum_{m=1}^{+\infty}\Bigl{(}\frac{(\log m)^{n}}{m}-\frac{(\log(m+1))^{n+1}-(\log m)^{n+1}}{n+1}\Bigr{)},$$

$$L^{\prime }(1,\chi )=-\sum _{a=1}^{q-1}\chi (a){\gamma }_1(a,q),$$

$$\gamma_{1}(a,q)=-\frac{1}{q}\Bigl{(}\frac{1}{2}(\log q)^{2}+\log q\ \psi\bigl{(}\frac{a}{q}\bigr{)}+\psi_{1}\bigl{(}\frac{a}{q}\bigr{)}\Bigr{)},$$

$$T(x):=\gamma_{1}+\psi_{1}(x)=-\frac{\log x}{x}-\sum_{m=1}^{+\infty}\Bigl{(}\frac{\log(x+m)}{x+m}-\frac{\log m}{m}\Bigr{)},$$

$$\gamma_{1}=\lim_{N\to+\infty}\Bigl{(}\sum_{j=1}^{N}\frac{\log j}{j}-\frac{(\log N)^{2}}{2}\Bigr{)}=-0.0728158454835\dotsc$$

$$L(1,\chi)=-\frac{1}{q}\sum_{a=1}^{q-1}\chi(a)\ \psi\bigl{(}\frac{a}{q}\bigr{)},$$

$$L^{\prime}(1,\chi)=-(\log q)L(1,\chi)+\frac{1}{q}\sum_{a=1}^{q-1}\chi(a)\ T\bigl{(}\frac{a}{q}\bigr{)},$$

$$\sum _{\chi \ne {\chi }_0}\frac{L^{\prime }}{L}(1,\chi )=-(q-2)\log q-\sum _{\chi \ne {\chi }_0}\frac{{\sum }_{a=1}^{q-1}\chi (a)T(a/q)}{{\sum }_{a=1}^{q-1}\chi (a)\psi (a/q)}.$$

$$L(s,\chi)=\frac{1}{\pi i}\Bigl{(}\frac{2\pi}{q}\Bigr{)}^{s}\Gamma(1-s)\frac{\tau(\chi)}{\sqrt{q}}\cos\Bigl{(}\frac{\pi s}{2}\Bigr{)}L(1-s,\overline{\chi})$$

$$\frac{L^{\prime}}{L}(s,\chi)=\log\Bigl{(}\frac{2\pi}{q}\Bigr{)}-\frac{\Gamma^{\prime}}{\Gamma}(1-s)-\frac{\pi}{2}\tan\Bigl{(}\frac{\pi s}{2}\Bigr{)}-\frac{L^{\prime}}{L}(1-s,\overline{\chi}),$$

$$\frac{L^{\prime}}{L}(0,\chi)=\log\Bigl{(}\frac{2\pi}{q}\Bigr{)}+\gamma-\frac{L^{\prime}}{L}(1,\overline{\chi}).$$

$$L^{\prime }(0,\chi )$$

$$\displaystyle=\frac{\log q}{q}\sum_{a=1}^{q-1}a\chi(a)+\sum_{a=1}^{q-1}\chi(a)\log\Bigl{(}\Gamma\bigl{(}\frac{a}{q}\bigr{)}\Bigr{)}$$

$$\displaystyle=-(\log q)L(0,\chi)+\sum_{a=1}^{q-1}\chi(a)\log\Bigl{(}\Gamma\bigl{(}\frac{a}{q}\bigr{)}\Bigr{)},$$

$$L(0,\chi )=-B_{1,\chi }:=-\frac{1}{q}\sum _{a=1}^{q-1}a\chi (a)\ne 0,$$

$$\sum _{\chi \text{odd}}\frac{L^{\prime }}{L}(1,\chi )$$

$$L^{\prime}(1,\chi)=(\gamma+\log(2\pi))L(1,\chi)+\frac{\tau(\chi)}{q}\sum_{a=1}^{q-1}\overline{\chi}(a)\ R\bigl{(}\frac{a}{q}\bigr{)},$$

$$R(x)$$

$$S(x)$$

$$S(x)=2\int_{0}^{+\infty}\Bigl{(}(x-1)e^{-t}+\frac{e^{-xt}-e^{-t}}{1-e^{-t}}\Bigr{)}\frac{\gamma+\log t}{t}\ \mathrm{d}t,\quad x>0.$$

$$\sum _{a=1}^{q-1}\overline{\chi }(a)R(a/q)=-\sum _{a=1}^{q-1}\overline{\chi }(a)S(a/q).$$

$$L(1,\chi)=2\frac{\tau(\chi)}{q}\sum_{a=1}^{q-1}\overline{\chi}(a)\log\Bigl{(}\Gamma\bigl{(}\frac{a}{q}\bigr{)}\Bigr{)},$$

$$\sum _{\begin{array}{c}\chi \ne {\chi }_0\\ \chi \text{even}\end{array}}\frac{L^{\prime }}{L}(1,\chi )$$

$$\mathfrak{G}_{q}^{+}=\frac{q-1}{2}\gamma+\frac{q-3}{2}\log(2\pi)-\frac{1}{2}\sum_{\begin{subarray}{c}\chi\neq\chi_{0}\\ \chi\,\textrm{even}\end{subarray}}\frac{\sum_{a=1}^{q-1}\overline{\chi}(a)\ S(a/q)}{\sum_{a=1}^{q-1}\overline{\chi}(a)\log\bigl{(}\Gamma(a/q)\bigr{)}}.$$

Full text

1. Introduction

Let $K$ be a number field and let $\zeta_{K}(s)$ be its Dedekind zeta-function. It is a well known fact that $\zeta_{K}(s)$ has a simple pole at $s=1$; writing the expansion of $\zeta_{K}(s)$ near $s=1$ as

[TABLE]

the Euler-Kronecker constant of $K$ is defined as

[TABLE]

In the special case in which $K=\mathbb{Q}(\zeta_{q})$ is a prime cyclotomic field, where $q$ is an odd prime and $\zeta_{q}$ is a primitive $q$-root of unity, we have that the Dedekind zeta-function satisfies $\zeta_{\mathbb{Q}(\zeta_{q})}(s)=\zeta(s)\prod_{\chi\neq\chi_{0}}L(s,\chi)$, where $\zeta(s)$ is the Riemann zeta-function, $L(s,\chi)$ are the Dirichlet $L$-functions, $\chi$ runs over the non trivial Dirichlet characters mod $q$ and $\chi_{0}$ is the trivial Dirichlet character mod $q$. By logarithmic differentiation, we immediately get that the Euler-Kronecker constant for the prime cyclotomic field $\mathbb{Q}(\zeta_{q})$ is

[TABLE]

where $\gamma$ is the Euler-Mascheroni constant. Sometimes the quantity $\mathfrak{G}_{q}$ is denoted as $\gamma_{q}$ but this conflicts with notations used in literature.

An extensive study of the properties of $\mathfrak{G}_{q}$ was recently started by Ihara [15, 16] and continued by many others; here we are mainly interested in computational problems involving $\mathfrak{G}_{q}$ and hence we just recall the paper by Ford-Luca-Moree [9]. We introduce a new method to compute the Euler-Kronecker constants of prime cyclotomic fields which is faster and uses less computing resources than the ones previously known. The new algorithm requires the values of the generalised gamma function $\Gamma_{1}$ at some rational arguments $a/q\in(0,1)$. Such a function has, for $q$ large and $a={o}\left(q\right)$, an order of magnitude exponentially smaller than the ones previously used to determine $\mathfrak{G}_{q}$, see Section 4 below. Moreover, the Fast Fourier Transform (FFT) used in this new approach allows a decimation in frequency strategy111We use here this nomenclature since it is standard in the literature on the Fast Fourier Transform, but it could be translated in number theoretic language using suitable properties of Dirichlet characters. that leads to gaining a factor $1/2$ in the quantity of needed precomputation operations, in the length of the involved transforms and in their memory usage. Our algorithm uses the formulae in Section 3 below that, according to Deninger [6] and Kanemitsu [14], were first proved in 1883 by Berger [1] and in 1929 by Gut [13].

Another interesting quantity related to $\mathfrak{G}_{q}$ is the Euler-Kronecker constant $\mathfrak{G}_{q}^{+}$ for $\mathbb{Q}(\zeta_{q}+\zeta_{q}^{-1})$, the maximal real subfield of $\mathbb{Q}(\zeta_{q})$. According to eq. (10) of Moree [20] it is defined as

[TABLE]

In Section 3.4 we will give a formula that leads a to a direct evaluation of $\mathfrak{G}_{q}^{+}$ in terms of some special functions values attained at some rationals $a/q\in(0,1)$. Moreover, in Section 3.4 we will use the previously proved relations to see why the quantity $\mathfrak{G}_{q}-\mathfrak{G}_{q}^{+}$ is much easier to compute than both its summands.

During such computations, as a by-product, we also evaluated the related quantity

[TABLE]

see Section 5. Other quantities related to $\mathfrak{G}_{q}$ are the generalised Euler constants in arithmetic progressions, sometimes also called Stieltjes constants in arithmetic progressions, denoted as $\gamma_{k}(a,q)$, $k\in\mathbb{N}$, $q\geq 1$, $1\leq a\leq q$, which are defined by

[TABLE]

see eq. (1.3)-(1.4) and (7.3) of Dilcher [7], where

[TABLE]

for $n\in\mathbb{N}$ and $z\in\mathbb{C}\setminus\{0,-1,-2,\dotsc\}$, $\psi_{n}(1)=-\gamma_{n}$, and the generalised Euler constants $\gamma_{n}$ are defined as

[TABLE]

by, e.g., eq. (3)-(4) of Bohman-Fröberg [3]. Remark that $\gamma_{0}=\gamma$. It is worth recalling that the functions $\psi_{n}(z)$ occur in Ramanujan’s second notebook, see [2, Chapter 8, Entry 22].

The quantities in (4) and, as we will see in sections 2-3 below, the one in (1), are hence connected with the values of $\psi_{n}$, $n\geq 1$, which is the logarithmic derivative of $\Gamma_{n}$, a generalised Gamma function, see Deninger [6], Dilcher [8] and Katayama [18], whose definition for $n=1$ is given in Section 3.2. In some sense we can say that the $\psi_{n}$-functions, $n\geq 1$, are the analogue of the usual digamma function. In the following we will denote as $\psi$ the standard digamma function $\Gamma^{\prime}/\Gamma$; we also remark that it can be represented as the function $\psi_{0}$ defined in (5).

The paper is organised as follows. In Sections 2-3 we will give the derivations of the main formulae we will need in the computations; such proofs are classical and are based on the functional equation for the Dirichlet $L$-functions, see Cohen’s books [4]-[5], for instance. Other useful references for this part are the papers of Deninger [6] and Dilcher [7, 8].

In Section 4 we will see how to implement the formulae of the previous two sections, starting form a straightforward application of the definitions (1)-(2) of $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$; then we will compare the Ford-Luca-Moree approach, based on the formulae of Section 2, with our new procedure, based on the formulae of Section 3. In particular we will see, in both cases, how to insert the Fast Fourier Transform and we will discuss their precisions, computational costs and memory usages. In Section 5 we describe how to compute $M_{q}$. Section 6 is devoted to provide more data on the generalised Euler constants in arithmetic progressions. Finally, in Section 7 we will provide several tables containing a comparison scheme of the different implementations and the computational results and running times of the different approaches. At the bottom of the paper we inserted two colored scatter plots for the normalised values of $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$ for every prime $q$, $3\leq q\leq 10^{6}$ and two scatter plots about $M_{q}$ and its normalised values for the same set of primes.

We finally remark that some of the ideas presented here will also be used in a joint work with Pieter Moree, Sumaia Saad Eddin and Alisa Sedunova on the computation of the Kummer ratio of the class number for prime cyclotomic fields, see [19].

Acknowledgements. Some of the calculations here described were performed using the University of Padova Strategic Research Infrastructure Grant 2017: “CAPRI: Calcolo ad Alte Prestazioni per la Ricerca e l’Innovazione”, http://capri.dei.unipd.it. I also wish to thank Karim Belabas and Bill Allombert (University of Bordeaux) for a couple of key suggestions about libpari and gp2c and Luca Righi (University of Padova) for his help in developing the quadruple precision versions of the fft-programs, in designing the parallelised precomputations and in organising the use of the cluster of the Dipartimento di Matematica “Tullio Levi-Civita”, http://computing.math.unipd.it/highpc, and the use of CAPRI.

2. The Ford-Luca-Moree method

Recall that $q$ is an odd prime. If we do not restrict to Dirichlet characters of a prescribed parity, we can use eq. (6.1) and (7.4) of Dilcher [7], as in Ford-Luca-Moree, see eq. (3.2) in [9]. In fact eq. (6.1) of [7] gives

[TABLE]

where $\gamma_{1}(a,q)$ is defined in (4) which, for $k=1$, becomes

[TABLE]

for any $q\geq 1$ and $1\leq a\leq q$, where $\psi,\psi_{1}$ are defined in (5). Again using (5), we define

[TABLE]

and, specialising (6), we also have

[TABLE]

To compute $\gamma_{1}$ and similar constants with a very large precision, see Section 6.3 below. We also remark here that the rate of convergence of the series in (7) is, roughly speaking, about $(\log m)/m^{2}$. Recalling now eq. (3.1) of [9], i.e.,

[TABLE]

by the orthogonality of Dirichlet characters and (8), we obtain eq. (3.2) of [9], i.e.,

[TABLE]

where $T(x)$ is defined in (7) (pay attention to the change of sign in (7) with respect to eq. (3.2) of [9]). Summarising, we finally get

[TABLE]

Formula (9) is the one used in the paper by Ford-Luca-Moree [9]. We will now explain how we can compute $\mathfrak{G}_{q}$ via (1) using the values of the two special functions $\psi$ and $T$, together with the values of the non trivial Dirichlet characters mod $q$.

From a computational point of view it is clear that in (9) we first have to evaluate $T(a/q)$ and $\psi(a/q)$ for every $1\leq a\leq q-1$. For the $\psi$-values we can rely on the PARI/GP function psi or, if less precision is sufficient, we can use the analogous function included in GSL, the gnu scientific library [11]222GSL provides just a double precision (in the sense of the C programming language precision) version of $\psi$; hence this is faster, but less accurate, than the computation of the $\log\Gamma$-values needed in the procedure described in the next section. If we use PARI/GP to precompute and store the $\psi$-values, then the costs of the precomputation and the input/output part of the FFT step have to be doubled, see Table 1.. For computing the $T$-values, a task for which there are no pre-defined functions in any software libraries we know, we can use the summing function sumnum of PARI/GP; this is the most time-consuming step of the procedure. Using the FFT algorithm to perform the sums over $a$, it is easy to see that computing $\mathfrak{G}_{q}$ via (9) has a computational cost of $\mathcal{O}\left(q\log q\right)$ arithmetical operations together with the cost of computing $q-1$ values of the $\psi$ and $T$ functions. For more details see Section 4.

3. Another method: distinguishing Dirichlet characters’ parities

3.1. Primitive odd Dirichlet character case.

Recall that $q$ is an odd prime, let $\chi\neq\chi_{0}$ be a primitive odd Dirichlet character mod $q$ and let $\tau(\chi):=\sum_{a=1}^{q}\chi(a)\,e(a/q)$, $e(x):=\exp(2\pi ix)$, be the Gauß sum associated with $\chi$. The functional equation for $L(s,\chi)$, see, e.g., the proof of Theorem 3.5 of Gun-Murty-Rath [12], gives

[TABLE]

and hence

[TABLE]

which, evaluated at $s=0$, gives

[TABLE]

By the Lerch identity about values of the Hurwitz zeta-function, see, e.g., Proposition 10.3.5 of Cohen [5], and the orthogonality of Dirichlet characters, we get

[TABLE]

since, see Proposition 9.5.12 and Corollary 10.3.2 of Cohen [5], we have

[TABLE]

where $B_{1,\chi}$ is the first $\chi$-Bernoulli number which is non-zero since $\chi$ is odd.

Summarising, by (10)-(12), we obtain

[TABLE]

From a computational point of view, in (13) we need to compute the $\log\Gamma$-values; to do so we can rely on an internal PARI/GP function or, if less precision is sufficient, we can use the analogous function included in the C programming language. We remark that, for $x\to 0^{+}$, \log\bigl{(}\Gamma(x)\bigr{)}\sim\log(1/x) and $\psi(x)\sim-1/x$; hence for $q$ large and $a={o}\left(q\right)$, the values of \log\bigl{(}\Gamma(a/q)\bigr{)} are exponentially smaller than the ones of $\psi(a/q)$. Moreover, to compute the first $\chi$-Bernoulli number $B_{1,\overline{\chi}}$, defined in (12), we just need an integral sequence.

3.2. Primitive even Dirichlet character case.

Recall that $q$ is an odd prime. Assume now that $\chi\neq\chi_{0}$ is a primitive even Dirichlet character mod $q$. We follow Deninger’s notation in [6] by writing $R(x)=-\frac{\partial^{2}}{\partial s^{2}}\zeta(s,x)|_{s=0}=\log(\Gamma_{1}(x))$, $x>0$, where $\zeta(s,x)$ is the Hurwitz zeta function, $s\in\mathbb{C}\setminus\{1\}$. By eq. (3.5)-(3.6) of [6] we have

[TABLE]

where, see eq. (2.3.2) of [6], the $R$-function can be expressed for every $x>0$ by

[TABLE]

It is worth recalling that comparing (15)-(16) with (7), we see that $\psi_{1}(x)=R^{\prime}(x)/2$ (note the different definition of $\gamma_{1}$ on page 174 of Deninger’s paper). Using (6) we have $S(1)=0$ and $R(1)=-\zeta^{\prime\prime}(0)$. An alternative definition of $S(x)$ for $x>0$, which will be useful during the computations, is implicitly contained in eq. (2.12) of Deninger [6]:

[TABLE]

By the orthogonality of the Dirichlet characters, we immediately get

[TABLE]

For $L(1,\chi)$, we use formula (2) of Proposition 10.3.5 of Cohen [5] and the parity of $\chi$ to get

[TABLE]

since $W(\chi)=\tau(\chi)/q^{1/2}$ for even Dirichlet characters, see Definition 2.2.25 of Cohen [4]. Summarising, using (14) and (18)-(19), if $\chi$ is an even Dirichlet character mod $q$, we finally get

[TABLE]

We remark that in (20) we can reuse the $\log\Gamma$-values already needed in (13). For computing the $S$-values, a task for which there are no pre-defined functions in any software libraries we know, we can use the PARI/GP functions sumnum and intnum; this is the most time-consuming step of the procedure. We also remark that, for $x\to 0^{+}$, $S(x)\sim(\log x)^{2}$ and $T(x)\sim\log(1/x)/x$; hence for $q$ large and $a={o}\left(q\right)$, the values of $S(a/q)$ are exponentially smaller than the corresponding ones of $T(a/q)$.

Finally we remark that for $S$ and $B_{1,\overline{\chi}}$ we just need the summation over half of the Dirichlet characters involved; hence in both cases in their computation using the Fast Fourier Transform we can implement the so-called decimation in frequency strategy that allow us to improve on both the speed and the memory usage of the actual computation, see Section 4.1 below. Using the FFT algorithm to perform the sums over $a$, it is easy to see that computing $\mathfrak{G}_{q}$ via (13) and (20) has a computational cost of $\mathcal{O}\left(q\log q\right)$ arithmetical operations together with the cost of computing $q-1$ values of the $\log\Gamma$-function and $(q-1)/2$ decimated in frequency values of the $S$-function333An explanation for this fact can be found at the end of Section 4.1..

3.3. On $\mathfrak{G}_{q}^{+}$: the constant attached to the maximal real subfield

of $\mathbb{Q}(\zeta_{q})$.

It is a consequence of the computations in this section that the Euler-Kronecker constant $\mathfrak{G}_{q}^{+}$ for $\mathbb{Q}(\zeta_{q}+\zeta_{q}^{-1})$, the maximal real subfield of $\mathbb{Q}(\zeta_{q})$, is directly connected with the $S$-function since, by (2) and (20), we have

[TABLE]

Hence in this case the relevant information is encoded in the $S$ and $\log\Gamma$ functions. Clearly $\mathfrak{G}_{q}^{+}$ can be obtained during the $\mathfrak{G}_{q}$-computation since it requires a subset of the data needed for getting $\mathfrak{G}_{q}$. In Figure 2 you can find a colored scatter plot of its values for every $q$ prime, $3\leq q\leq 10^{6}$.

Moreover, a direct evaluation of $\mathfrak{G}_{q}^{+}$ via (21) allow us to use a decimation in frequency strategy in the application of the FFT technique to evaluate the sums over $a$, see Sections 4.1-4.3.

3.4. Regarding $\mathfrak{G}_{q}-\mathfrak{G}_{q}^{+}$.

By (1)-(2), (13) and (21) it is trivial to get that

[TABLE]

This reveals that, from a practical point of view, $\mathfrak{G}_{q}-\mathfrak{G}_{q}^{+}$ is much easier to compute with respect to both $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$: this not just because, as for $\mathfrak{G}_{q}^{+}$, it requires a subset of the data needed for $\mathfrak{G}_{q}$ but also because it involves just one special function, $\log\Gamma$, which is directly available in many software libraries and in the C programming language.

In this case too, a direct evaluation of $\mathfrak{G}_{q}-\mathfrak{G}_{q}^{+}$ via (22) allow us to use a decimation in frequency strategy in the application of the FFT technique to evaluate the sums over $a$, see Sections 4.1-4.3. Some computational data about this quantity are also included in [19].

4. Comparison of methods, results and running times

First of all we notice that PARI/GP, v. 2.11.4, has the ability to generate the Dirichlet $L$-functions (and many other $L$-functions) and hence the computation of $\mathfrak{G}_{q}$, $\mathfrak{G}_{q}^{+}$ and $M_{q}$ can be performed using (1)-(3) with few instructions of the gp scripting language. This computation has a linear cost in the number of calls of the lfun function of PARI/GP and, at least on our Dell Optiplex desktop machine, it is slower than both the procedures we are about to describe.

Comparing (13) and (20) with (9), we see that in both cases we can rely on pre-defined functions to compute either the \log\bigl{(}\Gamma(a/q)\bigr{)}-values or the $\psi(a/q)$-values, $1\leq a\leq q-1$, and finally we have to evaluate the $T$ and $S$ functions respectively involved. We recall that, when taking $q$ very large, it is relevant to know their order of magnitude for $x\to 0^{+}$; it is easy to verify that \log\bigl{(}\Gamma(x)\bigr{)}\sim\log(1/x), $S(x)\sim(\log x)^{2}$, $\psi(x)\sim-1/x$ and $T(x)\sim\log(1/x)/x$. Hence for $x\to 0^{+}$, we have that \log\bigl{(}\Gamma(x)\bigr{)} and $S(x)$ are exponentially smaller than $\psi(x)$ and $T(x)$; a fact that will lead to a more accurate result when using a fixed precision in the final step of the computation. Another difference is that, for the odd Dirichlet characters, the first $\chi$-Bernoulli number in eq. (13) does not involve any special function, but just an integral sequence. So it seems reasonable to compare the following two approaches:

1. a)

use the $T$-series formulae and the $\psi$-values as in [9]; in this case we have two possible alternatives to evaluate the $\psi$-function: using GSL (gaining in speed but losing in precision) or using PARI/GP (with a much better precision, but doubling the needed hard disk storage and the number of input/output operations on the hard disk); 2. b)

use the $S$-function formulae for the even Dirichlet characters case and the first $\chi$-Bernoulli number for the odd one; remark that in both cases we have to evaluate a sum of the $\log\Gamma$-values.

This way we can extend the computation performed in [9], not only because we are developing a different implementation of the same formulae, but also because we can solve the problem in an alternative way which is faster, needs less computing resources, and uses functions having a much smaller order of magnitude, see Table 1 for a summary of these facts. In the computation we will use the PARI/GP scripting language to exploit its ability to accurately evaluate the series and integrals involved in the definition of the $T$ and $S$ functions, defined respectively in (7) and (16)-(17), via the functions sumnum or intnum.

4.1. Using the FFT algorithm

We also remark that the procedures a)-b) trivially require a quadratic number of arithmetical operations to perform the computations in (9), (13) and (20), but this can be improved by using the FFT algorithm and the following argument. Focusing on (9), (13) and (20), we remark that, since $q$ is prime, it is enough to get $g$, a primitive root of $q$, and $\chi_{1}$, the Dirichlet character mod $q$ given by $\chi_{1}(g)=e^{2\pi i/(q-1)}$, to see that the set of the non trivial characters mod $q$ is $\{\chi_{1}^{j}\colon j=1,\dotsc,q-2\}$. Hence, if, for every $k\in\{0,\dotsc,q-2\}$, we denote $g^{k}\equiv a_{k}\in\{1,\dotsc,q-1\}$, every summation in (9)-(13) and (20) is of the type

[TABLE]

where $e(x):=\exp(2\pi ix)$, $j\in\{1,\dotsc,q-2\}$, $\sigma=\pm 1$, and $f$ is a suitable function which assumes real values. As a consequence, such quantities are, depending on $\sigma$, the Discrete Fourier Transforms, or its inverse transformation, of the sequence $\{f(a_{k}/q)\colon k=0,\dotsc,q-2\}$. This idea was first formulated by Rader [22] and it was already used in [9] to speed-up the computation of these quantities via the use of FFT-dedicated software libraries.

For the approach b) we can also use the decimation in frequency strategy: assuming that in (23) one has to distinguish between the parity of $j$ (hence on the parity of the Dirichlet characters), letting $m=(q-1)/2$, for every $j=0,1,\dotsc,q-2$ we have that

[TABLE]

Let now $j=2t+\ell$, where $\ell\in\{0,1\}$ and $t\in\mathbb{Z}$. Then, the previous equation becomes

[TABLE]

where $t=0,\dotsc,m-1$, $\sigma=\pm 1$,

[TABLE]

Hence, if we just need the sum over the even, or odd, Dirichlet characters as in the procedure b) for $f(x)=S(x)$ or $f(x)=x$, instead of computing an FFT transform of length $q-1$ we can evaluate an FFT of length $(q-1)/2$, applied on a suitably modified sequence according to (24). Clearly this represents a gain in both the speed and the memory usage in running the actual computer program. Moreover, if the values of $f(a_{k}/q)$ have to be precomputed and stored, this also means that the quantity of information we have to save during the precomputation (which will be the most time consuming part), and to recall for the FFT algorithm, is reduced by a factor of $2$.

In Table 1 we give a summary of the main characteristics of both approaches for computing $\mathfrak{G}_{q}$; it is clear that the one using $T(x)$ beats the one which implements $S(x)$ only in the total number of the needed FFT transforms444In fact the FFT transforms can be independently performed and hence they can be executed in parallel; this eliminates the unique disadvantage in using the $S$-function method., but in any other aspect the latter is better. In particular the procedure b) is much faster in the precomputation part since its cost is $\leq 1/2$ than approach a)’s one.

4.2. Decimation in frequency for the even Dirichlet characters case

We make explicit the form that the sequence $b_{k}$ defined in (24), assumes in our cases.

It is useful to remark that from $\langle g\rangle=\mathbb{Z}^{*}_{q}$ it trivially follows that $g^{m}\equiv q-1\bmod{q}$, where $m=(q-1)/2$. Hence, recalling $a_{k}\equiv g^{k}\bmod q$, we obtain $a_{k+m}\equiv g^{k+m}\equiv a_{k}(q-1)\equiv q-a_{k}\bmod{q}$ and, as a consequence, we get

[TABLE]

So, inserting the reflection formula for $S(x)$, see eq. (3.3) of Dilcher [8]555Pay attention to the fact that the Deninger $S(x)$-function defined in (15)-(16) is equal to $-2\log(\Gamma_{1}(x))$ as defined in Proposition 1 of Dilcher [8]., into (24)-(25), for every $k=0,\dotsc,m-1$ and for $f(x)=S(x)$, using (16), the sequence $b_{k}$ becomes

[TABLE]

where $a_{k}\equiv g^{k}\bmod q$, while, using (17), we obtain

[TABLE]

in which we exploited the uniform convergence of the involved integrals. To optimise speed and precision, both equations (26)-(27) will be used during the actual computations; when possible we will exploit the exponential decay $e^{-ct}$, with $c=\min(a_{k}/q,1-a_{k}/q)$, of the integrand function in (27) using the PARI/GP function intnum. But when the parameter $c$ will become too small to give reliable results, we will switch to apply the PARI/GP function sumnum to eq. (26); in this case, roughly speaking, the decay order is $(\log n)/n^{2}$.

Hence, thanks to the previous formulae, the number of calls to the sumnum or intnum functions required in the precomputation of the $S$-values is reduced by a factor of $2$ with respect to the ones needed to precompute the $T$-values.

If we are just interested in the computation of $\mathfrak{G}_{q}^{+}$, we can directly use (21) in which we can embed (26)-(27) and the following remark about the needed $\log\Gamma$-values. Assuming $f(x)=\log\Gamma(x)$, using (25) and the well-known reflection formula $\Gamma(x)\Gamma(1-x)=\pi/\sin(\pi x)$, we obtain

[TABLE]

thus further simplifying the final computation by replacing the $\Gamma$-function with the $\sin$-function.

4.3. Decimation in frequency for the odd Dirichlet characters case

We make explicit the form that the sequence $c_{k}$ defined in (24), assumes in our cases.

If we are just interested in the computation of $\mathfrak{G}_{q}-\mathfrak{G}_{q}^{+}$, we can directly use (22); using the reflection formula $\Gamma(x)\Gamma(1-x)=\pi/\sin(\pi x)$ and arguing as in the previous paragraph, we obtain

[TABLE]

for every $k=0,\dotsc,m-1$, $m=(q-1)/2$, and hence $c_{k}$ is modified accordingly. In this case the gain of using the previous formula is that the number of needed evaluations of the $\log\Gamma$-function is reduced by a factor of $2$.

The case in which $f(x)=x$ is easier; using again $\langle g\rangle=\mathbb{Z}^{*}_{q}$, $a_{k}\equiv g^{k}\bmod q$ and $g^{m}\equiv q-1\bmod{q}$, we can write that $a_{k+m}\equiv q-a_{k}\bmod{q}$; hence

[TABLE]

so that in this case we obtain $c_{k}=e(\sigma k/(q-1))(2a_{k}/q-1)$ for every $k=0,\dotsc m-1$, $m=(q-1)/2$, $\sigma=\pm 1$.

4.4. Computations trivially summing over $a$ (slower but with more digits available).

Unfortunately in libpari the FFT-functions work only if $q=2^{\ell}+1$, for some $\ell\in\mathbb{N}$. So we had to trivially perform these summations and hence, in practice, this part is the most time consuming one in both the procedures a) and b) since it has a quadratic cost in $q$. Being aware of such limitations, we used PARI/GP (with the trivial way to compute the sum over $a$) to evaluate $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$ with these three approaches for every odd prime $q\leq 300$, on a Dell OptiPlex-3050 (Intel i5-7500 processor, 3.40GHz, 16 GB of RAM and running Ubuntu 18.04.2) using a precision of $30$ digits, see Table 2; we also inserted there the values of $M_{q}$, defined in (3), for the same set of primes. Such results largely extend the precision of the data in Table 1 on page 1472 of [9]. The computation of the values of Table 2 needed 19 seconds using the $S$-function, 33 seconds using the $T$-function and 51 seconds using PARI/GP lfun function. We also computed the values of $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$, with a precision of 30 digits, for $q=1009,2003,3001,4001,5003,$ $6007,7001,$ $8009,$ $9001,10007$, $20011,30011$, as you can see in Table 3. These numbers were chosen to heuristically evaluate how the computational cost depends on the size of $q$. In this case, in the fifth column of Table 3 we also reported the running time of the direct approach, i.e. using (1), the third and fourth columns are respectively the running times of the other two procedures. For these values of $q$ it became clear that the computation time spent in performing the sums over $a$ was the longest one. This means that inserting an FFT-algorithm is fundamental to further improve the performances of both the algorithms a)-b). We discuss this in more detail in the next paragraph.

4.5. Computations summing over $a$ via FFT (much faster but with less

digits available).

As we saw before, for large $q$ the time spent in summing over $a$ dominates the computational cost. So we implemented the use of FFT for this task. We first used the gp2c compiler tool to obtain suitable C programs to perform the precomputations of the needed $T$ and $S$-values with $38$ digits and save them to the hard disk666If we do not use the GSL to directly compute $\psi$, we need to insert its precomputation here.. Then we passed such values to the C programs which used the fftw [10] software library to perform the FFT step. In such a final stage the performance was thousand times faster than the one for the same stage trivially performed; as an example you can compare the running times for $q=10007$, $20011$, $30011$ in Tables 3 and 4. The running times for the approaches a) and b) reveal that the latter is faster, mainly because it requires less input operations to gain the stored precomputed information since the FFT works on a set of data of half the length than in the former case777If $\psi$ is precomputed using PARI/GP, then the gain ratio in the stored space and in the number of input/output operations is raised to $3/4$..

This way we computed the values of $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$ for $q=40009$, $42611$, $50021$, $60013$, $70001$, $80021$, $90001$, $100003$, $305741$, $1000003$, $4178771$, $6766811$, $10000019$, $28227761$, $75743411$ with the long double precision, see Table 4. These computations were performed with the Dell OptiPlex machine mentioned before.

Some of these $q$-values were chosen for their size and others with the help of $\mathcal{B}$, the “greedy sequence of prime offsets”, http://oeis.org/A135311, in the following way. We define $\mathcal{B}$ using induction, by $b(1)=0\in\mathcal{B}$ and $b(n)\in\mathcal{B}$ if it is the smallest integer exceeding $b(n-1)$ such that for every prime $r$ the set $\{b(i)\bmod r\colon 1\leq i\leq n\}$ has at most $r-1$ elements. An equivalent statement, assuming that the prime $k$-tuples conjecture holds, is that $b(n)$ is minimal such that $b(1)=0$ and there are infinitely many primes $q$ with $b(i)q+1$ prime for $2\leq i\leq n$, $n\geq 2$. Let now

[TABLE]

where $\mathcal{A}$ is an admissible set, i.e., $\mathcal{A}=\{a_{1},\dotsc,a_{s}\}$, $a_{i}\in\mathbb{N}$, $a_{i}\geq 1$, such that does not exist a prime $p$ such that $p\mid n\prod_{i=1}^{s}(a_{i}n+1)$ for every $n\geq 1$. Thanks to Theorem 2 of Moree [20], if the prime $k$-tuples conjecture holds and if $\mathcal{A}$ is an admissible set, then $\mathfrak{G}_{q}<(2-m(\mathcal{A})+{o}\left(1\right))\log q$ for $\gg x/(\log x)^{-|\mathcal{A}|-1}$ primes $q\leq x$. Moreover, by Theorem 6 of Moree [20], assuming both the Elliott-Halberstam and the prime $k$-tuples conjectures, if $\mathcal{A}$ is an admissible set then $\mathfrak{G}_{q}=(1-m(\mathcal{A})+{o}\left(1\right))\log q$ for $\gg x/(\log x)^{-|\mathcal{A}|-1}$ primes $q\leq x$.

The greedy sequence of prime offsets $\mathcal{B}$ has the property that any finite subsequence is an admissible set. With a PARI/GP script we computed the first $2089$ elements of $\mathcal{B}$ since for $\mathcal{C}:=\{b(2),\dotsc,b(2089)\}$ we get $m(\mathcal{C})>2$. So, if we are looking for negative values of $\mathfrak{G}_{q}$, it seems to be a good criterion to evaluate $\mathfrak{G}_{q}$ for a prime number $q$ such that $bq+1$ is prime for many elements $b\in\mathcal{C}$ (clearly it is better to start with the smaller available $b$’s). To be able to measure this fact, we define

[TABLE]

Some of the $q$-values written before in this paragraph are such that $v(q)>1.15$ so that, thanks to Moree’s results already cited, they are good candidates to have a negative Euler-Kronecker constant. The complete list of $q\leq 10^{10}$ such that $v(q)>1.2$ is towards the end of the PARI/GP script testseq that can be downloaded here: http://www.math.unipd.it/~languasc/EK-comput.html.

4.5.1. Data for the scatter plots.

After having evaluated the running times of the previous examples, we decided to provide the colored scatter plots, see Figures 1-2, of the normalised values of $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$ (both in long double precision) for every odd prime $q\leq 10^{6}$ thus enlarging the known range of the data on $\mathfrak{G}_{q}$ and $\mathfrak{G}_{q}^{+}$, see [9]. For performing the needed precomputations of the $S$-values, we used the cluster of the Department of Mathematics of the University of Padova; the cluster setting is described here: http://computing.math.unipd.it/highpc. The minimal value of $\mathfrak{G}_{q}/\log q$, $3\leq q\leq 10^{6}$, $q$ prime, is $0.13067\dotsc$ and it is attained at $q=305741$, as expected; the maximal value is $1.62693\dotsc$ and it is attained at $q=19$. The minimal value of $\mathfrak{G}^{+}_{q}/\log q$, $3\leq q\leq 10^{6}$, $q$ prime, is $0.451468\dotsc$ and it is attained at $q=918787$; the maximal value is $1.42626\dotsc$ and it is attained at $q=2053$. The points $(q,\mathfrak{G}_{q}/\log q)$ and $(q,\mathfrak{G}^{+}_{q}/\log q)$ in Figures 1-2 are colored in orange if $v(q)\leq 0.25$ ($65.65$% of the cases), in green if $0.25<v(q)\leq 0.5$ ($23.62$%), in blue if $0.5<v(q)\leq 0.75$ ($6.29$%), in black if $0.75<v(q)\leq 1$ ($4.21$%), and in red if $v(q)>1$ ($0.23$%). The behaviour of $\mathfrak{G}_{q}$ is the expected one since the red strip essentially corresponds with its minimal values, while the minima of $\mathfrak{G}^{+}_{q}$ seem to be less related to $v(q)$; we plan to investigate this phenomenon in the next future. The complete list of numerical results for $\mathfrak{G}_{q}$ and $\mathfrak{G}^{+}_{q}$ can be downloaded at the following web address: https://www.math.unipd.it/~languasc/EKcomput/results.

4.5.2. Computations for larger $q$.

For values of $q$ larger than $30$ millions the precomputation of $T$ and $S$, if performed on a single desktop computer, would require too much time; hence we parallelised them on the cluster previously mentioned. To check the correctness of such computations it is possible to use the following formulae; recalling that $\gamma=0.577215664901\dotsc$ and $\zeta^{\prime\prime}(0)=-2.006356455908\dotsc$, we have that

[TABLE]

Formula (29) is an immediate consequence of Theorem 2.5 of Deninger [6] and formula (30) follows from equation (7.10) of Dilcher [7].

Moreover, for being able to handle very large cases, we used a dedicated fftw interface888It is called the guru64 interface; see the user’s manual of fftw [10]. which is able to perform transforms whose length is greater than $2^{31}-1$.

In this way we were able to obtain an independent confirmation of Theorem 4 of [9] getting $\mathfrak{G}_{964477901}=-0.18237472563711916085\dotsc$, since we computed it using the quadruple precision. At the same time we also got $\mathfrak{G}^{+}_{964477901}=10.40222338242826353694\dotsc$ To do so we first split the computation, with a precision of 38 digits, of the needed decimated in frequency values of $S$ in $49$ subintervals $I_{j}$ of size $10^{7}$ each (for $T$ we would need $97$ intervals of such a length); the computation time required for each $I_{j}$ was on average about $1600$ minutes on one of the cluster’s machines. Then we passed such values to the programs that performed the FFT-step and got the final results. This last part needed about 23 minutes (long double precision) or 522 minutes (quadruple precision) of computation time on an Intel(R) Xeon(R) CPU E5-2650 v3 @ 2.30GHz, with 160 GB of RAM, running Ubuntu 16.04. A similar procedure let us to get analogous computation times for the long double precision evaluation of $\mathfrak{G}_{1217434451}=0.877596\dotsc$ and $\mathfrak{G}^{+}_{1217434451}=12.946690\dotsc$

We then looked for prime numbers $q$ such that $v(q)>v(964477901)=1.2369344\dotsc$ and we found that $v(2918643191)=1.2440460\dotsc$ In about 90 minutes of computation time for the FFT-step on the same machine mentioned before we got that $\mathfrak{G}_{2918643191}=0.302789\dotsc$ and $\mathfrak{G}^{+}_{2918643191}=12.573983\dotsc$, using the long double precision. In this case it seems that procedure in a) is much less stable than the one in b) probably because of the fact that $T(x)$ and $\psi(x)$ are much larger, for $x\to 0^{+}$, than $S(x)$ and $\log(\Gamma(x))$. Computations for further “good” candidates, in the sense that $v(q)>1.18$, like $q=193894451$, $212634221$, $251160191$ $538906601$, $1139803271$, $1217434451$, $1806830951$, $2488788101$, $2830676081$, $7079770931$ were also performed. The computations for these primes were performed on the cluster previously mentioned.

Moreover, for $q=9109334831$ we got that $\mathfrak{G}_{9109334831}=-0.248739\dotsc$, thus obtaining a new minimal value for $\mathfrak{G}_{q}$ and a new example of Theorem 4 of [9]; at the same time we also got $\mathfrak{G}^{+}_{9109334831}=12.128187\dotsc$ The precomputations for this case, performed with the same strategy used for the smaller primes $q$ mentioned in this paragraph, required about nine days on the cluster and the FFTs computation required about 1000 minutes on the Xeon machine mentioned before (this amount of time also depends on a runtime RAM swapping phenomenon) or 312 minutes on the new CAPRI infrastructure of the University of Padova (“Calcolo ad Alte Prestazioni per la Ricerca e l’Innovazione”; whose CPU is an Intel(R) Xeon(R) Gold 6130 CPU @ 2.10GHz, with 256 cores and equipped with 6TB of RAM). Such a result was then double-checked on CAPRI using the much slower algorithm a). A further new example of Theorem 4 of [9] we obtained is $\mathfrak{G}_{9854964401}=-0.096465\dotsc$ which required about ten days of time for the precomputations and 326 minutes for the FFT stage on CAPRI. As usual the result was double-checked using the approach a).

All the results mentioned in this paragraph are collected in Table 5. The PARI/GP scripts and the C programs used and the computational results obtained are available at the following web address: http://www.math.unipd.it/~languasc/EK-comput.html.

5. On the absolute value of the logarithmic derivative

of Dirichlet $L$-functions

Using (10)-(12), (14) and (18)-(19), for every odd prime $q$ we immediately get

[TABLE]

and

[TABLE]

Hence we can compute $M_{q}=\max_{\chi\neq\chi_{0}}|L^{\prime}/L(1,\chi)|=\max(M^{\textrm{odd}}_{q},M^{\textrm{even}}_{q})$ using the values of $\log\Gamma$ and $S$ obtained for the computation of $\mathfrak{G}_{q}$ and $\mathfrak{G}^{+}_{q}$. In Table 2 we give the values of $M_{q}$ for every odd prime up to $300$ computed, using PARI/GP, with a precision of $30$ digits. Using the data in Section 4.5.1 we also computed, on the Dell Optiplex machine previously mentioned, the values of $M_{q}$ and $M_{q}/\log\log q$ for every odd prime $q\leq 10^{6}$ and in Figures 3-4 we inserted their scatter plots that largely extend Figure 1 of Ihara-Murty-Shimura [17] (please remark that our $M_{q}$ is denoted as $Q_{m}$ there). Such data in Figures 3-4 also fit, for $q$ sufficiently large, with the estimate $M_{q}\leq(2+{o}\left(1\right))\log\log q$ as $q$ tends to infinity, proved, under the assumption of GRH, in Theorem 3 of [17]. We also remark that $M^{\textrm{odd}}_{q}>M^{\textrm{even}}_{q}$ for $62521$ cases over a total number of primes equal to $78497$ ($79.65\%$) and that $M^{\textrm{even}}_{q}>M^{\textrm{odd}}_{q}$ in the remaining $15976$ cases ($20.35\%$).

The complete list of numerical results for $M_{q}$ can be downloaded at the following web address: https://www.math.unipd.it/~languasc/EKcomput/results.

6. On the generalised Euler constants in arithmetic progressions $\gamma_{k}(a,q)$

Recall that $q$ is an odd prime. In the case we have to precompute $T(a/q)$ and we also need $\psi(a/q)$. Hence, as a by-product we can also obtain the values of the generalised Euler constants $\gamma_{0}(a,q)$ and $\gamma_{1}(a,q)$, see subsections 6.1-6.2. In practice this is done by activating an optional flag in the main gp script. The computation of $\gamma_{k}(a,q)$ for $k\geq 2$ is described in subsection 6.3.

6.1. Generalised Euler constants $\gamma_{0}(a,q)$

For $\gamma_{0}(a,q)$ with $1\leq a\leq q-1$, $q$ odd prime, by (4) we have

[TABLE]

Recalling that $\psi(1)=-\gamma$, we also have $\gamma_{0}(q,q)=(\gamma-\log q)/q.$

6.2. Generalised Euler constants $\gamma_{1}(a,q)$

For $\gamma_{1}(a,q)$ with $1\leq a\leq q-1$, $q$ odd prime, we can use (4) and (7). This way we get

[TABLE]

Moreover, since $\psi(1)=-\gamma$ and $T(1)=0$, we also have

[TABLE]

Using the formulae in the previous two paragraphs we computed $\gamma_{0}(a,q)$ and $\gamma_{1}(a,q)$ with $q$ prime, $3\leq q\leq 100$, $1\leq a\leq q$, in about 4 seconds of computation time with a precision of 30 digits.

Such results are listed towards the end of the gp-script file that can be downloaded here: http://www.math.unipd.it/~languasc/EK-comput.html.

6.3. The general case $\gamma_{k}(a,q)$, $k\geq 2$

The general case $\gamma_{k}(a,q)$, $k\in\mathbb{N}$, $k\geq 2$, $q\geq 1$, $1\leq a\leq q$, do not follow from the data already computed for the Euler-Kronecker constants since we need information about the values of $\psi_{n}(x)$, for every $2\leq n\leq k$. Such a direct computation of both $\psi_{n}(a/q)$ and $\gamma_{n}$ can be easily performed via eq. (4)-(5) using the PARI/GP summing function sumnum paying attention to submit a sufficiently fast convergent sum. For example, to compute $\gamma_{n}$, $n\in\mathbb{N}$, we used the formulae

[TABLE]

and

[TABLE]

which both easily follow from (6). We get, in less than 7 seconds of time and with a precision of at least $40$ digits, the results in Table 6; to be sure about the correctness of such results we computed them twice using the formulae (31)-(32) and then we compared the outcomes. These values are in agreement with the data on page 282 of Bohman-Fröberg [3] for $n=0,\dotsc,20$. For larger $n$’s the formulae in (31)-(32) seem to be not good enough to get precise results via the sumnum function with this precision level.

To compute $\psi_{n}(a/q)$ and, as a consequence, $\gamma_{k}(a,q)$, we can proceed in a similar way as we did for $T(a/q)$ and $\gamma_{1}(a,q)$, see the program Gen-Euler-constants.gp here http://www.math.unipd.it/~languasc/EK-comput.html. Towards the end of this program file you can find a large list (too long to be included here) of computed values of $\gamma_{k}(a,q)$ for $1\leq k\leq 20$, $1\leq q\leq 9$, $1\leq a\leq q$, with a precision of $20$ digits. In about 50 seconds of computation time we replicated Dilcher’s computations, since the values we got are in agreement with the data on pages S21-S24 of [7].

7. Tables