# Theta functions and optimal lattices for a grid cells model

Laurent Bétermin

Faculty of Mathematics, University of Vienna,
Oskar-Morgenstern-Platz 1, 1090 Vienna, Austria
<EMAIL_ADDRESS>ORCID id: 0000-0003-4070-3344

(date)

###### Abstract

Certain types of neurons, called “grid cells”, have been shown to fire on a
triangular grid when an animal is navigating on a two-dimensional environment,
whereas recent studies suggest that the face-centred-cubic (FCC) lattice is
the good candidate for the same phenomenon in three dimensions. The goal of
this paper is to give new evidences of these phenomena by considering a
infinite set of independent neurons (a module) with Poisson statistics and
periodic spread out Gaussian tuning curves. This question of the existence of
an optimal grid is transformed into a maximization problem among all possible
unit density lattices for a Fisher Information which measures the accuracy of
grid-cells representations in $\mathbb{R}^{d}$. This Fisher Information has
translated lattice theta functions as building blocks. We first derive
asymptotic and numerical results showing the (non-)maximality of the
triangular lattice with respect to the Gaussian parameter and the size of the
firing field. In a particular case where the size of the firing fields and the
lattice spacing match with experiments, we have numerically checked that it is
possible to find a value for the Gaussian parameter above which the triangular
lattice is always optimal. In the case of a radially symmetric distribution of
firing locations, we also characterize all the lattices that are critical
points for the Fisher Information at fixed scales belonging to an open
interval (we call these lattices “volume stationary”). It allows us to compare
the Fisher Information of a finite number of lattices in dimension 2 and 3 and
to give another evidences of the optimality of the triangular and FCC
lattices.

AMS Classification: Primary 49N20; Secondary 62P10.
Keywords: Grid cells, Lattices, Optimization, Fisher Information, Theta
functions.

###### Contents

  1. 1 Introduction and setting
    1. 1.1 Presentation of the problem
    2. 1.2 Some known results on optimal lattices
    3. 1.3 Numerical and asymptotic results
    4. 1.4 Characterization of volume-stationary lattices
    5. 1.5 Conclusion and open problems
  2. 2 Description of the model, derivation of formula (1.3) and scaling
  3. 3 Alternative formula and proof of Theorem 1.2
  4. 4 Volume-stationary lattices - Proof of Theorem 1.4
  5. 5 Numerical investigation
    1. 5.1 Parametrization of lattices in dimension 2 and 3
    2. 5.2 The pure discrete case: combinations of $\mathcal{Q}_{L}(y)$
    3. 5.3 Numerical results for $\mathcal{F}_{\mu}^{\alpha}$ in dimension 2
    4. 5.4 Numerical results for $\mathcal{F}_{\mu}^{\alpha}$ in dimension 3

##  1 Introduction and setting

###  1.1 Presentation of the problem

A challenging mathematical problem is to justify rigorously why periodic
patterns arise in nature and experiments: densest packing [52, 20], atoms in a
solid [31, 39], triangular lattices of Ginzburg-Landau vortices in type II
superconductors [1, 43], rich polymorphic behavior for systems with Coulombian
interactions [3, 33], etc. This type of problem, also called “Crystallization
problem” (see [17]), is often considered as variational, which means that
optimal periodic structures can be seen as extrema of certain functionals
usually derived from simplified models.

A typical example in Neurobiology is the existence of grid cells in the medial
entorhinal cortex (MEC) of the brain discovered by Hafting et al. [29] and
that have been brought to light for many mammals (rats, humans, etc.), see
e.g. [42]. Each of these neurons is tuned to the position of the animal and
fires when it crosses the sites of an approximate triangular lattice (also
called “hexagonal lattice”, see (1.7) for a precise formula and Figure 2 for a
representation) during a two-dimensional navigation. Furthermore, several
evidences have been found in the brains of flying bats and humans concerning
face-centred-cubic (FCC) lattices encoding the spatial representation in
three-dimensional displacements (see e.g. [54, 32]). It also appears that grid
cells can be split into modules where the firing periodic patterns share the
same scale and orientation but are shifted around an average position
belonging to a firing field (i.e. the corresponding spatial region that evokes
firing). Also, it is important to notice that this phenomenon occurs at
different scales: the lattice spacings of the grids variy in a discrete way
throughout the MEC as shown in [47]. In this paper, we aim to design a model
supporting the fact that one module is tuned to a triangular or a FCC pattern.

Several Mathematical attempts to justify the emergence of the triangular grid
for one module have been made (see e.g. [45, 51, 53, 35, 24, 46, 41, 2]), and
the present work is inspired by the one of Mathis et al. [35] where an
Information Theory point of view has been chosen. The goal is to maximize the
Fisher Information’s trace measuring the accuracy of grid cells representation
(also called “resolution”). The main novelty in our work compared to [35] is
that our model does not simplify into a best packing problem leading to the
triangular and FCC lattices as optimizers. Furthermore, we are considering
more spread out Gaussian tuning curves instead of short-range bump functions.
All these new assumptions transform the short-range packing problem derived in
[35] into an energy maximization problem pretty much the same way as the hard-
sphere crystallization question [30] and the best packing problems [50, 52,
20] were transformed into a minimization question among lattices for Gaussian
interactions in [36, 21].

We now briefly introduce our setting. The reader can refer to Section 2 for a
complete derivation and explanation of our main formula based on several
assumptions that are only partially recalled in this introduction. The main
statistical notions are also well explained in [23]. We assume that neurons
(grid cells) fire independently following a Poisson distribution. The average
firing number of a neuron is given by a lattice periodic tuning curve
$\Omega^{\alpha}_{L}$ which tunes the position $x\in\mathbb{R}^{d}$ of the
animal to the positive real number $\Omega^{\alpha}_{L}(x)\in\mathbb{R}_{+}$,
where

$\Omega_{L}^{\alpha}(x)=\theta_{L+x}(\alpha):=\sum_{p\in
L}e^{-\pi\alpha|p+x|^{2}},\quad\alpha>0,\quad\quad
L=\bigoplus_{i=1}^{d}\mathbb{Z}u_{i},$ (1.1)

the Euclidean norm being denoted by $|.|$ and $\\{u_{i}\\}_{1\leq i\leq d}$
being a basis of $\mathbb{R}^{d}$. Such infinite discrete set
$L\subset\mathbb{R}^{d}$ is called a $d$-dimensional lattice and
$\theta_{L+x}(\alpha)$ is called the translated lattice theta function (see
also [15]) which quantifies a Gaussian interaction between a point
$x\in\mathbb{R}^{d}$ and a lattice $L$ with Gaussian parameter $\alpha$.
Notice that the constant $\pi$ is only here for technical reasons and to fit
perfectly with the usual definition of theta functions. We call
$\mathcal{L}_{d}(V)$ the set of lattices with co-volume $V$, i.e.

$\mathcal{L}_{d}(V)=\left\\{L=\bigoplus_{i=1}^{d}\mathbb{Z}u_{i}:\\{u_{i}\\}_{i}\textnormal{
basis of $\mathbb{R}^{d}$},|\det(u_{1},...,u_{d})|=V\right\\}.$

Notice that the choice of a tuning curve of Gaussian type describes well
(based on the central limit theorem) the spatial positions at which neurons
fire, at the population level (see e.g. [24, 27]). Furthermore, since the
translated theta function is also traditionally seen as a building blocks for
a large class of lattice energies (see e.g. [19, 7, 16]), we have chosen to
focus on this important Gaussian case in this paper. However, our main formula
will be written in Section 2 in terms of a general function $f$ instead of the
Gaussian function $r\mapsto e^{-\pi\alpha r^{2}}$.

Since a module is given by a set of translated copies of the lattice
$\\{L+y_{j}\\}_{j}$ and tuning curves $\Omega^{\alpha}_{L+y_{j}}(x)$, we are
assuming that the empirical measure $\delta_{Y}$ associated with the set of
shifts $Y=\\{y_{j}\\}_{j}$ converges to a Radon measure $\mu$ (i.e. the firing
distribution) with support $\Sigma$ (called “firing field” in this context).
We will write $\mu\in\mathcal{M}(\Sigma)$. It is well-known that $\Sigma$
scales as the lattices in the experiments (see e.g. [18]), in the sense that
multiplying the lattice distances by $\lambda>0$ implies that we have to
replace the firing field by $\lambda\Sigma$. In order to take this scaling
effect into consideration, we define, for any $\lambda>0$, the rescaled
measure $\mu_{\lambda}$ on the Borelian set $\mathcal{B}_{d}$ of
$\mathbb{R}^{d}$ by

$\mu_{\lambda}(F)=\mu(\lambda F),\quad\forall F\in\mathcal{B}_{d}.$ (1.2)

Figure 1 depicts the situation in $\mathbb{R}^{2}$. Our goal is to decode the
position (e.g. $x=0$) of the animal given the number of spikes of a neuron
population. Therefore, as in [35], we are considering the Fisher information
$J_{M}(0)$ per neuron associated to the module
$M:=(L,\Omega_{L}^{\alpha},\mu,\Sigma)$ whose inverse of the trace is a lower
bound for the error made in the decoding process (this is the so-called
Cramer-Rao bound, see e.g. [23, Thm. 11.10.1]). In our setting, the trace of
$J_{M}(0)$ is shown (see Section 2) to be equal to (see also (2.14) for an
expanded formula)

$\mathcal{F}^{\alpha}_{\mu}(L):=\int_{\Sigma}\mathcal{Q}_{L}^{\alpha}(y)d\mu(y),\quad\textnormal{where}\quad\mathcal{Q}_{L}^{\alpha}(y):=\left|\nabla_{y}\sqrt{\theta_{L+y}(\alpha)}\right|^{2},$
(1.3)

and the question we are investigating is to find the lattice for which the
decoding error is the smallest possible., i.e. the lattice $L$ maximizing
$\mathcal{F}_{\mu}^{\alpha}$.

Figure 1: General configuration with arbitrary lattice $L$, firing field
$\Sigma=B_{R}$ (repeated periodically) and a generic point $y\in\Sigma$.
Figure 2: Triangular configuration with lattice $\mathsf{A}_{2}$ defined by
(1.7), firing field $\Sigma=B_{R}$ (repeated periodically) and a generic point
$y\in\Sigma$.

According to the following scaling formula (see Proposition 2.2) which hols
for any $\alpha>0$, $\lambda>0$, $\Sigma\subset\mathbb{R}^{d}$,
$\mu\in\mathcal{M}(\Sigma)$ and $L\in\mathcal{L}_{d}(1)$,

$\mathcal{F}_{\mu_{\lambda}}^{\alpha}(\lambda
L)=\lambda^{-2}\mathcal{F}_{\mu}^{\lambda^{2}\alpha}(L),$ (1.4)

it is therefore enough to consider unit density lattices in our problem stated
as follows.

Maximization Problem: Given $\alpha>0$ and $\mu$, what is the maximizer of
$\mathcal{F}^{\alpha}_{\mu}$ in $\mathcal{L}_{d}(1)$?

######  Remark 1.1 (Size of firing fields and lattice spacing).

It is implicit that $\Sigma$, i.e. the support of $\mu$, is also fixed in the
above problem. Furthermore, it has been observed that firing fields of a
single grid cell do not seem to overlap (see e.g. [29, Fig. 1]). In
particular, the ratio between the diameter of $\Sigma$ and the lattice spacing
of the grid appears to be approximately constant along the MEC – equal to
$\approx 0.3$ for Brun et al. [18] – and very recent experimental results by
Nagele et al. [37] also suggest that this ratio might be actually slightly
smaller (the firing fields are narrower). Our numerics in Section 5 (see
Figure 10) will support the fact that too much overlapping implies the non-
existence of a maximizer for $\mathcal{F}_{\mu}^{\alpha}$ in
$\mathcal{L}_{2}(1)$.

In this paper, we give several types of results concerning this maximization
problem: asymptotic, numerical and constrained. Indeed, it appears that
finding a rigorous proof of any optimality result for
$\mathcal{F}_{\mu}^{\alpha}$, where $\alpha,\mu$ are fixed, is out of reach
(see for instance Section 5.2). However, it is possible to study the behaviour
of the functional $\mathcal{F}_{\mu}^{\alpha}$ as $\alpha\to 0$, or
numerically when the parameters are fixed. Furthermore, staying in the class
of lattices that are critical points of $\mathcal{F}_{\mu}^{\alpha}$ in an
open interval of scales leads to comparing only a finite number of highly
symmetric lattices.

###  1.2 Some known results on optimal lattices

The type of optimality problem investigated in this paper has recently
attracted a lot of attention (see e.g. [17, Sect. 2.5]), especially in
Mathematical Physics and for energies per point of type

$E_{f}[L]:=\sum_{p\in
L}f(|p|^{2}),\quad\textnormal{where}\quad|f(r)|=O(r^{-\frac{d}{2}-\eta}),\eta>0\quad\textnormal{as
$r\to+\infty$},$ (1.5)

and where $f$ is interpreted as an interaction potential between particles.
The exponential case – i.e. where $f(r^{2})=e^{-\pi\alpha r^{2}}$, $\alpha>0$,
is a Gaussian function – appears to be of fundamental importance. In this
case, the energy defined in (1.5) is called the lattice theta function

$\theta_{L}(\alpha):=\sum_{p\in L}e^{-\pi\alpha|p|^{2}},$ (1.6)

which appears to be a building block for any energy of type $E_{f}$ where
$f=\mathcal{L}[\mu_{f}]$ is the Laplace transform of a measure $\mu_{f}$ on
$\mathbb{R}_{+}$ (see e.g. [19, 7, 16]). When $\mu_{f}$ is a nonnegative
measure, then $f$ is called completely monotone, which means that for all
$r>0$ and all $k\in\mathbb{N}\cup\\{0\\}$, $(-1)^{n}f^{(n)}(r)\geq 0$ (this is
a consequence of Hausdorff-Bernstein-Widder Theorem [6]). Therefore, the
optimality of a lattice in $\mathcal{L}_{d}(1)$ for
$L\mapsto\theta_{L}(\alpha)$ and for all $\alpha>0$ ensures the same result
for $E_{f}$ when $f$ is completely monotone. If we allow our structures $L$ to
be more general periodic configurations, then such optimality result is called
universal optimality and has been recently solved in dimension
$d\in\\{8,24\\}$ in [21] where the best (densest) packings of spheres (see
[52, 20]) $\mathsf{E}_{8}$ and $\Lambda_{24}$ are the unique minimizers in the
respective dimensions. However, this problem is still open in dimension $d=2$
where the triangular lattice

$\mathsf{A}_{2}:=\sqrt{\frac{2}{\sqrt{3}}}\left[\mathbb{Z}(1,0)\oplus\mathbb{Z}\left(\frac{1}{2},\frac{\sqrt{3}}{2}\right)\right],$
(1.7)

which is the two-dimensional best packing (see [50]), is conjectured to be the
unique minimizer. When restricted to the set of lattices, the lattice theta
function has been shown by Montgomery [36] to be minimized in
$\mathcal{L}_{2}(1)$ by $\mathsf{A}_{2}$ for any $\alpha>0$. Furthermore, such
result cannot be true in dimension $d=3$ as explained in [44] since the best
candidates, which are the Face-Centred-Cubic (FCC) and Body-Centred-Cubic
(BCC) lattices respectively defined by

$\displaystyle\mathsf{D}_{3}:=2^{-\frac{1}{3}}\left[\mathbb{Z}(1,0,1)\oplus\mathbb{Z}(0,1,1)\oplus\mathbb{Z}(1,1,0)\right]$
$\displaystyle\mathsf{D}_{3}^{*}:=2^{\frac{1}{3}}\left[\mathbb{Z}(1,0,0)\oplus\mathbb{Z}(0,1,0)\oplus\mathbb{Z}\left(\frac{1}{2},\frac{1}{2},\frac{1}{2}\right)\right]$

are dual of each other but not unimodular, i.e.
$\mathsf{D}_{3}^{*}\neq\mathsf{D}_{3}$. We recall that the dual lattice of $L$
is defined by $L^{*}=\\{y\in\mathbb{R}^{d}:y\cdot p\in\mathbb{Z},\forall p\in
L\\}$. Let us also recall the results of [15, 9] where $\mathsf{D}_{3}$ (resp.
$\mathsf{D}_{3}^{*}$) has been shown to be a strict local minimizer of
$L\mapsto\theta_{L}(\alpha)$ in $\mathcal{L}_{3}(1)$ for sufficiently large
(resp. small) $\alpha$ and a saddle point for sufficiently small (resp. large)
$\alpha$.

The extrema of the translated theta function
$(L,y)\mapsto\theta_{L+y}(\alpha)$ have also been studied. On the one hand,
all the critical points of $y\mapsto\theta_{L+y}(\alpha)$ as well as their
nature are known when $L$ is an orthorhombic (see [15, Sect. III.3]) or
triangular lattice (see [5]) and for all $\alpha>0$, see also Section 5.2. For
general lattices, only the trivial critical points are known (center of the
unit cell and midpoints). Therefore, all the zeros of
$\mathcal{Q}_{L}^{\alpha}$ are also known in these cases. On the other hand,
only few results are available concerning the extrema of
$L\mapsto\theta_{L+y}(\alpha)$ for fixed $y\neq 0$. In the case where
$y=c_{L}$ is the center of the unit cell of $L$,
$L\mapsto\theta_{L+c_{L}}(\alpha)$ does not have any minimizer (see [15, Prop.
1.3]) and the maximizer among two-dimensional lattices (resp. among
$d$-dimensional orthorhombic lattices) with the same density is the triangular
lattice $\mathsf{A}_{2}$ (resp. the simple cubic lattice $\mathbb{Z}^{d}$) as
we have shown in [12] (resp. in [15, Thm. 1.4] by mainly using tools from
[26]). Notice that other two-dimensional results of this type might be derived
from new discoveries concerning combinations of theta functions by Luo and Wei
in [34].

Finally, concerning the averaging of the shifts for the lattice theta
functions, i.e. for energies of the form
$\int_{\Sigma}\theta_{L+y}(\alpha)d\mu(y)$, recent advances have been made in
[14, 11] when $\mu$ is a radially symmetric measure sufficiently rescaled
around the origin or generated by a completely monotone kernel $\rho$ as
$d\mu(x)=\rho(|x|^{2})dx$. In particular, the above mentioned results [36, 15,
9, 21] for the lattice theta function $L\mapsto\theta_{L}(\alpha)$ on
$\mathcal{L}_{d}(1)$ still hold in these particular cases for
$L\mapsto\int_{\Sigma}\theta_{L+y}(\alpha)d\mu(y)$. The asymptotic results
presented in the next section are actually based on these facts combined with
the work of Regev and Stephens-Davidowitz [40] on the growth of (translated)
lattice theta functions.

###  1.3 Numerical and asymptotic results

Whereas $\mathcal{F}_{\mu}^{\alpha}$ defined by (1.3) has the translated
lattice theta function $\theta_{L+y}(\alpha)$ as a building block, it cannot a
priori be considered like a lattice energy of type $E_{f}$ defined by (1.5).
It is also clear that finding the global maximizer of
$\mathcal{F}_{\mu}^{\alpha}$ using an analytic or variational method appears
to be out of reach but numerically accessible. Moreover, some asymptotic
results can be shown.

It is indeed usual to study the asymptotic behaviour with respect to $\alpha$
of a problem involving lattice theta functions (see e.g. [15, Thm. 1.6]). When
$\alpha\to 0$, an alternative formula (see (3.1)) using a recent result by
Regev and Stephens-Davidowitz in [40] combined with new results on the soft
lattice theta function [14, 11] lead to the following asymptotic result.

######  Theorem 1.2 (The small $\alpha$ case - Non-maximality results).

Let $\mu\in\mathcal{M}(\Sigma)$ be radially symmetric with support
$\Sigma=B_{R}\subset\mathbb{R}^{2}$, for some $R>0$. Then:

  1. 1.

For all $L\in\mathcal{L}_{2}(1)$, $\exists\alpha_{0}>0$,
$\exists\lambda_{0}>0$, such that for all $\alpha\in(0,\alpha_{0})$, for all
$\lambda\in(0,\lambda_{0})$,

$\mathcal{F}^{\alpha}_{\mu_{\lambda}}(L)>\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\mathsf{A}_{2}).$

  2. 2.

For $\mu$ such that $d\mu(x)=\rho(|x|^{2})dx$, where $\rho$ is a completely
monotone function and all $L\in\mathcal{L}_{2}(1)$ there exists $\alpha_{0}$
such that for all $\alpha\in(0,\alpha_{0})$,
$\mathcal{F}^{\alpha}_{\mu}(L)>\mathcal{F}^{\alpha}_{\mu}(\mathsf{A}_{2})$.

  3. 3.

There exists $\alpha_{1}>0$ and $\lambda_{1}>0$ such that for all
$\alpha\in(0,\alpha_{1})$ and $\lambda\in(0,\lambda_{1})$,
$\mathcal{F}^{\alpha}_{\mu_{\lambda}}$ does not have any maximizer in
$\mathcal{L}_{2}(1)$.

In dimension $d=3$, the two first points hold for $L$ in an open ball of
$\mathcal{L}_{3}(1)$ centred at $\mathsf{D}_{3}^{*}$ and the third point also
holds.

This result agrees with the numerics we have performed for small $\alpha$ and
where the triangular lattice minimizes $\mathcal{F}_{\mu}^{\alpha}$ (see
Figure 9). Furthermore, our main numerical findings are the following (see
Section 5), choosing $\Sigma=B_{R}$, $\mu=\sigma_{R}$ being the uniform
measure on $B_{R}$ and $\alpha=\frac{10}{\pi}$:

  * •

Dimension 2. Maximality of the triangular lattice.

    * –

Global optimality of the triangular lattice. The triangular lattice
$\mathsf{A}_{2}$ is the unique maximizer of $\mathcal{F}_{\mu}^{\alpha}$ in
$\mathcal{L}_{2}(1)$ for
$R\in\left\\{0.1,0.2,0.3,0.4,0.5,\frac{1}{2}\sqrt{\frac{2}{\sqrt{3}}}\right\\}$.
The last value of $R$ corresponds to half of the side length of
$\mathsf{A}_{2}$.

    * –

Non-existence of a maximizer for large $R$. If $R\geq 0.59$, then
$\mathcal{F}_{\mu}^{\alpha}$ does not have a maximizer in $\mathcal{L}_{2}(1)$
and $\mathsf{A}_{2}$ is its unique minimizer in $\mathcal{L}_{2}(1)$.

    * –

Comparison with experimental values from [18]. If $R=0.16$, which corresponds
to the case where the ratio between the diameter of $\Sigma$ and the lattice
spacing of the triangular lattice is $0.3$ (as suggested by [18]), we observe
that $\mathsf{A}_{2}$ is the unique maximizer of
$\mathcal{F}_{\mu}^{\bar{\alpha}}$ in $\mathcal{L}_{2}(1)$ if
$\bar{\alpha}>\bar{\alpha}_{0}\approx 1.25$. For
$\bar{\alpha}\leq\bar{\alpha}_{0}$, $\mathcal{F}_{\mu}^{\bar{\alpha}}$ does
not seem to have a maximizer in $\mathcal{L}_{2}(1)$ and $\mathsf{A}_{2}$ is a
minimizer in $\mathcal{L}_{2}(1)$.

  * •

Dimension 3. Maximality of the FCC lattice

    * –

Comparison of cubic lattices. There exists $R_{0}\approx 0.57$ such that for
all $R\leq R_{0}$,

$\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3})>\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3}^{*})>\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{3}).$

    * –

Local maximality of the FCC lattice. $\mathsf{D}_{3}$ is a local maximizer of
$\mathcal{F}_{\mu}^{\alpha}$ in $\mathcal{L}_{3}(1)$ for
$R\in\\{0.1,0.2,0.3,0.4,0.5,2^{-5/6}\\}$. The last value of $R$ corresponds to
half of the side length of $\mathsf{D}_{3}$.

According to the scaling formula (1.4), the same qualitative behaviour holds
for any other value of $\alpha$ with rescaled quantities.

###  1.4 Characterization of volume-stationary lattices

In [47], it has been observed that rats have grid cells with grid spacing from
35.2cm to 171.7cm, from dorsal to ventral MEC. The set of possible scales has
been shown to be discrete for all the studied species, with different values.
If we believe that grid cells are universal among a lot of animals, in such a
way that grid scales are covering all the possible values of a certain (open)
interval, we can focus on lattices $L$ that are critical points (or maximizer)
of $\mathcal{F}_{\mu}$ in $\mathcal{L}_{d}(V)$ for all $V$ in an open
interval. We show a result which is similar to the one obtained in [10, Sect.
III] for energies of type $E_{f}$ defined by (1.5), using tools from [25, 28].

We first recall the notion of strongly eutactic layer exactly as we have done
it in [10, Sect. III].

######  Definition 1.1 (Strongly eutactic layer).

Let $L\in\mathcal{L}_{d}(1)$. We say that a layer $\mathfrak{m}=\\{p\in
L;|p|=\lambda\\}$, for some $\lambda>0$, of $L$ is strongly eutactic if
$\sharp\mathfrak{m}=2k$ for some $k\in\mathbb{N}$ and, for any
$x\in\mathbb{R}^{d}$,

$\sum_{p\in\mathfrak{m}}\frac{(p\cdot x)^{2}}{|p|^{2}}=\frac{2k}{d}|x|^{2}.$

######  Remark 1.3.

After a suitable renormalization, $\mathfrak{m}$ is also called a spherical
$2$-design (see e.g. [4]).

######  Theorem 1.4 (Volume-stationary lattices for
$\mathcal{F}^{\alpha}_{\mu}$).

We assume that $\Sigma\subset\mathbb{R}^{d}$ is compact and
$\mu\in\mathcal{M}(\Sigma)$ is radially symmetric. Then,
$L\subset\mathbb{R}^{d}$ is a critical point of
$L\mapsto\mathcal{F}_{\mu_{\lambda}}(\lambda L)$ in $\mathcal{L}_{d}(1)$ for
all $\lambda$ in an open interval $I$ if and only if all the layers of $L$ are
strongly eutactic and $I=(0,\infty)$. In particular,
$L\in\\{\mathsf{A}_{2},\mathbb{Z}^{2},\mathbb{Z}^{3},\mathsf{D}_{3},\mathsf{D}_{3}^{*}\\}$
in dimension $d\in\\{2,3\\}$.

Therefore, restricting our maximization problem to volume-stationary lattices
is very simple since one have only to compare the lattices belonging to
$\\{\mathsf{A}_{2},\mathbb{Z}^{2},\mathbb{Z}^{3},\mathsf{D}_{3},\mathsf{D}_{3}^{*}\\}$.
Furthermore, Theorem 1.4 tells us that
$\mathsf{A}_{2},\mathbb{Z}^{2},\mathbb{Z}^{3},\mathsf{D}_{3}$ and
$\mathsf{D_{3}}^{*}$ are critical points of $\mathcal{F}_{\mu}^{\alpha}$ for
any $\alpha>0$ and $\mu$ radially symmetric on a ball $\Sigma=B_{R}$ for any
$R>0$. In this case, combined with our numerical findings stated above (see
also Section 5), we obtain new evidences that $\mathsf{A}_{2}$ and
$\mathsf{D}_{3}$ are optimal for the grid cells problem if we restrict our
study to volume-stationary lattices and radially symmetric firing field and
measure $\mu$. More precisely, our numerics suggest that, for $\Sigma=B_{R}$
and $\mu=\sigma_{R}$ being the uniform measure on $B_{R}$, we have to
following.

  * •

For all fixed $R>0$, there exists $\alpha_{R}=\alpha(R,d)$ such that

    * –

in dimension $d=2$,
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{A}_{2})>\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{2})$
(resp.
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{A}_{2})<\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{2})$)
if $\alpha>\alpha_{R}$ (resp. if $\alpha<\alpha_{R}$);

    * –

in dimension $d=3$, for all $\alpha>\alpha_{R}$,
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3})>\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3}^{*})>\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{3})$.
Furthermore, there exists $\tilde{\alpha}_{R}=\tilde{\alpha}(R)$ such that for
all $\alpha<\alpha_{R}$,
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3})<\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3}^{*})<\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{3})$.

  * •

For all $\alpha>\alpha_{d}$ for some $\alpha_{d}>0$, there exists
$R_{\alpha}=R(\alpha,d)$ such that

    * –

in dimension $d=2$,
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{A}_{2})>\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{2})$
(resp.
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{A}_{2})<\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{2})$)
if $R<R_{\alpha}$ (resp. if $R>R_{\alpha}$);

    * –

in dimension $d=3$, for all $R<R_{\alpha}$,
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3})>\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3}^{*})>\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{3})$.
Furthermore, there exists $\tilde{R}_{\alpha}=\tilde{R}(\alpha)>0$ such that,
for all $R>\tilde{R}_{\alpha}$,
$\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3})<\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3}^{*})<\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{3})$.

###  1.5 Conclusion and open problems

We conclude that, whereas the rigorous study of $\mathcal{F}_{\mu}^{\alpha}$
stays a very challenging problem, our asymptotic and numerical investigations
show new evidences of the optimality of the triangular and face-centred-cubic
lattices for the grid cells problem in dimension 2 and 3. In particular, our
results suggest the following conjecture that also covers the experimental
case where $R=0.16$ (see Remark 1.1 on the size of the firing fields).

######  Conjecture 1.5 (Maximality of the best packing in dimension 2 and 3).

Let $\sigma_{R}$ be the uniform measure on $B_{R}$ and
$\alpha=\frac{10}{\pi}$. Furthermore, let
$R_{2}:=\frac{1}{2}\sqrt{\frac{2}{\sqrt{3}}}$ and $R_{3}:=2^{-\frac{5}{6}}$.
Then we have:

  1. 1.

For all $R\leq R_{2}$, $\mathsf{A}_{2}$ is the unique maximizer of
$\mathcal{F}_{\sigma_{R}}^{\alpha}$ in $\mathcal{L}_{2}(1)$.

  2. 2.

For all $R\leq R_{3}$, $\mathsf{D}_{3}$ is the unique maximizer of
$\mathcal{F}_{\sigma_{R}}^{\alpha}$ in $\mathcal{L}_{3}(1)$.

This problem by itself turns out to be a very interesting mathematical
question with the need to derive new results for translated lattice theta
functions. For instance, one may want to push further the properties of the
Heat Kernel and to use the connection with $\mathcal{F}_{\mu}^{\alpha}$
explained in Remark 2.1.

Moreover, we have also compared our model to the results found by Brun et al.
in [18] about the size of the firing field $\Sigma$ related to the lattice
spacing. It appears that the optimality properties of our model support the
experimental observations of Brun et al. for a wide range of parameters
$\alpha>1.25$. Also notice that in their recent work [37], Nagele et al.
suggest that the ratio between the diameter of $\Sigma$ (i.e. $2R$ if
$\Sigma=B_{R}$) and the lattice spacing should be slightly smaller than the
one found in [18]. Nevertheless, it seems that one could find again a range of
$\alpha$ such that the triangular lattice again maximizes the Fisher
Information’s trace.

For a future work, it would also be interesting to study the $\mu$-dependence
of the maximizer of $\mathcal{F}_{\mu}^{\alpha}$, in particular when $\mu$ is
not uniform nor radially symmetric. In the latter case, one could try to
explain the existence of sheared grids observed for instance in [48], also for
grid cells when for example the animal visits the boundary of its displacement
zone. Another direction or research would also be to replace the Gaussian
function by another tuning curve and to study the same type of optimality
problem. Finally, we think that our results and observation might be useful
for other biological or physical models where periodic Fisher Information with
Poisson distribution are involved.

Plan of the paper. In Section 2, we explain how we obtain formula (1.3) and on
which mathematical assumptions we have built our model. Furthermore, the
scaling formula (1.4) is shown in Proposition 2.2. Section 3 is devoted to the
proof of Theorem 1.2 and to the alternative formula 3.1. The characterization
of volume-stationary lattices is done in Section 4 where the proof of Theorem
1.4 is given. Finally, our numerical investigations are presented in Section 5
where we first discuss the parametrization of lattices as well as the
properties of $\mathcal{Q}_{L}^{\alpha}(y)$.

##  2 Description of the model, derivation of formula (1.3) and scaling

In this section, we are justifying formula (1.3) by listing all our
assumptions. As explained in the introduction, we are following Mathis et al.
[35] for building our model, with only small modifications.

We consider neurons (called grid cells) that are firing (spiking) according to
the position (the stimulus) $x\in\mathbb{R}^{d}$ of an animal. We are also
considering an arbitrary fixed time interval $[0,T]$ where the neurons are
spiking. We are assuming the following hypothesis (A1)-(A6) where the Gaussian
is replaced for the moment by an arbitrary function $f$ as in (1.5).

  * (A1)

Independence of neurons. The neurons we consider are independently firing.

Given $n$ independent neurons such that the $i^{th}$ neuron fires $k_{i}$
times in the time interval $[0,T]$, the probability of having
$K=(k_{1},...,k_{n})$ spikes when an animal is at position
$x\in\mathbb{R}^{d}$ is

$P(K|x)=\prod_{i=1}^{n}P_{i}(k_{i}|x),$ (2.1)

where $P_{i}(k_{i}|x)$ is the probability of firing the $i^{th}$ neuron
$k_{i}$ times when the animal is at $x$.

  * (A2)

Firing Poisson statistics. The neuron’s firing follows a Poisson distribution.
More precisely, we assume that

$P_{i}(k_{i}|x)=\frac{(\Omega_{i}(x))^{k_{i}}e^{-\Omega_{i}(x)}}{k_{i}!},$
(2.2)

where the tuning curve $\Omega_{i}(x)$ is the average firing number of the
$i^{th}$ neuron at position $x$ in the time interval $[0,T]$.

  * (A3)

Spread out lattice-periodic tuning curve shape. All the tuning curves of the
neurons are identically equal to

$\Omega_{i}(x)=\Omega_{L}(x)=E_{f}[L+x]:=\sum_{p\in
L}f(|p+x|^{2}),\quad\forall i\in\\{1,...,n\\},$ (2.3)

for a lattice $L\in\mathcal{L}_{d}(1)$ and a function
$f:[0,\infty)\to\mathbb{R}$ as in (1.5), i.e.

$|f(r)|=O(r^{-\frac{d}{2}-\eta}),\eta>0\quad\textnormal{as $r\to+\infty$},$

in order for $\Omega_{L}(x)$ to be absolutely convergent. We will call
$\Omega_{L}:\mathbb{R}^{d}\to\mathbb{R}_{+}$ the grid cell’s tuning curve.
When $f$ has its maximum at $r=0$, it means that there is more spikes when the
animal cross a lattice site, which is the usual assumptions made on $f$. This
tuning curve’s shape is slightly different than the periodic boundary
conditions assumed by Mathis et al. in [35].

We now define a grid module as an ensemble of $N\in\mathbb{N}$ shifted grid
cells $\\{\Omega_{L+y_{j}}\\}_{j=1}^{N}$ where the vectors
$Y=(y_{j})_{j}\subset\mathbb{R}^{d}$ are called the spatial phases of the
module. Furthermore, we have for all $x\in\mathbb{R}^{d}$, all lattice
$L\in\mathcal{L}_{d}(1)$ and all shifting vector $y_{j}\in Y$,

$\Omega_{L+y_{j}}(x)=E_{f}[L+y_{j}+x]:=\sum_{p\in L}f(|p+x+y_{j}|^{2}).$ (2.4)

Therefore, for any spatial phase with shifting vector $y_{j}$, we write the
probability in (2.1) as

$P_{j}(K|x)=e^{-n\Omega_{L+y_{j}}(x)}\prod_{i=1}^{n}\frac{(\Omega_{L+y_{j}}(x))^{k_{i}}}{k_{i}!},\quad
j\in\\{1,...,N\\}.$ (2.5)

As explained in [35], we are interested in solving the following mathematical
question: given a spike count vector $K=(k_{1},...,k_{n})$, where is the
animal? The estimation of this position is written $\widehat{x}(K)$. We
therefore want to minimize the square of the error in the decoding process
which is given by

$\varepsilon(x|\widehat{x}):=\mathbb{E}_{P(K|x)}\left(\|x-\widehat{x}(K)\|^{2}\right).$
(2.6)

Minimizing the error $\varepsilon(x|\widehat{x})$ is the same as minimizing
the trace of the covariance matrix $\textnormal{Cov}(x,\widehat{x})$.
Furthermore, the Fisher Information (see e.g. [23, Eq. (11.291)]) associated
to the probability $P$ with Poisson distribution is defined by the matrix
$J(x)=(J_{\ell,m}(x))_{\ell,m}$ where

$J_{\ell,m}(x):=\int\left[\partial_{x_{\ell}}\log P(K|x)\partial_{x_{m}}\log
P(K|x)\right]P(K|x)dK.$ (2.7)

  * (A4)

Unbiased estimator assumption. We assume that the Cramer-Rao lower bound holds
(see e.g. [23, Thm. 11.10.1]), i.e. $\textnormal{Cov}(x,\widehat{x})\geq
J(x)^{-1}$. This is actually automatically satisfied in the independent
Poisson process case.

For the grid module, (A1) ensures that the Fisher Information is the sum of
Fisher Information for all the $N$ spatially shifted phases $Y$, i.e.

$J_{\ell,m}(x)=\sum_{j=1}^{N}\int\left[\partial_{x_{\ell}}\log
P_{j}(K|x)\partial_{x_{m}}\log P_{j}(K|x)\right]P_{j}(K|x)dK.$ (2.8)

A straightforward computation gives

$\mathrm{Tr}J(x)=\sum_{\ell=1}^{d}\sum_{j=1}^{N}\left(\partial_{x_{\ell}}\Omega_{L+y_{j}}(x)\right)^{2}\int\left(-n+\frac{\sum_{i=1}^{n}k_{i}}{\displaystyle\Omega_{L+y_{j}}(x)}\right)^{2}P_{j}(K|x)dK.$
(2.9)

The integral is in fact the Fisher Information of $n$ independent Poisson
processes with the same parameter $\lambda=\Omega_{L+y_{j}}(x)$, which means
that it is $n$ times the Fisher Information of one single Poisson process with
the same parameter (see e.g. [23, p. 395]). Therefore, for all $y_{j}\in Y$,

$\int\left(-n+\frac{\sum_{i=1}^{n}k_{i}}{\displaystyle\Omega_{L+y_{j}}(x)}\right)^{2}P_{j}(K|x)dK=\frac{n}{\Omega_{L+y_{j}}(x)}.$
(2.10)

We therefore obtain that

$\mathrm{Tr}J(x)=n\sum_{j=1}^{N}\mathcal{Q}_{L}(y_{j}+x),\quad\textnormal{where}\quad\mathcal{Q}_{L}(y):=\frac{\left|\nabla_{x}E_{f}[L+y]\right|^{2}}{E_{f}[L+y]}.$
(2.11)

Optimizing such lattice energy is a huge challenge since the maximizer varies
a lot with $x$ and $Y$ (see Section 5.2). That is why it is not appropriate to
our study. Therefore, we are again following [35] in order to transform the
discrete set of shifts $Y$ into a continuous one written $\Sigma$.

  * (A5)

Aggregation of spatial phases. We assume that
$\displaystyle\frac{1}{N}\sum_{j=1}^{N}\delta_{y_{j}}\rightharpoonup\mu\in\mathcal{M}(\Sigma)$
as $N\to\infty$ in the weak-$\ast$ sense, with support in a compact set
$\Sigma\subset\mathbb{R}^{2}$. This is suggested by experiments (see [29]).
The set $\Sigma$ can be considered as the support of the grid cells module,
also called the firing field of the module. Actually, there is a scale
dependence between $\Sigma$ and $L$ in the sense that, if the distances in the
lattice are multiplied by a real number $\lambda$, then it is the same for
$\Sigma$ (see e.g. [38, 47] and Proposition 2.2). Furthermore, it has been
observed in [18] that the smallest width of $\Sigma$ is 30cm-50cm for a
lattice spacing (corresponding to the triangular lattice) of 50cm-1m (see also
[37] for other data with smaller firing fields width compared to lattice
spacing), i.e. firing fields are not overlapping (as in Figure 2). We also
notice that, presented as above, $\mu$ should be a probability measure.
However, in the following we are considering a more general Radon measure.

Therefore, we want to maximize the following Fisher Information per neuron for
the module $M=(L,\Omega_{L}^{\alpha},\mu,\Sigma)$, where the limit is
justified by the fact that $y\mapsto\mathcal{Q}_{L}(y)$ is continuous and
bounded on the compact set $\Sigma$,

$J_{M}(x):=\frac{1}{n}\lim_{N\to\infty}\frac{\mathrm{Tr}J(x)}{N}=\int_{\Sigma}\frac{\left|\nabla_{x}\Omega_{L+y}(x)\right|^{2}}{\Omega_{L+y}(x)}d\mu(y)=\int_{\Sigma}\frac{\left|\nabla_{x}E_{f}[L+y+x]\right|^{2}}{E_{f}[L+y+x]}d\mu(y).$
(2.12)

  * (A6)

Re-centering of the problem. We assume that $x=0$. Thus, the problem is
simplified to maximize the resolution at $x=0$ of the decoding process.

Therefore, given $f$ satisfying assumptions (1.5) and
$\mu\in\mathcal{M}(\Sigma)$, we want to maximize among lattices $L$ the
following functional:

$\mathcal{F}_{\mu}(L)=J_{M}(0)=\int_{\Sigma}\frac{\left|\nabla_{y}E_{f}[L+y]\right|^{2}}{E_{f}[L+y]}d\mu(y).$
(2.13)

Furthermore, if $f\geq 0$, then we have the following more compact form

$\mathcal{F}_{\mu}(L)=\int_{\Sigma}\left|\nabla_{y}\sqrt{E_{f}[L+y]}\right|^{2}d\mu(y)$
(2.14)

which gives us (1.3) by choosing $f(r)=e^{-\pi\alpha r}$ and renaming the
functional $\mathcal{F}_{\mu}^{\alpha}$ to show the $\alpha$-dependence.

For the reader’s convenience, the elementary form of
$\mathcal{F}_{\mu}^{\alpha}(L)$ is

$\mathcal{F}_{\mu}^{\alpha}(L)=4\pi^{2}\alpha^{2}\displaystyle\bigints_{\Sigma}\frac{\displaystyle\sum_{i=1}^{d}\left(\sum_{p\in
L}(p_{i}+y_{i})e^{-\pi\alpha|p+y|^{2}}\right)^{2}}{\displaystyle\sum_{p\in
L}e^{-\pi\alpha|p+y|^{2}}}d\mu(y).$

######  Remark 2.1 (Fisher Information in terms of the Heat Kernel).

It has to be noticed that $\mathcal{F}_{\mu}^{\alpha}$ can be expressed in
terms of the Heat Kernel on $L$. More precisely, and as we already pointed out
in [5, 13], let $u_{L}(y,t)$ be the temperature at point $y$ and at time $t$,
if at initial time $t=0$ a heat source of unit strength is placed at each
point of $L$, i.e. $u_{L}$ is defined, for any lattice
$L\subset\mathbb{R}^{d}$, any $y\in\mathbb{R}^{d}$ and any $t>0$ as the
solution of

$\displaystyle\left\\{\begin{array}[]{ll}\partial_{t}u_{L}(y,t)=\Delta_{y}u_{L}&\text{for
$(y,t)\in\mathbb{R}^{d}\times(0,\infty)$}\vspace{0.3ex}\\\
u_{L}(y,0)=\sum_{p\in L}\delta_{p}&\text{for
$y\in\mathbb{R}^{d}$},\end{array}\right.$

where $\delta_{p}$ is the Dirac measure at $p\in\mathbb{R}^{d}$. Then, for any
$\alpha>0$, any lattice $L\subset\mathbb{R}^{d}$ and any measure
$\mu\in\mathcal{M}(\Sigma)$,

$\mathcal{F}_{\mu}^{\alpha}(L)=\frac{1}{\alpha^{\frac{d}{2}}}\int_{\Sigma}\left|\nabla_{y}\sqrt{u_{L}\left(y,\frac{1}{4\pi\alpha}\right)}\right|^{2}d\mu(y).$
(2.15)

We now show the scaling formula already known by Mathis et al. in [35, Eq
(22)] and applied to the Gaussian periodic tuning curve. It is crucial to
remember that, as recalled above and according to [38, 47], changing the
density of the lattice $L$ also changes the size of the firing field $\Sigma$
accordingly.

######  Proposition 2.2 (Scaling formula).

For any $\lambda>0$, any $\mu\in\mathcal{M}(\Sigma)$, any
$L\in\mathcal{L}_{d}(1)$ and any $\alpha>0$, we have

$\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\lambda
L)=\lambda^{-2}\mathcal{F}^{\lambda^{2}\alpha}_{\mu}(L).$ (2.16)

###### Proof.

By a simple change of variable, we have

$\displaystyle\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\lambda L)$
$\displaystyle=\int_{\lambda\Sigma}\left|\nabla_{y}\sqrt{\theta_{\lambda
L+y}(\alpha)}\right|^{2}d\mu_{\lambda}(y)$
$\displaystyle=\frac{1}{\lambda^{2}}\int_{\Sigma}\left|\nabla_{y}\sqrt{\theta_{\lambda
L+\lambda y}(\alpha)}\right|^{2}d\mu(y)$
$\displaystyle=\frac{1}{\lambda^{2}}\int_{\Sigma}\left|\nabla_{y}\sqrt{\theta_{L+y}(\lambda^{2}\alpha)}\right|^{2}d\mu(y)=\lambda^{-2}\mathcal{F}^{\lambda^{2}\alpha}_{\mu}(L).$

∎

Since for any $\alpha>0$, any $\mu\in\mathcal{M}(\Sigma)$ and any
$L\in\mathcal{L}_{d}(1)$, $\lim_{\lambda\to
0}\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\lambda L)=\lim_{\lambda\to
0}\lambda^{-2}\mathcal{F}^{\lambda^{2}\alpha}_{\mu}(L)=+\infty$, it is
sufficient to look for the maximizer of $\mathcal{F}_{\mu}^{\alpha}$ for fixed
$\alpha$, $\mu\in\mathcal{M}(\Sigma)$ and $L\in\mathcal{L}_{d}(1)$. If the
density is not fixed, then it is then sufficient to take the high density
limit to find $+\infty$ as the maximum of $\mathcal{F}_{\mu}^{\alpha}$. Notice
that we also have, for any $\alpha>0$, any $\mu\in\mathcal{M}(\Sigma)$ and
$L\in\mathcal{L}_{d}(1)$,
$\lim_{\lambda\to\infty}\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\lambda
L)=\lim_{\lambda\to\infty}\lambda^{-2}\mathcal{F}^{\lambda^{2}\alpha}_{\mu}(L)=0$.

##  3 Alternative formula and proof of Theorem 1.2

We are proving Theorem 1.2 in this part, using an alternative formula given in
[40].

######  Proposition 3.1 (Alternative formula for
$\mathcal{F}_{\mu}^{\alpha}$).

For any $\alpha>0$, any $\mu\in\mathcal{M}(\Sigma)$ and any
$L\in\mathcal{L}_{d}(1)$,

$\mathcal{F}^{\alpha}_{\mu}(L)=4\pi^{2}\alpha^{2}\int_{\Sigma}\theta_{L+y}(\alpha)\left|\underset{w\sim
D_{L+y,\alpha}}{\mathbb{E}[w]}\right|^{2}d\mu(y),$ (3.1)

where $D_{L+x,\alpha}$ is the discrete Gaussian distribution over $L+x$ with
parameter $\alpha$, i.e. the probability distribution over $L+x$ that assigns
probability proportional to $e^{-\pi\alpha|w|^{2}}$ to each vector $w\in L+x$.

######  Remark 3.2.

Let us specify that the above expected value is actually the following vector:

$\underset{w\sim
D_{L+y,\alpha}}{\mathbb{E}[w]}=\frac{1}{\theta_{L+y}(\alpha)}\sum_{p\in
L}(p+y)e^{-\pi\alpha|p+y|^{2}}.$ (3.2)

###### Proof.

Following [40, p. 6] (with $s=\alpha^{-1/2}$), a straightforward computation
gives (see (3.2))

$\frac{\nabla_{x}\theta_{L+x}(\alpha)}{\theta_{L+x}(\alpha)}=-2\pi\alpha\underset{w\sim
D_{L+x,\alpha}}{\mathbb{E}[w]}.$ (3.3)

It follows that

$\frac{|\nabla_{y}\theta_{L+y}(\alpha)|^{2}}{\theta_{L+y}(\alpha)}=4\pi^{2}\alpha^{2}\theta_{L+y}(\alpha)\left|\underset{w\sim
D_{L+y,\alpha}}{\mathbb{E}[w]}\right|^{2},$ (3.4)

which shows (3.1). ∎

We now show Theorem 1.2 by directly applying (3.1).

###### Proof of Theorem 1.2.

A straightforward computation shows that, for any $\mu\in\mathcal{M}(\Sigma)$
and any $\lambda>0$, as $\alpha\to 0$,

$\mathcal{F}^{\alpha}_{\mu_{\lambda}}(L)=4\pi^{2}\alpha^{2}\int_{\Sigma}\theta_{L+y}(\alpha)\left|\underset{w\sim
D_{L+y,\alpha}}{\mathbb{E}[w]}\right|^{2}d\mu_{\lambda}(y)\sim
4\pi^{2}\alpha^{2}\int_{\Sigma}\theta_{L+y}(\alpha)d\mu_{\lambda}(y),$ (3.5)

which is uniquely minimized by $\mathsf{A}_{2}$ (resp. locally minimized by
$\mathsf{D}_{3}^{*}$) by a direct application of [14, Thm. 2 and 3] (resp.
[11, Thm. 4.3]) as $\lambda$ is small enough or if $\mu$ has the form
$d\mu(x)=\rho(|x|^{2})dx$ and $\rho$ is a completely monotone function. The
fact that $\mathcal{F}^{\alpha}_{\mu_{\lambda}}$ does not have any maximizer
as $\alpha,\lambda$ are small enough is again a simple consequence of our work
in [11] and the fact that the lattice theta function
$L\mapsto\theta_{L}(\alpha)$ does not have a maximizer in
$\mathcal{L}_{d}(1)$. Indeed, it is enough to take a sequence of orthorhombic
lattices with degenerate shapes for showing that the theta functions goes to
infinity. ∎

######  Remark 3.3.

The same is true in dimension $d\in\\{8,24\\}$ replacing $\mathsf{A}_{2}$ by
$\mathsf{E}_{8}$ or $\Lambda_{24}$, again by an application of [11, Thm 4.5]
and the universal optimality of these lattices shown in [21].

##  4 Volume-stationary lattices - Proof of Theorem 1.4

We now show Theorem 1.4 by using a combination of Proposition 2.2, formula
(3.1) and our works in [14, 10].

###### Proof of Theorem 1.4.

We consider $\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\lambda L)$ for $\lambda>0$,
then by assumption we have that
$\nabla_{L}\mathcal{F}^{\alpha}_{\mu_{\lambda}}(\lambda L)=0$ for all
$\lambda\in I$, where $I$ is an open interval of $\mathbb{R}_{+}$. We
therefore have, for all $\lambda\in I$,

$\nabla_{L}\int_{\Sigma}\left|\nabla_{y}\sqrt{\theta_{\lambda
L+y}(\alpha)}\right|^{2}d\mu(y)=0,$ (4.1)

the gradient $\nabla_{L}$ being taken in $\mathcal{L}_{d}(1)$. Since
$\lambda\mapsto\nabla_{L}\mathcal{F}^{\alpha}_{\mu}(\lambda L)$ is analytic by
integration on a compact set and absolute convergence, its zero are therefore
isolated and it follows that necessarily $I=(0,\infty)$. We now use formula
(3.1) combined with Proposition 2.2 and we get

$\displaystyle\mathcal{F}_{\mu_{\lambda}}^{\alpha}(\lambda
L)=\lambda^{-2}\mathcal{F}_{\mu}^{\lambda^{2}\alpha}(L)=4\pi^{2}\lambda^{2}\alpha^{2}\int_{\Sigma}\theta_{L+y}(\lambda^{2}\alpha)\left|\underset{w\sim
D_{L+y,\lambda^{2}\alpha}}{\mathbb{E}[w]}\right|^{2}d\mu(y).$

Therefore, since the above squared expectation converges to a constant for
small $\lambda^{2}\alpha$, we obtain for $\lambda<\lambda_{0}$ sufficiently
small,

$\nabla_{L}\left\\{\int_{\Sigma}\theta_{L+y}(\lambda^{2}\alpha)d\mu(y)\right\\}=0.$
(4.2)

By Poisson Summation Formula (see e.g. [22, Eq. (43)]), we obtain as in [14,
11] that, for all $\lambda,\alpha>0$,

$\displaystyle\int_{\Sigma}\theta_{L+y}(\lambda^{2}\alpha)d\mu(y)=\int_{\Sigma}\frac{1}{(\lambda^{2}\alpha)^{\frac{d}{2}}}\sum_{q\in
L^{*}}e^{-\frac{\pi|q|^{2}}{\lambda^{2}\alpha}}e^{2i\pi q\cdot y}d\mu(y)$
$\displaystyle=\frac{1}{(\lambda^{2}\alpha)^{\frac{d}{2}}}\sum_{q\in
L^{*}}e^{-\frac{\pi|q|^{2}}{\lambda^{2}\alpha}}\int_{\Sigma}e^{2i\pi q\cdot
y}d\mu(y)$
$\displaystyle=\frac{1}{(\lambda^{2}\alpha)^{\frac{d}{2}}}\sum_{q\in
L^{*}}e^{-\frac{\pi|q|^{2}}{\lambda^{2}\alpha}}\hat{\mu}(q),$

where $\hat{\mu}$ is the Fourier transform of the measure $\mu$. By
analyticity of
$\lambda\mapsto\nabla_{L}\int_{\Sigma}\theta_{L+y}(\lambda^{2}\alpha)d\mu(y)$,
we therefore get again that, for all $\lambda>0$,

$\nabla_{L}\sum_{q\in
L^{*}}e^{-\frac{\pi}{\alpha\lambda^{2}}|q|^{2}}\widehat{\mu}(q)=0.$ (4.3)

Since $\mu$ is radially symmetric, then its Fourier transform is radially
symmetric, i.e. $\widehat{\mu}(q)=h(|q|)$ for some function
$h:\mathbb{R}_{+}\to\mathbb{R}$. Now, the same argument can be used as in [25,
28, 10] to conclude that the first layer of $L^{*}$ must be strongly eutactic.
Indeed, (4.3) implies that, for all $t>0$ (which could be written
$t=\lambda^{2}\alpha$),

$\nabla_{L^{*}}\theta_{L^{*}}\left(\frac{1}{t}\right)=0,$

which implies that

$\nabla_{L^{*}}\zeta_{L^{*}}(s)=0$

for all $s>d$ for the Epstein zeta function defined by

$\zeta_{L^{*}}(s):=\sum_{q\in
L^{*}\backslash\\{0\\}}\frac{1}{|q|^{s}}=\int_{0}^{\infty}\left(\theta_{L^{*}}\left(\frac{t}{\pi}\right)-1\right)\frac{t^{\frac{s}{2}-1}}{\Gamma(s/2)}dt,$

where the second equality follows from a simple Laplace transform argument
(see e.g. [7]). We therefore apply [28, Corollary 1] to conclude that all the
layers of $L^{*}$ are necessarily strongly eutactic. In dimension
$d\in\\{2,3\\}$, all these lattices are characterized (see for instance [28,
Corollary 2]) and they are the one listed in the theorem’s statement for which
the list of duals stays unchanged. This completes the proof. ∎

######  Remark 4.1.

The proof still holds if the Gaussian is replaced by a function $f$ which is
the Laplace transform of a measure. We also suspect that the result is also
true for some non-radially symmetric measure $\mu$.

##  5 Numerical investigation

We are using the software Scilab to perform our numerical investigation. All
the values are computed by using a sufficient number of points for estimating
our infinite sums according to the value of $\alpha$ we choose (mostly
$\alpha=\frac{10}{\pi}\approx 3.18$).

###  5.1 Parametrization of lattices in dimension 2 and 3

We recall here how to parametrize a $d$-dimensional lattice in order to
perform our numerical computations. We basically follow the lines of [8, 9] by
stating that a two-dimensional (resp. three-dimensional) lattice of unit
density can be parametrized by two real numbers $(x,y)$ (resp. five real
numbers $(u,v,x,y,z)$) belonging to a fundamental domain
$\mathcal{D}_{2}\subset\mathbb{R}^{2}$ (resp.
$\mathcal{D}_{3}\subset\mathbb{R}^{5}$) where only one copy of each lattice
appears (see e.g. [49, Sect. 1.4]). More precisely, two-dimensional and three-
dimensional lattices of unit density can be written respectively as

$\displaystyle
L=\mathbb{Z}\left(\frac{1}{\sqrt{y}},0\right)\oplus\mathbb{Z}\left(\frac{x}{\sqrt{y}},\sqrt{y}\right)\in\mathcal{L}_{2}(1)$
$\displaystyle
L=2^{\frac{1}{6}}\left[\mathbb{Z}\left(\frac{1}{\sqrt{u}},0,0\right)\oplus\mathbb{Z}\left(\frac{x}{\sqrt{u}},\frac{v}{\sqrt{u}},0\right)\oplus\mathbb{Z}\left(\frac{y}{\sqrt{u}},\frac{vz}{\sqrt{u}},\frac{u}{v\sqrt{2}}\right)\right]\in\mathcal{L}_{3}(1).$

In dimension $d=2$, this fundamental domain $\mathcal{D}_{2}$ is easy to
describe whereas it is more complicated in dimension $d=3$. In particular, we
have

$\mathcal{D}_{2}:=\left\\{(x,y)\in\mathbb{R}^{d}:0\leq
x\leq\frac{1}{2},y>0,x^{2}+y^{2}\geq 1\right\\}.$

The square lattice $\mathbb{Z}^{2}$ and triangular lattice $\mathsf{A}_{2}$
are respectively represented by the points $(0,1)\in\mathcal{D}_{2}$ and
$\left(\frac{1}{2},\frac{\sqrt{3}}{2}\right)\in\mathcal{D}_{2}$. We omit to
write what an example of fundamental domain of $\mathcal{L}_{3}(1)$ could be,
but we recall that the simple cubic lattice $\mathbb{Z}^{3}$, the FCC lattice
$\mathsf{D}_{3}$ and the BCC lattice $\mathsf{D}_{3}^{*}$ are respectively
represented by $(2^{1/3},1,0,0,0)\in\mathcal{D}_{3}$,
$(1,1,0,1/2,1/2)\in\mathcal{D}_{3}$ and
$(2^{-1/3},1,0,1/2,1/2)\in\mathcal{D}_{3}$. We recall that, since the three-
dimensional case is too demanding in terms of computational time, we will only
compare certain lattices and compute the Hessian matrices at the FCC lattice.

###  5.2 The pure discrete case: combinations of $\mathcal{Q}_{L}(y)$

We first explain why the problem of maximizing the discrete version of the
Fisher Information’s trace (2.11) (choosing again $x=0$), i.e.

$\mathrm{Tr}J(0)=n\sum_{j=1}^{N}\mathcal{Q}_{L}(y_{j}),\quad\mathcal{Q}_{L}(y):=\frac{\left|\nabla_{y}\theta_{L+y}(\alpha)\right|^{2}}{\theta_{L+y}(\alpha)},$

is very difficult and therefore why it is convenient to consider – as
suggested by experiments [29, 42, 37] – a continuous firing field $\Sigma$
instead of a set of points $Y$. We only consider the two-dimensional case, but
the same remarks hold for higher dimensions.

For any $L=\bigoplus_{i=1}^{d}\mathbb{Z}u_{i}$, if for all
$j\in\\{1,...,N\\}$,
$y_{j}\in\\{u_{i},\frac{u_{i}}{2},\frac{1}{2}\sum_{i}u_{i}\\}$ for any
$i\in\\{1,...,d\\}$, then $\mathrm{Tr}J(0)=0$ since these points $y_{j}$ are
known to be critical points for $y\mapsto\theta_{L+y}(\alpha)$ – and more
generally $y\mapsto E_{f}[L+Y]$ for any function $f$ satisfying (1.5). They
are the trivial zeros of $\mathrm{Tr}J(0)$. In the particular case
$L=\mathsf{A}_{2}$, the barycenters of the primitive triangles are also
critical points of $y\mapsto\theta_{L+y}(\alpha)$ for all $\alpha>0$ as shown
in [5] (see Figure 4). Furthermore, in Figure 3, we show that the maximizer of
$\mathrm{Tr}J(0)$ varies a lot with $\\{y_{j}\\}_{j}$.

Figure 3: Plot of
$t\mapsto\mathcal{Q}^{\alpha}_{\mathsf{A}_{2}}(0.1,0.1)+\mathcal{Q}^{\alpha}_{\mathsf{A}_{2}}(t,0.2)-\left(\mathcal{Q}^{\alpha}_{\mathbb{Z}^{2}}(0.1,0.1)+\mathcal{Q}^{\alpha}_{\mathbb{Z}^{2}}(t,0.2)\right)$
for $\alpha=\frac{10}{\pi}$ and $t\in[0,1]$, showing that the maximizer of
$\mathrm{Tr}J_{M}(0)$ is not always triangular. Figure 4: Zeros of
$x\mapsto\nabla_{x}\theta_{\mathsf{A}_{2}+x}(\alpha)$ in a patch of
$\mathsf{A}_{2}$. Following [5], $\bullet$ are the maximizers, $\circ$ the
minimizers and $\times$ the saddle points of the function
$x\mapsto\theta_{\mathsf{A}_{2}+x}(\alpha)$. These are the only critical
points of this function.

It has also to be noticed that
$L\mapsto\displaystyle\mathcal{Q}^{\alpha}_{L}(y)$ alone is not always
maximized by the same lattice as $y$ varies in $\Sigma=B_{R}$. Indeed, writing
$y=(r\cos\theta,r\sin\theta)$ and fixing $r$, we compare the energies of the
square and triangular lattices by plotting

$\theta\mapsto
g_{r}(\theta):=\mathcal{Q}_{\mathbb{Z}^{2}}^{\alpha}(r\cos\theta,r\sin\theta)-\mathcal{Q}_{\mathsf{A}_{2}}^{\alpha}(r\cos\theta,r\sin\theta),\quad\alpha=\frac{10}{\pi}.$
(5.1)

It is easy to give an explicit formula for $g_{r}(\theta)$ since this is the
difference of

$\displaystyle\mathcal{Q}_{\mathbb{Z}^{2}}^{\alpha}(r\cos\theta,r\sin\theta)=$
$\displaystyle
4\pi^{2}\alpha^{2}\frac{\displaystyle\left(\sum_{m,n\in\mathbb{Z}}(m+r\cos\theta)e^{-\pi\alpha[(m+r\cos\theta)^{2}+(n+r\sin\theta)^{2}]}\right)^{2}}{\displaystyle\sum_{m,n\in\mathbb{Z}}e^{-\pi\alpha[(m+r\cos\theta)^{2}+(n+r\sin\theta)^{2}]}}$
$\displaystyle\quad\quad\quad+4\pi^{2}\alpha^{2}\frac{\displaystyle\left(\sum_{m,n\in\mathbb{Z}}(n+r\sin\theta)e^{-\pi\alpha[(m+r\cos\theta)^{2}+(n+r\sin\theta)^{2}]}\right)^{2}}{\displaystyle\sum_{m,n\in\mathbb{Z}}e^{-\pi\alpha[(m+r\cos\theta)^{2}+(n+r\sin\theta)^{2}]}}$

and

$\displaystyle\mathcal{Q}_{\mathsf{A}_{2}}^{\alpha}(r\cos\theta,r\sin\theta)=$
$\displaystyle
4\pi^{2}\alpha^{2}\frac{\displaystyle\left(\sum_{m,n\in\mathbb{Z}}\left(\lambda(m+\frac{n}{2})+r\cos\theta\right)e^{-\pi\alpha[(\lambda(m+\frac{n}{2})+r\cos\theta)^{2}+(\frac{\lambda\sqrt{3}}{2}n+r\sin\theta)^{2}]}\right)^{2}}{\displaystyle\sum_{m,n\in\mathbb{Z}}e^{-\pi\alpha[(\lambda(m+\frac{n}{2})+r\cos\theta)^{2}+(\frac{\lambda\sqrt{3}}{2}n+r\sin\theta)^{2}]}}$
$\displaystyle\quad\quad\quad+4\pi^{2}\alpha^{2}\frac{\displaystyle\left(\sum_{m,n\in\mathbb{Z}}\left(\frac{\lambda\sqrt{3}}{2}n+r\sin\theta\right)e^{-\pi\alpha[(\lambda(m+\frac{n}{2})+r\cos\theta)^{2}+(\frac{\lambda\sqrt{3}}{2}n+r\sin\theta)^{2}]}\right)^{2}}{\displaystyle\sum_{m,n\in\mathbb{Z}}e^{-\pi\alpha[(\lambda(m+\frac{n}{2})+r\cos\theta)^{2}+(\frac{\lambda\sqrt{3}}{2}n+r\sin\theta)^{2}]}}.$

Our observations are the following (see Figure 5):

  * •

For small $r$, we observe that $g_{r}(\theta)<0$ for all $\theta\in[0,2\pi]$.
The triangular lattice has higher energy than the square lattice.

  * •

For $r>0.24$, there are alternation of signs for $g_{r}$ on some intervals of
$[0,2\pi]$.

  * •

For $r=0.6$, the positive part of $g_{r}$ is more important than the negative
one.

Therefore, showing the maximality of $\mathsf{A}_{2}$ for the integral of
$\mathcal{Q}_{L}^{\alpha}$ on a ball of radius larger than $R_{0}=0.24$ is
necessarily tricky. Furthermore, it is not surprising that the triangular
lattice is almost never a critical point of the integrand in
$\mathcal{L}_{2}(1)$. To illustrate this fact, we have plotted in Figure 6 the
function
$\theta\mapsto|\nabla_{L}\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)|^{2}$
on $[0,2\pi]$ for $L=\mathsf{A}_{2}$, $\alpha=10/\pi$ and fixed
$r\in\\{0.2,0.3\\}$. Notice that it is straightforward to compute
$|\nabla_{L}\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)|^{2}$, recalling
that $\mathsf{A}_{2}$ corresponds to the point $(x,y)=(1/2,\sqrt{3}/2)$, since

$\displaystyle|\nabla_{L}\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)|_{L=\mathsf{A}_{2}}|^{2}$
$\displaystyle=\left(\partial_{x}\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)|_{(x,y)=(1/2,\sqrt{3}/2)}\right)^{2}+\left(\partial_{y}\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)|_{(x,y)=(1/2,\sqrt{3}/2)}\right)^{2},$

where

$\displaystyle\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)$
$\displaystyle=4\pi^{2}\alpha^{2}\frac{\displaystyle\left(\sum_{m,n\in\mathbb{Z}}\left(\frac{(m+xn)}{\sqrt{y}}+r\cos\theta\right)e^{-\pi\alpha\left[\left(\frac{(m+xn)}{\sqrt{y}}+r\cos\theta\right)^{2}+(n\sqrt{y}+r\sin\theta)^{2}\right]}\right)^{2}}{\displaystyle\sum_{m,n\in\mathbb{Z}}e^{-\pi\alpha\left[\left(\frac{(m+xn)}{\sqrt{y}}+r\cos\theta\right)^{2}+(n\sqrt{y}+r\sin\theta)^{2}\right]}}$
$\displaystyle+4\pi^{2}\alpha^{2}\frac{\displaystyle\left(\sum_{m,n\in\mathbb{Z}}(n\sqrt{y}+r\sin\theta)e^{-\pi\alpha\left[\left(\frac{(m+xn)}{\sqrt{y}}+r\cos\theta\right)^{2}+(n\sqrt{y}+r\sin\theta)^{2}\right]}\right)^{2}}{\displaystyle\sum_{m,n\in\mathbb{Z}}e^{-\pi\alpha\left[\left(\frac{(m+xn)}{\sqrt{y}}+r\cos\theta\right)^{2}+(n\sqrt{y}+r\sin\theta)^{2}\right]}}.$

Figure 5: Plots of $g_{r}$ (see (5.1)) on $[0,2\pi]$ for
$\alpha=\frac{10}{\pi}$, $r\in\\{0.1,0.3,0.5,0.6\\}$ in reading direction.

Figure 6: Plot of
$\theta\mapsto|\nabla_{L}\mathcal{Q}^{\alpha}_{L}(r\cos\theta,r\sin\theta)|^{2}$
on $[0,2\pi]$ for $L=\mathsf{A}_{2}$, $\alpha=\frac{10}{\pi}$ and fixed
$r\in\\{0.2,0.3\\}$.

###  5.3 Numerical results for $\mathcal{F}_{\mu}^{\alpha}$ in dimension 2

We are considering the case $\alpha=\frac{10}{\pi}$ (see Figure 7), in such a
way that

$\theta_{L+y}(\alpha)=\sum_{p\in L}e^{-10|p+y|^{2}}.$

We are also assuming that $\mu=\sigma_{R}$ is the uniform measure on the disk
$\Sigma=B_{R}$ for some $R>0$. For simplicity, we omit to renormalize the
measure, in such a way that $\mu$ is almost never a probability measure. It
does not matter in our study since we are considering optimizers at fixed $R$.

Figure 7: Plot of $r\mapsto f(r^{2})=e^{-10r^{2}}$.

We observe that for $R\in\\{0.1,0.2,0.3,0.4,0.5\\}$, the triangular lattice
$\mathsf{A}_{2}$ is the unique maximizer of $\mathcal{F}_{\mu}^{\alpha}$ in
$\mathcal{L}_{2}(1)$. In Figure 8, we give one example of our numerical
outputs as a contour plot for the particular case $R=0.5$. We have also
checked that the same maximality result holds for
$R=\frac{1}{2}\sqrt{\frac{2}{\sqrt{3}}}$ being the half side length of
$\mathsf{A}_{2}$, i.e. the radius of the balls reaching the densest lattice
packing of unit density in $\mathbb{R}^{2}$. If $R\geq 0.59$, the square
lattice has a higher energy than the triangular one and
$\mathcal{F}_{\mu}^{\alpha}$ does not have any maximizer.

In order to illustrate Theorem 1.2, we choose $\alpha=\frac{2}{\pi}$,
$\Sigma=B_{R}$ $\mu=\sigma_{R}$, and we observe in Figure 9) that for $R=0.1$,
$\mathsf{A}_{2}$ is the unique minimizer of $\mathcal{F}_{\mu}^{\alpha}$ in
$\mathcal{L}_{2}(1)$ and the functional does not have any maximizer.

Figure 8: Case $\mu=\sigma_{R}$ uniform on $\Sigma=B_{R}$, $R=0.5$ and
$\alpha=\frac{10}{\pi}$. The triangular lattice $\mathsf{A}_{2}$,
corresponding to the point $(1/2,\sqrt{3}/2)$, is the unique maximizer of
$\mathcal{F}_{\mu}^{\alpha}$ in $\mathcal{L}_{2}(1)$.

Figure 9: Illustration of Theorem 1.2. Case $\mu=\sigma_{R}$ uniform on
$\Sigma=B_{R}$, $R=0.1$ and $\alpha=\frac{2}{\pi}$. There is no maximizer
(left). Furthermore, the triangular lattice $\mathsf{A}_{2}$, corresponding to
the point $(1/2,\sqrt{3}/2)$, is the unique minimizer of
$\mathcal{F}_{\mu}^{\alpha}$ in $\mathcal{L}_{2}(1)$.

If we compare the Fisher Information of the two-dimensional volume-stationary
lattices, i.e. $\mathsf{A}_{2}$ and $\mathbb{Z}^{2}$ for fixed $\alpha=10/\pi$
and increasing $R$, we obtain the graph of Figure 10. We observe that, for
$R\in[0.1,0.57]$,
$\mathcal{F}_{\sigma_{R}}^{\alpha}(\mathsf{A}_{2})>\mathcal{F}_{\sigma_{R}}^{\alpha}(\mathbb{Z}^{2})$
whereas the square lattice has a higher Fisher Information for
$R\in[0.58,0.7]$. We also remark that the firing fields start to overlap when
$R>\frac{1}{2}\sqrt{\frac{2}{\sqrt{3}}}\approx 0.5373$.

Figure 10: Plot of
$R\mapsto\mathcal{F}_{\sigma_{R}}^{\alpha}(\mathsf{A}_{2})-\mathcal{F}_{\sigma_{R}}^{\alpha}(\mathbb{Z}^{2})$
where $\alpha=\frac{10}{\pi}$ and $\sigma_{R}$ is the uniform measure on
$\Sigma=B_{R}$, for $R\in[0.1,0.7]$.

We have also numerically checked the case where the ratio between $2R$ and the
lattice spacing of the triangular lattice is approximately $0.3$. Applied to
the triangular lattice $\mathsf{A}_{2}$ of side length $\sqrt{2/\sqrt{3}}$, it
means that $R=0.16$. For this value of $R$, we have found that the triangular
lattice is always maximal for $\alpha\in(1.25,+\infty)$ (see Figure 11).
Notice that even though a large zone of these figures seems to reach the
maximum around the triangular lattice, we have checked that the minimum is
really given by the latter.

Figure 11: Plot of $\mathcal{F}_{\sigma_{R}}^{\alpha}$ where $R=0.16$ and
$\alpha\in\\{1.3,3,5\\}$. The triangular lattice appears to be a maximizer in
all these cases.

###  5.4 Numerical results for $\mathcal{F}_{\mu}^{\alpha}$ in dimension 3

For our investigation in dimension $d=3$, we have again chosen
$\alpha=\frac{10}{\pi}$. Since the numerical computations are very demanding,
we are only giving the following type of results: some local optimality
results and some comparisons of volume-stationary lattices, i.e.
$\mathbb{Z}^{3},\mathsf{D}_{3}$ and $\mathsf{D}_{3}^{*}$.

As in the previous section on two-dimensional lattices, it is clear that
$\mathsf{D}_{3}$, $\mathsf{D}_{3}^{*}$ and $\mathbb{Z}^{3}$ are almost never
critical points for $L\mapsto\mathcal{Q}_{L}^{\alpha}(y)$ for a given $y$. It
is only true for some points $y$ being the holes or deep holes of $L$ (center
of cells, midpoints, etc.) as we have already explained for $\mathsf{A}_{2}$
(see Figure 3). It is therefore not sufficient to study
$L\mapsto\mathcal{Q}_{L}^{\alpha}(y)$ in order to investigate the extrema of
$\mathcal{F}_{\mu}^{\alpha}$.

By comparing $\mathcal{Q}_{L}^{\alpha}(y)$ for $y\in\Sigma=B_{R}$, we observe
that there exists $R_{0}\approx 0.1$ such that

$\mathcal{Q}_{\mathsf{D}_{3}}(y)>\mathcal{Q}_{\mathsf{D}_{3}^{*}}(y)>\mathcal{Q}_{\mathbb{Z}^{3}}(y),\quad\forall
y\in B_{R_{0}},$

which directly means that, for $\mu=\sigma_{R_{0}}$,

$\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3})>\mathcal{F}_{\mu}^{\alpha}(\mathsf{D}_{3}^{*})>\mathcal{F}_{\mu}^{\alpha}(\mathbb{Z}^{3}).$
(5.2)

We have also checked on a set of discrete values of $R$ that (5.2) still holds
for $R\in[0.1,0.57]$. Notice that $0.57>2^{-5/6}=:R_{\mathsf{D}_{3}}$ which is
the radius of spheres giving the densest packing of unit density in
$\mathbb{R}^{3}$, on a FCC lattice.

Furthermore, we have also numerically checked the local maximality of
$\mathsf{D}_{3}$ for radii $R\in\\{0.1,0.2,0.3,0.4,0.5,2^{-\frac{5}{6}}\\}$ by
computing the Hessian of $\mathcal{F}_{\mu}^{\alpha}$ at $L=\mathsf{D}_{3}$
which appears to be definite positive in all these cases.

Acknowledgments: I would like to thank the Austrian Science Fund (FWF) for its
financial support through the project F65 as well as the DFG-FWF international
joint project FR 4083/3-1/I 4354. I would also like to thank Martin Stemmler
for helpful discussions about the size of the firing fields. I am also
grateful to the anonymous referees for their suggestions.

## References

  * [1] A. Abrikosov.  The Magnetic Properties of Superconducting Alloys.  Journal of Physics and Chemistry of Solids, 2:199–208, 1957.
  * [2] F. Anselmi, M.M. Murray, and B. Franceschiello.  A computationalmodel for grid maps in neural populations.  Journal of Computational Neuroscience, 48:149–159, 2020.
  * [3] M. Antlanger, G. Kahl, M. Mazars, L. Samaj, and E. Trizac.  Rich polymorphic behavior of Wigner bilayers.  Physical Review Letters, 117(11):118002, 2016.
  * [4] C. Bachoc and B. Venkov.  Modular Forms, Lattices and Spherical Designs.  Réseaux euclidiens, designs sphériques et groupes, L’Ens. Math.(Monographie 37, J. Martinet, ed.):87–111, 2001.
  * [5] A. Baernstein II.  A minimum problem for heat kernels of flat tori.  Contemp. Math., 201:227–243, 1997.
  * [6] S. Bernstein.  Sur les fonctions absolument monotones.  Acta Math., 52:1–66, 1929.
  * [7] L. Bétermin.  Two-dimensional Theta Functions and Crystallization among Bravais Lattices.  SIAM J. Math. Anal., 48(5):3236–3269, 2016.
  * [8] L. Bétermin.  Local variational study of 2d lattice energies and application to Lennard-Jones type interactions.  Nonlinearity, 31(9):3973–4005, 2018.
  * [9] L. Bétermin.  Local optimality of cubic lattices for interaction energies.  Anal. Math. Phys., 9(1):403–426, 2019.
  * [10] L. Bétermin.  Minimizing lattice structures for Morse potential energy in two and three dimensions.  J. Math. Phys., 60(10):102901, 2019.
  * [11] L. Bétermin.  Minimal Soft Lattice Theta Functions.  Constr. Approx., 52(1):115–138, 2020.
  * [12] L. Bétermin and M. Faulhuber.  Maximal theta functions - Universal optimality of the hexagonal lattice for Madelung-like lattice energies.  Preprint. arXiv:2007.15977, 2020.
  * [13] L. Bétermin and H. Knüpfer.  On Born’s conjecture about optimal distribution of charges for an infinite ionic crystal.  J. Nonlinear Sci., 28(5):1629–1656, 2018.
  * [14] L. Bétermin and H. Knüpfer.  Optimal lattice configurations for interacting spatially extended particles.  Lett. Math. Phys., 108(10):2213–2228, 2018.
  * [15] L. Bétermin and M. Petrache.  Dimension reduction techniques for the minimization of theta functions on lattices.  J. Math. Phys., 58:071902, 2017.
  * [16] L. Bétermin and M. Petrache.  Optimal and non-optimal lattices for non-completely monotone interaction potentials.  Anal. Math. Phys., 9(4):2033–2073, 2019.
  * [17] X. Blanc and M. Lewin.  The Crystallization Conjecture: A Review.  EMS Surv. in Math. Sci., 2:255–306, 2015.
  * [18] V.H. Brun, T. Solstad, K.B. Kjelstrup, M. Fyhn, M.P. Witter, E.I. Moser, and M.B. Moser.  Progressive increase in grid scale from dorsal to ventral medial entorhinal cortex.  Hippocampus, 18:1200–1212, 2008.
  * [19] H. Cohn and A. Kumar.  Universally optimal distribution of points on spheres.  J. Amer. Math. Soc., 20(1):99–148, 2007.
  * [20] H. Cohn, A. Kumar, S. D. Miller, D. Radchenko, and M. Viazovska.  The sphere packing problem in dimension 24.  Ann. of Math., 185(3):1017–1033, 2017.
  * [21] H. Cohn, A. Kumar, S. D. Miller, D. Radchenko, and M. Viazovska.  Universal optimality of the $E_{8}$ and Leech lattices and interpolation formulas.  to appear in Ann. of Math., 2019.  arXiv:1902:05438.
  * [22] J. H. Conway and N. J. A. Sloane.  Sphere Packings, Lattices and Groups, volume 290.  Springer, 1999.
  * [23] T. M. Cover and J. A. Thomas.  Elements of Information Theory.  Wiley, New York, second edition edition, 2006.
  * [24] T. D’Albis and R. Kempter.  A single-cell spiking model for the origin of grid-cell patterns.  PLoS Comput Biol, 13(10):e1005782, 2017.
  * [25] B. N. Delone and S. S. Ryshkov.  A Contribution to the Theory of the Extrema of a Multidimensional zeta-function.  Dokl. Akad. Nauk SSSR, 173(4):991–994, 1967.
  * [26] M. Faulhuber and S. Steinerberger.  Optimal Gabor frame bounds for separable lattices and estimates for Jacobi theta functions.  J. Math. Anal. Appl., 445(1):407–422, 2017.
  * [27] K. Gerlei, J. Passlack, I. Hawes, B. Vandrey, H. Stevens, I. Papastathopoulos, and M.F. Nolan.  Grid cells are modulated by local head direction.  Nature Communications, 11(4228), 2020.
  * [28] P. M. Gruber.  Application of an Idea of Voronoi to Lattice Zeta Functions.  Proc. Steklov Inst. Math., 276:103–124, 2012.
  * [29] T. Hafting, M. Fyhn, S. Molden, M.B. Moser, and E.I. Moser.  Microstructure of a spatial map in the entorhinal cortex.  Nature, 436:801–6, 2005.
  * [30] R. C. Heitmann and C. Radin.  The Ground State for Sticky Disks.  J. Stat. Phys., 22:281–287, 1980.
  * [31] L. Kantorovich.  Quantum Theory of the Solid State.  Fundamental Theories of Physics. Springer Netherlands, 1st edition edition, 2004.
  * [32] M. Kim and E.A. Maguire.  Can we study 3D grid codes non-invasively in the human brain? Methodological considerations and fMRI findings.  Neuroimage, 186:667–678, 2019.
  * [33] S. Luo, X. Ren, and J. Wei.  Non-hexagonal lattices from a two species interacting system.  SIAM J. Math. Anal., 52(2):1903–1942, 2020.  Preprint. arXiv:1902.09611.
  * [34] S. Luo and J. Wei.  On minima of sum of theta functions and Mueller-Ho Conjecture.  Preprint. arXiv:2004.13882, 2020.
  * [35] A. Mathis, M.B. Stemmler, and A.V.M. Herz.  Probable nature of higher-dimensional symmetries underlying mammalian grid-cell activity patterns.  eLife 2015.4, pages 1–29, 2015.
  * [36] H. L. Montgomery.  Minimal Theta Functions.  Glasg. Math. J., 30(1):75–85, 1988.
  * [37] J. Nagele, A.V.M. Herz, and M.B. Stemmler.  Untethered firing fields and intermittent silences: Why grid-cell discharge is so variable.  Hippocampus, DOI: 10.1002/hipo.23191, 2020.
  * [38] G.J. Quirk, R.U. Muller, J.L. Kubie, and J.B. Ranck Jr.  The positional firing properties of medial entorhinal neurons: Description and comparison with hippocampal place cells.  J Neurosci, 12:1945–1963, 1992.
  * [39] C. Radin.  Low temperature and the origin of crystalline symmetry.  International Journal of Modern Physics B, 1(5 and 6):1157–1191, 1987.
  * [40] O. Regev and N. Stephens-Davidowitz.  An Inequality for Gaussians on Lattices.  SIAM J. Discrete Math., 31(2):749–757, 02 2017.  Preprint. arXiv:1502.04796.
  * [41] S. Rosay, S. Weber, and M. Mulas.  Modeling grid fields instead of modeling grid cells.  Journal of Computational Neuroscience, 47:43–60, 2019.
  * [42] D.C. Rowland, Y. Roudi, M.B. Moser, and E.I. Moser.  Ten years of grid cells.  Annual Review of Neuroscience, 39(19-40), 2016.
  * [43] E. Sandier and S. Serfaty.  From the Ginzburg-Landau Model to Vortex Lattice Problems.  Comm. Math. Phys., 313(3):635–743, 2012.
  * [44] P. Sarnak and A. Strömbergsson.  Minima of Epstein’s Zeta Function and Heights of Flat Tori.  Invent. Math., 165:115–151, 2006.
  * [45] T. Solstad, E.I. Moser, and G.T. Einevoll.  From grid cells to place cells: A mathematical model.  Hippocampus, 16:1026–1031, 2006.
  * [46] B. Sorscher, G.C. Mel, S. Ganguli, and S.A. Ocko.  A unified theory for the origin of grid cells through the lens of pattern formation.  In Advances in Neural Information Processing Systems 32, pages 10003–10013. Curran Associates, Inc., 2019.
  * [47] H. Stensola, T. Stensola, T. Solstad, K. Frøland, M.B. Moser, and E.I. Moser.  The entorhinal grid map is discretized.  Nature, 492:72–78, 2012.
  * [48] T. Stensola, H. Stensola, M.B. Moser, and E.I. Moser.  Shearing-induced asymmetry in entorhinal grid cells.  Nature, 518:207–212, 2015.
  * [49] A. Terras.  Harmonic Analysis on Symmetric Spaces and Applications II.  Springer New York, 1988.
  * [50] L. Fejes Tóth.  Über die dichteste Kugellagerung.  Math. Z., 48(676–684), 1943.
  * [51] B.W. Towse, C. Barry, D. Bush, and N. Burgess.  Optimal configurations of spatial scale for grid cell firing under noise and uncertainty.  Phil. Trans. R. Soc. B, 369:20130290, 2014.
  * [52] M. Viazovska.  The sphere packing problem in dimension 8.  Ann. of Math., 185(3):991–1015, 2017.
  * [53] X.-X Wei, J. Prentice, and V. Balasubramanian.  A principle of economy predicts the functional architecture of grid cells.  Elife, 4(e08362), 2015.
  * [54] M.M. Yartsev and N. Ulanovsky.  Representation of Three-Dimensional Space in the Hippocampus of Flying Bats.  Science, 340(6130):367–372, 2013.