This is the first official “research” thread of the Polymath15 project to upper bound the de Bruijn-Newman constant . Discussion of the project of a non-research nature can continue for now in the existing proposal thread. Progress will be summarised at this Polymath wiki page.

The proposal naturally splits into at least three separate (but loosely related) topics:

- Numerical computation of the entire functions , with the ultimate aim of establishing zero-free regions of the form for various .
- Improved understanding of the dynamics of the zeroes of .
- Establishing the zero-free nature of when and is sufficiently large depending on and .

Below the fold, I will present each of these topics in turn, to initiate further discussion in each of them. (I thought about splitting this post into three to have three separate discussions, but given the current volume of comments, I think we should be able to manage for now having all the comments in a single post. If this changes then of course we can split up some of the discussion later.)

To begin with, let me present some formulae for computing (inspired by similar computations in the Ki-Kim-Lee paper) which may be useful. The initial definition of is

where

is a variant of the Jacobi theta function. We observe that in fact extends analytically to the strip

as has positive real part on this strip. One can use the Poisson summation formula to verify that is even, (see this previous post for details). This lets us obtain a number of other formulae for . Most obviously, one can unfold the integral to obtain

In my previous paper with Brad, we used this representation, combined with Fubini’s theorem to swap the sum and integral, to obtain a useful series representation for in the case because expressions such as diverge as approaches . Nevertheless we can still perform the following contour integration manipulation. Let be fixed. The function decays super-exponentially fast (much faster than , in particular) as with ; as is even, we also have this decay as with (this is despite each of the summands in having much slower decay in this direction – there is considerable cancellation!). Hence by the Cauchy integral formula we have

Splitting the horizontal line from to at and using the even nature of , we thus have

Using the functional equation , we thus have the representation

where is the oscillatory integral

The formula (2) is valid for any . Naively one would think that it would be simplest to take ; however, when and is large (with bounded), it seems asymptotically better to take closer to , in particular something like seems to be a reasonably good choice. This is because the integrand in (3) becomes significantly less oscillatory and also much lower in amplitude; the term in (3) now generates a factor roughly comparable to (which, as we will see below, is the main term in the decay asymptotics for ), while the term still exhibits a reasonable amount of decay as . We will use the representation (2) in the asymptotic analysis of below, but it may also be a useful representation to use for numerical purposes.

** — 1. Numerical investigation of — **

Python and Matlab code to compute for medium-sized values of are now available in the Polymath15 github repository. It is expected that all the zeroes of these functions for are real (this is equivalent to the Riemann hypothesis). For , is zero precisely when is a non-trivial zero of the Riemann zeta function, so the first few zeroes of occur at approximately

As increases, we have observed that the zeroes drift slightly to the left. This is consistent with theory, for instance the number of zeroes of of real part between and is known to be asymptotically

so as increases we should expect a few more zeroes in this region. (I had incorrectly omitted the denominator in a previous version of (4).) It seems like a reasonable near-term aim to improve the numerics to the point where we can confirm this asymptotic.

A related theoretical result is that the gaps between zeroes should behave locally like an arithmetic progression for large , in the sense that

This would also be nice to confirm numerically.

Theory also gives that the functions decay roughly like as ; see later sections for more precise asymptotics. To see this decay numerically for large , it may be necessary to switch over to a representation such as (2) with close to , otherwise the effect of numerical roundoff error may become unmanageably large.

** — 2. Dynamics of zeroes — **

Let denote the zeroes of . For sake of discussion let us suppose that the zeroes are always simple (this is in fact predicted by theory assuming RH; see Corollary 2 of Csordas-Smith-Varga. It may in fact be possible to prove this claim unconditionally, but we may not need this claim for the present project). If we implicitly differentiate the equation

in time, we (formally, at least) obtain the equation

The functions obey the backwards heat equation

and thus we have

If we Taylor expand around the zero as

for some coefficients with non-zero (because we are assuming the zeroes to be simple) then we have after a brief calculation

and also

On the other hand, from Weierstrass factorisation we expect (formally at least) that

and thus we should have

Putting all this together, we should obtain the dynamics

This is not rigorous for a number of reasons, most notably that the sum here is not absolutely convergent, but these equations should hold in a principal value sense at least. In the regime this was established in Lemma 2.4 of Csordas-Smith-Varga; it may be possible to establish this in the entire region .

If we write , we obtain the dynamics

Informally, the real parts repel each other, while the imaginary parts attract each other. In particular, once a zero is real, it should stay real.

If a zero is not real, then it has a complex conjugate . Isolating the attraction that the imaginary part feels to its complex conjugate , we obtain

Suppose that is a zero with maximal imaginary part (it is a result of Ki, Kim, and Lee that there are only finitely many non-real zeroes for with ). Then all the summands in (7) are non-positive, hence we have the differential inequality

Hence if denotes the maximal imaginary part of a zero of , we (formally, at least) have

or equivalently

whenever is positive. Thus the zeroes will be forced into the real axis in finite time, and in fact we can establish the bound

from this reasoning (this is a result of de Bruijn).

One could presumably do better by a more careful analysis of the sum in (7). The only way it seems that the inequality could be close to sharp is if the offending non-real zero is somehow isolated far away from all other zeroes except for its complex conjugate . For large , this is presumably in contradiction with the asymptotics (4); for small , perhaps one could use numerics to exclude this possibility?

At time , we know that the first zeroes are real (a result of Gourdon). Thus any non-real zero will initially have a very large value of . It would be nice to somehow be able to say that these zeroes continue to have very large real part for positive values of as well, but unfortunately the velocity in (5) could be large and negative if the zero is just to the left of another zero that is repelling it towards the origin. Is there some way to stop this? One may have to work with “clusters” of zeroes and study their centre of mass (or some similar statistic) as this may behave in a better way than an individual position function . Perhaps some of the identities in this previous post could be adapted for this purpose?

** — 3. Asymptotics of — **

Standard calculations give the gaussian integral identity

whenever are complex numbers with , where we integrate along a horizontal contour such as and we use the standard branch of the complex logarithm. More generally, Laplace’s method suggests that one has the approximation

whenever is an oscillating phase that has a single stationary point (with ) and is a slowly varying amplitude, and the integral is along a contour that does not start or end too close to . (Here one uses the standard branch of the square root function.) There are several ways to make this approximation rigorous, such as the method of steepest descent, the saddle point method, or the method of stationary phase, but for this discussion let us work completely informally. One can apply this method to analyse the integrals (3). For the highest accuracy, one should use the phase ; this is the approach taken for instance in my paper with Brad (in the case), but has the drawback that the stationary point has to be defined using the Lambert -function. A more “quick and dirty” approach, which seems to give worse error bounds but still gives broadly correct qualitative conclusions, is only take the phase and treat the factor as an amplitude. In this case, the stationary point occurs at

(where we use the branch of the logarithm that makes lie in the strip (1)), with

and

This suggests that we have the approximation

(I think that in order for this approximation to be rigorous, one needs to take to be close to the imaginary part of , in particular close to .) Now, from (2) one has

Here we consider the regime where is large and is positive but fairly small. Let’s just look at the term

Taking and , we approximately have

and so (after some calculation, and dropping terms of size ) we have the somewhat crude approximation

where

and

In particular, we expect the magnitude to behave like

Similarly for the other three expressions that appear in (9). If and is large, this suggests that the terms in (9) dominate, and furthermore of the four terms, it is the third term which dominates the other three. Thus we expect to have

where

Dividing the phase by and using the argument principle, we now see where the asymptotic (4) is supposed to come from. Taking magnitudes we also expect

in particular should be non-zero for fixed and large .

Hopefully we will be able to make these asymptotics more precise, and also they can be confirmed by numerics (in particular there may be some sign errors or other numerical inaccuracies in my calculations above which numerics might be able to detect).

## 92 comments

Comments feed for this article

27 January, 2018 at 6:58 pm

Sujit NairTerry, do you have any particular in mind when you say numerically verify ?

I am trying to understand how the numerical verification is going to be used?

27 January, 2018 at 8:49 pm

Terence TaoNot yet; but I would like to get an idea of how large needs to be compared with for these asymptotics to become visible. This would give us some intuition as to how large the numerical zero free region would have to be for various choices of , though perhaps we would not be directly using these asymptotics to bound the size of the zero free region.

Also, I think using the numerics to confirm a theoretical prediction would help root out any sign errors or other bugs in either the numerics or theory.

27 January, 2018 at 8:18 pm

AnonymousIn the statement of the third topic, can also be allowed (since for each , has only finitely many zeros).

It seems that (6) implies that the term in (7) and (8) should be .

It is not clear where is formula (9).

27 January, 2018 at 8:44 pm

Terence TaoThe equation (9) is now visible.

As for (7) and (8): if is the complex conjugate of , then will equal .

27 January, 2018 at 8:30 pm

AnonymousIt seems that in (4), the remainder term should be , and in the next displayed formula should be

27 January, 2018 at 8:47 pm

Terence TaoGood point; the estimates for for a fixed are in many ways much better than the estimates for the (far more intensively studied) (and the zeroes for the latter certainly do not arrange themselves in an arithmetic progression asymptotically), but there is no contradiction between the two, because the error terms for the former blow up in the limit .

28 January, 2018 at 6:12 am

KMThere are some updates.

There is now a function which can take complex initial guesses as input Ht_complex_root_finder(z,t) and try to find the closest complex root (for small to medium z values). Tried near the first few zeroes and it always converges near the actual zero location.

For eg.

guess_z, t, estimated_root_z

28.2 + 0.1j, 0, 28.26945 + 8.79635e-09j

28.2 + 0.1j, 0.5, 27.98018 – 4.55597e-10j

39 + 0.4j, 0, 42.04407 + 3.96173e-15j

58 + 0.3j, 0.35, 60.51891 – 1.37051e-05j

There are convergence problems as the imaginary part in the initial guess is increased too much, but the imaginary parts of potential roots are expected to be on the smaller side, so we can start with smaller guesses.

The approx function in the last section of the article, for computing H_t for large values of z has also been put in, but the calculations have to be verified.

28 January, 2018 at 8:26 am

AnonymousIs there an error-bound on the imaginare part? that means, how can you relate the numerical error to the for the error-free region in with ?

28 January, 2018 at 8:47 am

AnonymousWell the numerical error is quite influenced by roundoff errors, so I think unless epsilon is much greater than the computed imaginary parts, it will be difficult to confidently relate the two.

28 January, 2018 at 8:55 am

KBut isn’t this the goal of such calculations? finding the zeros, identifying the imaginary part of them to get ? (and improving the precision of the calculation to decrease ) I might have misunderstood this, sorry. in this case, how else would one find zero-free regions?

28 January, 2018 at 9:28 am

AnonymousYeah we will just have to keep improving the precision while making the calculations as we keep decreasing epsilon. For a given T and epsilon, if we compute all the zeroes in the bounded region and they turn out to have x<T and y<<epsilon, then we can be confident that the x epsilon region is zero free.

28 January, 2018 at 9:44 am

Anonymous..The last line should be xepsilon should then be zero free..

28 January, 2018 at 10:00 am

KThanks, then I understood it correctly. I am just wondering about a quantitative way to get a lower-bound , or more concretely concerning your answer, how to quantify ? Sorry if this is obvious for the experts, i am just curious.

28 January, 2018 at 7:26 am

dhjpolymathHere is some input values for which zeros have been calculated using the version of utility, h_t_probable…, h_t_real… which was available after the first few commits on github. In what follows “zeros a-b” correspond to replacing “enumerate(first_1000_H0_zeros)” with “enumerate(first_1000_H0_zero[600:700])” and the step in t is .01.

t=.20-.42, zeros 300-400,

t=.01-.19 zeros 400-500,

t=.15-.27 zeros 300-400,

t=.01-.09 zeros 500-600,

t=.05-.38 zeros 100-200,

t=.05-.26 zeros 200-300,

t=.15-.27 zeros 400-500,

t=.10-.28 zeros 300-400,

t=.15-.38 zeros 200-300.

These computations overlap some in places, are ongoing, and are set to terminate after t reaches .50 ( currently ~5 hours in on multiple 2.5 Ghz processors).

28 January, 2018 at 12:00 pm

Anonymousenumerate(first_1000_H0_zero[600:700]) would be replaced by enumerate(first_1000_H0_zero[a:b]) from your description.

28 January, 2018 at 7:27 am

AnonymousThere is a right parenthesis missing in the formula just before (2).

[Corrected, thanks – T.]28 January, 2018 at 8:12 am

AnonymousAfter ” using the even nature”, has a missing closing bracket.

[Corrected, thanks – T.]28 January, 2018 at 9:07 am

Sujit NairThis comment is targetted to those who plan to work on the numerical part of this project

=============================================

Hey All, I checked out the most recent commit and it has broken some functions. For example, “python python/H0_real_compute.py” doesn’t work anymore. And I don’t even know if it is supposed to work or not :(

I have some suggestions to make it easy for

(a) coders to understand the goals and objective of the code

(b) coders to contribute code to this project.

(c) non-coders to set up the code on their machines and run the files with minimal handholding

I foresee more people joining this project, and I recommend the following to build and maintain momentum.

1. Maintain a coding guideline — style (PEP-8 for python), appropriate comments, and docstrings.

2. Make clear what the expectations are for sending a PR. Having a good code structure and style will go a long way.

3. Only merge to master after a PR has been done and approved by at least one peer.

4. Migrate the code to Python 3+. The only minor conflicts I see are in the print statements and it is quite easy to port to 3+. The industry is moving in that direction. Moreover, Numpy is only going to support 2.7 until the end of the year.

Since polymath15 is just taking off, I thought I should put it out there before it becomes too big and unwieldy. I am happy to contribute in any way I can.

Thoughts?

28 January, 2018 at 9:42 am

AnonymousThanks for the guidelines. 2) and 3) are especially important and its possible I made some changes without realizing it. Have also added you as a collaborator there so you can make the necessary changes. Do let me know how I can help in cleaning up things. – KM

28 January, 2018 at 9:47 am

Terence TaoOne possible way to limit the number of complex zeroes of in the range is to get a lower bound on the number of real zeroes by numerically observing sign changes in , and then also obtain an upper bound on the total number of zeroes by establishing an effective version of the formula (4). If we believe in the Riemann hypothesis, these two numbers should always only differ by for any fixed . So even without believing in the RH, this sort of computation should show that there are only complex zeroes at most of real part at most , whenever is large enough that we can numerically find all the real zeroes through sign changes. If we are _really_ lucky, we might be able to get the bound of complex zeroes below 2 in some cases, in which case there are no complex zeroes whatsoever in this region since they have to come in pairs.

One could imagine estimating by the argument principle on a rectangle with vertices for some moderately large . For large enough, the contribution of the two long horizontal sides should be manageable by theory, and one just has to numerically integrate in the short vertical line segment from to , which might be doable even for quite large values of . Presumably tricks like this are already used in the literature on numerical verification of the RH; we might want to take a look at that literature for inspiration.

28 January, 2018 at 10:12 am

Terence TaoA variant strategy: thanks to numerical verification of RH, we know that there are no complex zeroes of in the rectangle for some fairly large (of the order of ). If one can show that there are no zeroes of for and , then I think that this implies that there are no complex zeroes of in the region , because any complex zeroes of have no way to enter this region as time advances without becoming real, at which point they stay real indefinitely (note that all the zeroes of are known to have imaginary part at most ). The point here is that the region that we have to keep zero-free has size rather than and so might be numerically feasible to work with even for quite large . The one issue though is that we may need a fairly fine mesh of 's, because there could be a complex zero that crosses this boundary at very high speed that might be somehow missed using a coarse mesh of 's. Actually I'm not sure how to do a rigorous error analysis at present (but presumably theoretical estimates such as (4) will play a role).

28 January, 2018 at 10:05 am

AulaI would prefer “these equations should hold” instead of “it should be that these equations hold” which seems too clumsy to me.

[Edited – T.]30 January, 2018 at 10:33 am

AulaI was quoting the second line below the first display above (5) which looks unchanged.

[Hopefully the edit is more permanent now – T.]28 January, 2018 at 10:10 am

AnonymousFor , the assumed dynamics for the simple zeros of can be proved as follows:

Since is entire in both and and also satisfies the backwards heat equation, it follows (as shown above) that for each real , any simple zero satisfies

(*)

(dropping the “time” index for simplicity)

The main idea of the proof is to use the fact that the RHS involve only derivatives with respect to (thanks to the heat equation) and some analytic properties of for each to construct a sequence of polynomials converging locally uniformly on to with the following properties:

(i) Each zero of each is also a zero with the same multiplicity

of both and .

(ii) Each zero of is a zero of some

Note that for each simple(!) zero of , the RHS of (*) is the limit of its values for the corresponding sequence of polynomials (converging locally uniformly to ). And the desired result follows (from the above ODE established above for polynomial zeros) since the LHS of (*) describes the dynamics for the same simple zeros (shared by both and its approximating polynomial sequence ).

To construct the polynomials, we use the fact ( KKL paper) that for each , is an entire even function of order 1, so by Hadamard’s factorization theorem, it has the form (after pairing its zeros and , which are distinct from zero, using any order with possibly finitely many repetitions – according to their multiplicities)

(**)

Where are the zeros of with positive real parts (and also, possibly, purely imaginary with positive imaginary parts), such that they are distinct from zero, ordered in any desired order (because (**) is absolutely locally uniformly convergent!) with possible repetitions (according to their multiplicities).

It remains only to define the polynomials

(***)

It follows that these polynomials have all the required properties above.

28 January, 2018 at 3:17 pm

AnonymousThe proof I gave above has the issue that the approximating polynomial sequence not necessarily satisfy the heat equation (as does), but fortunately, I found a much simpler (and direct) proof: one needs only(!) to repeat your derivation above, using that the product representation (**) is converging locally uniformly to , thereby allowing termwise logarithmic derivative in (**) with the following result for each simple(!) zero of :

where is the multiplicity of at , and the convergent(!) summation is over the zeros of (including multiplicities) which are distinct from and having positive real parts or are purely imaginary with positive imaginary parts.

28 January, 2018 at 1:16 pm

NickI recently added sinh-sinh, exp-sinh, and trapezoidal quadrature into Boost.Math. They work in arbitrary precision using boost.multiprecision, but do not work in the complex plane, so these integrals will need to be split into real and complex parts to make them work in this case. I haven’t tried using it with interval arithmetic, but that might help establish more rigorous bounds.

28 January, 2018 at 3:46 pm

Terence TaoHere’s a small question that arose from trying to understand the dynamics of a simple zero under backwards heat flow. Suppose we replace by a monomial for some natural number. Then the backwards heat flow (or ) gives the family of polynomials

so for instance

for ,

for ,

for , and so forth.

Experimentally, it seems that in all of these cases, the zeroes of are purely imaginary for negative times and purely real for positive times. I wonder if this can be proven rigorously? It suggests that the only way zeroes can collide is if they approach each other vertically, and then repel each other horizontally.

28 January, 2018 at 4:23 pm

AnonymousOne may try to build a “Newtonian kinematic model” for the zeros distribution dynamics involving “acceleration” and “forces” applied on each zero by the others (assuming that the “mass measure” for each zero is its multiplicity) obviously we have the mass conservation law (i.e. multiplicity conservation law for clusters of zeros), it seems likely that “total momentum” of a cluster of zeros is also conserved with “elastic collisions” among zeros.

28 January, 2018 at 10:00 pm

Terence TaoNewtonian mechanics certainly provides a lot of intuition, but it should be taken with a grain of salt, because Newton’s second law is a second-order equation, whereas the ODE here is first-order, with the repulsive “forces” directly affecting velocity rather than acceleration. Related to this is that these equations have a definite arrow of time: real zeroes repel each other going forward in time, but attract each other going backwards in time.

In my previous paper with Brad, we interpreted these equations as a gradient flow rather than a Newtonian system, but I am not sure if that is a useful perspective to take for the current project.

29 January, 2018 at 8:30 am

AnonymousLet be the “velocity” of a simple zero of . Its “acceleration” should be

So (by using the expressions for the velocities in terms of the zeros locations) the acceleration can be expressed as a (complicated) function of the distances among the zeros.

This “acceleration function” seem to decay according to a “inverse cube law” with the distances among the zeros, but this second order ODE does not seem to have any apparent time arrow.

29 January, 2018 at 9:01 am

Sujit NairYeah, I am not sure how useful a first order flow vs. a second order flow perspective is. Both, the gradient flow and a Newtonian flow (with friction) brings down the potential function over time. The rate at which the potential decreases are different for these flows. For example, in a Newtonian flow with small friction (relative to inertial mass), there will be lots of oscillations before it settles down.

But, what we want in this context is an understanding of the contours of the potential itself, independent of where it is plugged in (gradient flow or a Newtonian flow with friction).

The purely imaginary for negative times and purely imaginary for positive times also makes it delicate. What stops the following: purely imaginary for , purely real for and bouncing around in intermediate range?

Terry, you are missing a power of two in the derivative in the exponential sign above.

[Not sure what you are referring to here, but I am using as a synonym for . -T.]29 January, 2018 at 9:27 am

AnonymousIt may also be useful to note that the zeros dynamics seems to be invariant with respect to time translations, zeros locations translations and rotations. Does these symmetries indicate the possibility that this ODE zero dynamics has (still unknown) corresponding conserved invariants (analogous to “energy”, “linear momentum” and “angular momentum”) as in Newtonian mechanics

29 January, 2018 at 10:21 am

Terence TaoIt’s not invariant under rotations (except by 180 degrees), because the reciprocal function is not rotation invariant. It’s not a Hamiltonian system, so symmetries don’t correspond with conservation laws; instead, it is a gradient flow, so I think symmetries have a chance instead of corresponding with monotonicity formulae. In the case of finitely many real zeroes, in this previous blog post there were some interesting monotonicity formulae arising from the Hamiltonian , the energy , the centre of mass , and virial interaction . Don’t know if any of these continue to be useful in our situation though, with infinitely many real zeroes and a finite number of additional complex zeroes.

29 January, 2018 at 9:53 am

Terence TaoI think I am able to verify my previous claim that all the zeroes of the polynomials are real for and imaginary for . By replacing with (and with ) it suffices to prove the first claim. Using the analysis in Section 2 of the above post we know that the zeroes of are eventually real and simple. On the other hand, they are self-similar: the zeroes of are times the zeroes of . So the zeroes must be all real and simple.

By Rouche's theorem, this says more generally that when is a family of entire functions evolving by heat flow, at at some time they have a repeated zero of some order at some location , then shortly after the zeroes will lie near the horizontal line through , and shortly before they will lie near the vertical line through . If we have the functional equation , and is real, this means in particular that shortly after the zeroes will be purely real, and shortly before they will not be real (except perhaps for one zero if is odd).

So this seems to tell us how collisions of zeroes work for our functions : it is possible for one or more pairs of complex zeroes fall into the real line (possibly hitting an existing real zero when doing so) and create a real zero of order two or more, which then immediately splits into real simple zeroes. Zeroes that are real stay real for all time, and can't collide with each other (but can collide with an incoming pair of complex zeroes as mentioned above). Finally it is in principle possible for two or more complex zeroes to collide with each other if they are approaching in an asymptotically vertical direction, in which case they will split away from each other in a horizontal direction. As there are only finitely many complex zeroes after any positive time, I think a compactness argument shows that there are only finitely many collisions after any positive time, though I don't see any way to rule them out completely.

29 January, 2018 at 3:37 pm

AnonymousI think that in case of a multiple zero of , the dynamics of the zeros can be more complex, depending on the Puiseux expansion of each zero trajectory as a power series in some fractional power of . These Poiseux expansions of the trajectories of all the zeros (which can be determined by the Newton polygon method) can start with various directions (for instance, if the trajectory of some zero has a power expansion in in a neighborhood of , then there are two corresponding (or “conjugate”) trajectories of two other zeros such that all the three “conjugating” zero trajectories are “exploding” from the location of the multiple zero at with initial directions apart (so they can’t initially start at directions parallel to the horizontal or vertical axes.)

29 January, 2018 at 4:49 pm

Terence TaoActually, it looks like the only Puiseux series dynamics that is consistent with the equation is if the zeroes move away from the fixed point at a rate of . More precisely, if has a zero of order at , then for near , it appears that the zeroes of near take the form as , where are the (real, simple) roots of the degree polynomial . Indeed, one can normalise and , in which case is equal to , and the claim follows by rescaling by and using Rouche’s theorem.

In response to a previous blog post, there was a tweet displaying some pictures of dynamics of zeroes under heat flow. I wonder if that sort of code could be placed in the repository, it might be useful to play with to get some intuition for this dynamics.

29 January, 2018 at 6:45 pm

AnonymousGabriel’s code is now available in the matlab folder of the repo.

29 January, 2018 at 7:45 pm

Terence TaoThanks! I see there’s a lot of new python code in the repo, I wonder if someone can summarise what’s been going on there in the last few days.

29 January, 2018 at 8:44 pm

Sujit NairTerry, regarding the code, here are some updates.

1. There was some amount of refactoring along with the addition of code of conduct and how to contribute.

2. The outputs are moved out of the main code base. The plots and numbers are in output/plots and output/numbers folders respectively.

3. There is some discussion happening on the GitHub /issues page. I think some of it should be moved to this blog, or maybe a separate page dedicated to just the numerical part.

4. We have some additional data for the first 1000 roots of for various values of .

30 January, 2018 at 2:32 am

AnonymousAlso worth mentioning that only the data for the first 100 or so zeroes for a given t should be treated as accurate. Ht computation is currently relying on cos(zu) computation (theta=0), which is not very reliable for numerics after a point.

We are attempting to code the I_t_theta and K_t_theta functions for better estimates.

30 January, 2018 at 5:03 am

AnonymousIt means that if there is a zero of with multiplicity , there is a neighborhood of such that for any sufficiently close to , all zeros of in this neighborhood are simple and can be either

(i) analytic trajectory passing through , or

(ii) singular trajectory which is analytic in (i.e. has a converging power series in ) and consists of an “incoming singular trajectory” for .

If there are analytic trajectories and singular trajectories, then (clearly) .

It is interesting to observe that must be even! (this is because each singular power series in has also a distinct(!) “conjugate” series induces by changing the argument in the series to – thereby giving four(!) singular semi-branches – two incoming and two outgoing).

It seems that any partition of the multiplicity into analytic trajectories and singular trajectories (with corresponding incoming branches and outgoing branches) is allowed provided that m_2 is even.

This seems to imply that the velocity of any (incoming or outgoing) zero near any multiple zero of is bounded by .

29 January, 2018 at 3:46 pm

Fabien FriedliI think another way to prove this fact is to use Sturm theorem which relates the number of real roots of a polynomial to the number of sign changes in the corresponding Sturm sequence. If I’m not mistaken, in this particular case the n-th polynomial in the Sturm sequence is given by

where if is odd and if is even.

(I hope the latex code appears correctly, sorry if that is not the case.)

In particular, if , the main coefficient is always positive, so there is no sign changes in the sequence evaluated at and there is sign changes in the sequence evaluated at . Since the polynomial is of degree , we deduce that all the roots are real and simple.

This proof (if correct) is a bit longer and not as nice but it doesn’t use the analysis of Section 2. Also maybe it’s possible to have more information about the location of the zeros by using Sturm theorem on some finite interval instead of , since we have an explicit formula.

28 January, 2018 at 9:30 pm

Sujit NairThere is a typo in the third line after “In my previous paper with Brad,”.

You meant “for in the case”, right?

[Not quite; I’ve fixed the text for now, but unfortunately due to the HTML parser eating text between < and > signs, the problem may recur with my next edit of the post if I’m not careful. -T]29 January, 2018 at 1:58 am

Simo RyuAlways appreciate your good explanation on such interesting topic sir.

30 January, 2018 at 8:44 am

Terence TaoHere is an example of how some further analysis of the dynamics of zeroes can make it easier to numerically establish a bound of the form for some .

Let be a zero with positive imaginary part , then by equation (7) in the post we have the dynamics

where is the conjugate zero to . The first term is of course negative, reflecting the downward pull of the conjugate zero. Any real zero (with ), or more generally any zero lower than or equal to the current zero (so also contributes a non-positive term. A zero lying above the current zero (so ) contributes a positive term of , but its conjugate zero contributes a negative term of . So the net contribution of this zero and its conjugate is only positive if

which after some tedious calculation simplifies to

(This probably relates to some classical fact about the geometry of hyperbolae, but I don’t recognise it immediately.)

In particular, if there are no zeroes in the hyperbola above , then the net velocity is negative.

This has the following consequence:

PropositionSuppose are such that1. All the zeroes of in the strip are real;

2. There are no zeroes of in the hyperbolic wedge above for any . (Actually one only needs to check this for , as we know that there are no zeroes above this height.)

Then there are no zeroes of in the region for any .

ProofSuppose not, then let be the first time where there is a zero entering the forbidden region. By hypothesis 1, is positive. But by hypothesis 2, is negative, and so the zero in fact was in the forbidden region previous to time , a contradiction. (To be fully rigorous, one may also have to treat separately the case of a repeated zero.)This proposition is similar to a statement I made in a previous comment, in which and the wedge was replaced by a ray . But this new formulation has the advantage that the wedge in which one has to prevent zeroes is a positive distance away from the real line, which should make it easier to numerically exclude zeroes from this region, since we do expect there to be zeros traveling across the real axis.

30 January, 2018 at 10:36 am

AnonymousDoes condition 2 above should apply for all ?

30 January, 2018 at 12:11 pm

Terence TaoI had omitted a condition , although it is not terribly important. The region that one needs to keep zero free sits above and has diameter , so it is hopefully numerically feasible to actually verify this condition even for quite large values of .

30 January, 2018 at 1:37 pm

AnonymousIn particular, the simple condition should be sufficient (since for , each point of the hyperbolic wedge has imaginary part greater than 1). So it seems that what is required now is an effective(!) upper bound on the horizontal velocities of zeros with real parts .

30 January, 2018 at 1:15 pm

CraigSo, I’ve been playing around with some deformations of various formulas; I think I may have found one of interest.

Define . We have . We can also propagate this deformation back through the definition of to get deformed zeta functions: . (I think — may have messed up a bit, please check me). One benefit of this deformation is that you can now (hopefully) get a useful series expression for even when . Unfortunately, I’m not sure how useful the difference-differential equation is.

30 January, 2018 at 1:37 pm

CraigIf you additionally turn the 4 in into a deformation parameter, you lose the deformed zeta functions, keep the difference-differential equation, and gain a 3rd-order PDE for in and the new parameter. That last one may be of interest, giving a new flow for the zeros.

30 January, 2018 at 5:17 pm

AnonymousI think that in considering the zeros dynamics, one may assume without loss of generality that all the zeros of are simple, because the set (say) of “singular times” for which has a multiple zero is countable (each such is an intersection time of two trajectories from a countably many trajectories of countably many zeros of ) and therefore has measure zero, and the integrability of the zeros velocities (thanks to the continuity of the zeros) shows that the “singular times” set can be ignored for the purpose of estimating the zeros dynamics.

30 January, 2018 at 6:02 pm

SHI was looking through the data of the first 1000 zeros of H_.2(z), and the asymptotic (4) above appears to be very tightly matched. However, the asymptotic spacing of zeros seems way off, i.e. numerical errors seem unable to confirm this asymptotic. Also, some zeros are very close to each other. Does this imply a double root occurs? For example, I see the entry:

0.2 977 2790.7960430143326

0.2 978 2790.7962358872796

Can anyone confirm my observations?

30 January, 2018 at 6:27 pm

AnonymousPerhaps the root finding algorithm somehow converged to the same root of (hence the step size of should be decreased).

31 January, 2018 at 6:19 am

AnonymousCurrently please treat only the first 100 or so zeroes for a given t as reliable. We are implementing less oscillatory formulas for calculating H_t for larger values of z which should make the estimates much better.

31 January, 2018 at 4:08 pm

meditationataeI used the PARI/gp software package to numerically estimate the 10th zero of H_{0.2} and I obtain z_10 ~= 99.4597665 . The file: out_01_1000_20_22.txt on the git repo gives the value: z_10 ~= 99.547664956, or twice 49.773832478, which is the same as the tenth zeta zero according to Odlyzko:

http://www.dtc.umn.edu/~odlyzko/zeta_tables/zeros1

It might be useful at some point to add an estimated precision to the computed zeroes, for example when comparing various calculation methods.

1 February, 2018 at 7:20 am

SHBy “implementing less oscillatory formulas” do you mean you are using equation (3) above numerically, instead of the original equation for H_t ?

3 February, 2018 at 4:26 am

KMYes, there are some functions in the repo Ittheta and Ktheta implemented (but not present in the main branch yet), where we try to do this. These still have to be sorted out completely.

On the other hand, now that we have a approx functional eqn, and if we also get expansions of the first few error terms soon, then it could be much faster and still accurate enough compared to the expensive quadrature based approaches.

31 January, 2018 at 2:46 am

Lior SilbermanVague idea: when numerically find zeroes of , a key feature is that while the error term in the counting formula (Riemann–von Mangoldt, generalized in (4) above) can somewhat large, it is oscillatory and has quite a bit of cancellation. This observation (due to Turing, improving on Littlewood) makes it easy to tell if a count of zeroes is complete.

Perhaps something similar works for the deformations , in that the error term of equation (4) above is small on average over .

31 January, 2018 at 9:47 am

Lior SilbermanApologies for the terrible editing — I should have been more careful late at night.

[Unfortunately this implementation of wordpress does not permit users to edit their comments, but you can send corrections either in email to me or in an additional comment and I will be able to insert them -T.]31 January, 2018 at 3:15 am

SMTypo: “{x} is sufficiently large depending on {x} and {\varepsilon}” Twice {x}.

[Corrected, thanks – T.]31 January, 2018 at 11:02 am

Terence TaoSome back of the envelope calculations as to the point at which the single term dominates the expression (9) (thus producing a zero-free region):

let’s say and . Using the predicted asymptotics in the post, the summands decay in like , so the terms should dominate as long as . (There are some problems with these asymptotics for extremely large , but the contribution of those terms should be very small from some other estimates that I have not yet worked out in detail.)

The terms involving are about as large as the terms involving , so these should be negligible when .

The term involving is about as large as that of , so that should be negligible when .

So it seems that we will need to be in a regime in which before we can use these asymptotics to guarantee a zero free region for for . Given that we’re not likely to get to be much larger than , this suggests that we can realistically hope to take to be something like , which would lead to an upper bound for which is something like . Somewhat discouragingly, further gains only improve logarithmically with (roughly speaking, it seems like we will be able to prove a bound of the shape where is the largest real part for which we can numerically evaluate ), so there may be a limit as to how dramatic an improvement one would get from serious computer power.

In the converse direction, suppose there was a serious violation of the Riemann hypothesis at some scale , so that for some and . From the Riemann-von Mangoldt formula one would expect this zero to move at velocity comparable to (from the attraction of all the real zeroes at distance O(1) from this zero), so that it should take time about for this zero to hit the real axis. So this bound of is unlikely to be improvable other than in the size of the constant.

31 January, 2018 at 1:55 pm

AnonymousIn line 13, it seems that ” should be corrected ( instead of ).

[Corrected, thanks -T.]31 January, 2018 at 2:19 pm

AnonymousIs it expected that the upper bound on is optimized for optimal about the same size ?

31 January, 2018 at 4:42 pm

Terence TaoIt seems from the back of the envelope calculation that and will both end up being comparable to , which suggests that in the inequality , the first term will dominate. To put it another way, lowering the parameter looks like it is going to be more important than lowering the parameter (I’m assuming here that the computational difficulty of evaluating for and depends primarily on and not so much on or ).

1 February, 2018 at 10:33 am

KMSuppose for a start we try to tackle the case of t=0.49. c*exp(1/t) could be quite reasonable unless the constant c is large. Within this range one could also run root finding algorithms in the complex plane to get an estimate of sigma_max (the imaginary values will likely come out to be << sqrt(1-2t)). Do you think this is a feasible strategy with the current numerical algos (which aren't behaving well for large T yet)?

1 February, 2018 at 3:09 am

AnonymousIs it possible to get an effective upper bounds on the (finite!) number of the non-real zeros of and for their real parts for any given ?

1 February, 2018 at 9:50 am

Terence TaoQuite possibly. What one can certainly do is give an effective upper bound for the real part of a zero with ; it seems the KKL asymptotics give a bound of shape , which when combined with the Riemann-von Mangoldt formula for gives an exponential type upper bound for the number of zeroes of real part larger than . Presumably one can do better than this by importing some of the technology used to prove zero density theorems for the zeta function.

Ki and Kim (see Theorem 2.2 of this paper) also have some argument that shows that if all but finitely many zeroes have real part less than at some time , then at some time shortly afterwards, all but finitely many zeroes will in fact be real. I have not looked at this argument in detail and don’t know how effective it can be made.

1 February, 2018 at 3:56 am

AnonymousIt is interesting to observe that (by symmetry) the dynamics of each simple non-real zero of is closely related (in fact a reflection with respect to the real and imaginary axes) to the dynamics of its three “reflections” (except when is purely imaginary, i.e. , but in this case any purely imaginary zero must have even multiplicity since is even). Therefore, non-real zeros of are “appearing in quartets”.

1 February, 2018 at 8:56 am

Terence TaoI’m recording here a heuristic relationship between the asymptotics in the last part of the blog post and the Riemann-Siegel formula, which is a key tool in computing the Riemann zeta function numerically at large heights. Roughly speaking, this formula asserts that

for high in the critical strip, where . After some manipulation using Stirling’s formula, this eventually gives

where (assuming I have not made any sign errors etc.)

where and

.

In contrast, the heuristics in the above blog post seem to lead to the more general approximation

where

and and are as before. (I’ll write these calculations more carefully on the wiki soon.) The -dependent factor in contains a term of the form , which begins to induce significant decay in the summation once . So whereas computing accurately at height requires about operations, the number of operations drops down to once is significantly larger than . Indeed, for much larger than , the first term in the summation for dominates, and (assuming also) the second term is significantly smaller, which is why has no zeroes in this regime.

It is probably going to be important to develop some rigorous analogue of the Riemann-Siegel formula (with explicit error term) for for all between 0 and some target in order to do numerics accurately for large values of , and also to create effective zero-free regions for . This is probably what I’ll be focusing my efforts on.

1 February, 2018 at 9:15 am

AnonymousPerhaps the identification of the main error term (i.e. the second dominant term in some asymptotic expansion) may help to find a good effective upper bound on the error term.

1 February, 2018 at 4:52 pm

Terence TaoOK, I’ve written some computations at http://michaelnielsen.org/polymath1/index.php?title=Asymptotics_of_H_t which suggest that one has a generalised Riemann-Siegel approximation of the form

when is large, where

It would be interesting to see if this approximation is already detectable in the numerics. If so, we could as a first pass pretend that this approximation is exact, and see if we are able to conclude a good zero free region for this approximation to . The error term should be reasonably small – the usual Riemann-Siegel approximate functional equation has a relative error of , so if we can get into the range of or so then we should be in fairly good shape. (Of course, we will then have to do the messy work of actually getting an explicit error term in the above approximation, but this looks doable, albeit tedious. Also, it may be possible to add more terms to the Riemann-Siegel expansion to improve the error beyond .)

1 February, 2018 at 9:18 pm

meditationataeI tried to calculate explicitly the first zero of the above approximation to H , for t=0.

For t=0, I find using PARI/gp :

H(29.967325) = 1.087455808883869198146106527380154275744 E-11

as an approximation to the first zero of H_0, using your

approximate formulae.

This uses the defined functions H and F, with t a parameter set to zero:

F(Z) = Pi^(-Z/2)*gamma((Z+4)/2)*sum(Y=1, floor(sqrt(imag(Z)/(2*Pi))), exp((t/16)*(log((Z+4)/(2*Pi*Y*Y))^2))*Y^(-Z))

H(X,Y)= ( F((1+I*X-Y)/2) + conj(F((1+I*X+Y)/2)) )/4

2 February, 2018 at 8:50 am

Terence TaoThat’s not too bad, considering that the actual zero is near 28 or so and the approximation is only supposed to get good for large values of X! Hopefully the agreement gets better as X increases.

2 February, 2018 at 12:16 pm

KMOn estimating the first 2000 zeroes for t=0, using corresponding actual 2*zeta zero as a starting guess, we do find the median error decreasing with increasing T (average error seems to decrease much more slowly)

n in z_n avg z_n median absolute error

1-400 771.185 0.131

401-800 1874.222 0.076

801-1200 2837.853 0.072

1201-1600 3741.913 0.061

1601-2000 4608.312 0.053

Also the function is much more well behaved near large T values and computes fairly rapidly, even with arbitrary precision libraries. We could now try integrating this into the Odlyzko–Schönhage framework.

4 February, 2018 at 9:29 am

David Bernier (@doubledeckerpot)Up to x=500.0, I find 106 sign changes in H_0(x) , at a fine mesh/grid. According to tables of zeta zeros, there are 108 zeros with imaginary part below 250.0 .

One possibility is that the approximation has two sign changes (two zeros) missing, due to having the wrong sign. Otherwise, one would need a finer grid to capture the 2 missing zeros.

4 February, 2018 at 10:14 am

David Bernier (@doubledeckerpot)According to my computations, of the first 108 zeros of , only two are missed by counting sign changes.

The approximation to has constant sign in the interval [440.0, 444.0]; but has zeros with imaginary parts 220.71491884, and 221.43070555 , according to tables of zeta zeros.

2 February, 2018 at 3:32 am

AnonymousA possible way to verify this approximation (as well as to derive an effective upper bound on the error term) is to use a good known asymptotic approximation for and to propagate (in the “time” parameter ) according to the backwards heat equation (satisfied by ) to be an asymptotic approximation for .

The propagation of the initial approximation error should give the approximation error . Moreover, since the heat kernel is positive, it may be used to derive an upper bound on the approximation error by propagating a known upper bound on the the initial approximation error .

2 February, 2018 at 8:42 am

Terence TaoThis is a promising idea! One speed bump to this is that the fundamental solution to the backwards heat equation is only positive definite if one propagates in an imaginary direction rather than a real one. Let’s illustrate this using the toy function , which is very roughly of the form that shows up in the actual approximation. To propagate to negative times is quite easy: in terms of heat kernels, one has , and one also has the explicit integral representation which is quite tractable (this is basically the representation I use in my paper with Brad). But for positive times one instead has

which is problematic because can get rather close to the poles of on the upper imaginary axis. One may be able to get around this by shifting the contour, but I don't immediately see a clean way to do it. Certainly there is a big problem with the convergence of the formal integral as when is positive. In the case of the actual functions , one can get around this by using the functional equations to work only on the side of the integral, but this leads to messier formulae. But a good model problem would be to find a useful description of the forward propagation of by the backwards heat equation. (Equivalently: how does the Gamma function behave under the forward heat equation?)

2 February, 2018 at 9:40 am

AnonymousIs it possible to propagate the Gammma function by (at least formally) propagating the integrand in its integral representation?

2 February, 2018 at 9:49 am

Terence TaoYes, but the formal integral is quite divergent, and I don’t see a quick way to renormalise it to be convergent again. If for instance one propagates

formally by the heat flow, one obtains

which has a rather nasty singularity as or . There may be a way to formally shift contours to get a more tractable integral, though I don’t immediately see how to do it: to make the term manageable one would like to move away from the real axis at an angle of at least 45 degrees, but on the other hand the term becomes extremely problematic as soon as the imaginary part of reaches .

ADDED LATER: actually, it looks like the term only becomes a problem on the right-half plane, not on the left. So one could potentially get a good representation by using for instance a contour consisting of the upper imaginary axis and the right real axis (assuming for sake of argument that ). I’ll have to think more carefully about this…

2 February, 2018 at 12:12 pm

AnonymousA similar idea is to find a (sufficiently simple) integral representation of a good initial approximation of and than propagate (initially at least formally, and later perhaps by analytic continuation) to by propagating the integrand in its representing integral (i.e. a kind of formal “propagation under the integral sign”).

2 February, 2018 at 9:04 am

Terence TaoThis approximation also helps explain some of the Ki-Kim-Lee results. IF one inspects what happens to near some large real number by making the ansatz , and throws away some lower order terms (anything of size of the main terms), and using the Stirling approximation when is large and is small, it looks like one eventually gets the approximation

for , where are the polar coordinates of the quantity

In particular, one expects the zeroes to be real and equally spaced, with a spacing of about (in -coordinates) or (in -coordinates).

1 February, 2018 at 11:40 am

Lior SilbermanThese formulae also go under the title “approximate functional equation” (though in modern versions there is nothing approximate about them). A good reference is Theorem 5.3 of Iwaniec–Kowalski’s book; I also have notes from a short course by Andy Booker on computing many values of L-functions through the AFE and Poisson sum, with a view toward locating zeroes.

1 February, 2018 at 10:24 pm

meditationataeThe errors in the Riemann-Siegel approximation could be taken (i.e. divided by) relative to { \exp ( – \pi x/8 ) } , although this factor seems a bit too small.

2 February, 2018 at 2:19 pm

meditationataeAccording to Odlyzko’s zeta zeros tables, zeta has a zero at rho ~= 1/2 + i*49998.500133391 . So H_0 has a zero at 2 Im(rho) ~= 99997.000 . Your Riemann-Siegel type approximation to H_0 has a zero at x = 99997.02, giving an absolute error of about 0.02 .

For very large T, I have no tables, and PARI/gp’s own zeta(.) takes too long to finish.

3 February, 2018 at 4:21 am

KMOdlyzko’s page also has a link with first 2 million zeroes. Will be trying that out.

http://www.dtc.umn.edu/~odlyzko/zeta_tables/zeros6.gz

2 February, 2018 at 8:55 pm

Polymath15, second thread: generalising the Riemann-Siegel approximate functional equation | What's new[…] thread of the Polymath15 project to upper bound the de Bruijn-Newman constant , continuing this previous thread. Discussion of the project of a non-research nature can continue for now in the existing proposal […]

2 February, 2018 at 8:58 pm

Terence TaoAs this thread is getting rather long, I am rolling over to a new thread: https://terrytao.wordpress.com/2018/02/02/polymath15-second-thread-generalising-the-riemann-siegel-approximate-functional-equation/

4 February, 2018 at 3:15 pm

David Bernier (@doubledeckerpot)Something is a bit off in my calculation using formulas (2) and (3), because it finds a zero of close to 30.185 …

Broken over several lines, the input to PARI/gp is:

z=30.185;

t=0.0;

thet=Pi/8 – 1/(4*z);

K=0.0;

for(Y=1, 40,

b=z-9*I;

beta=Pi*Y^2;

K=K+2*Pi*Pi*(Y^4)*sum(J=1,64, intnumgauss(X=J/16, (J+1)/16, exp( t*(X+I*thet)^2 – beta*exp(4*(X+I*thet)) + I*b*(X+I*thet)), ) ) ;

b = z-5*I; K=K – 3*Pi*(Y^2)*sum(J=1,64, intnumgauss(X=J/16, (J+1)/16, exp( t*(X+I*thet)^2 – beta*exp(4*(X+I*thet)) + I*b*(X+I*thet)), ) ) );

print( real((K+conj(K)))/2)

output:

2.034481 E-8

4 February, 2018 at 3:52 pm

David Bernier (@doubledeckerpot)Changing the range of J to: 0, 64 improves things a lot. The formula evaluates to ~= 0 for 28.2694… (using Gaussian quadrature) :

z=2*z1;

t=0.0;

thet=Pi/8 – 1/(4*z);

K=0.0;

for(Y=1, 40, b=z-9*I;

beta=Pi*Y^2;

K=K+2*Pi*Pi*(Y^4)*sum(J=0,64, intnumgauss(X=J/16, (J+1)/16, exp( t*(X+I*thet)^2 – beta*exp(4*(X+I*thet)) + I*b*(X+I*thet)), ) ) ;

b = z-5*I;

K=K – 3*Pi*(Y^2)*sum(J=0,64, intnumgauss(X=J/16, (J+1)/16, exp( t*(X+I*thet)^2 – beta*exp(4*(X+I*thet)) + I*b*(X+I*thet)), ) ) );

print( real((K+conj(K)))/2)

output using real-precision = 70 digits:

-1.767 E-15

I can get closer to an output of zero by, e.g., using real-precision = 150 digits.

12 February, 2018 at 11:32 am

AnonymousIt seems that because of the horizontal repulsion among zeros and clusters of zeros of and because of the increasingly regular zero gaps as increases, it is expected that for sufficiently large , the horizontal zero velocities should become (asymptotically for large ) more uniform with a (slow) trend of drifting away(!) from the origin. Is it possible to show such asymptotic trend of the horizontal velocities?

If true, it may be exploited to show that a certain zero-free domain for some should remain so for any (i.e. by showing that zeros are not entering to this domain.)

14 February, 2018 at 3:13 am

Vassilis PapanicolaouDear Professor Tao,

thank you for presenting all these stuff to us in such a nice way.

This is my first comment to “Polymath.”

I would like to share few rough ideas regarding the dynamics of the zeros of :

It seems clear that is entire in both and . Hence the zeros are the branches of a global analytic function, say , defined on an infinite-sheeted Riemann surface . The projections of the branch points of on the complex plane are the values of for which we will have multiple zeros of , and these are countably many points, say , (possibly without accumulation points in the complex plane). Then, by analytic continuation one expects that the equations giving the dynamics of the zeros will continue to hold (in some sense?) for any complex , except at .