Deep Learning Gauss-Manin Connections

The Gauss-Manin connection of a family of hypersurfaces governs the change of the period matrix along the family. This connection can be complicated even when the equations defining the family look simple. When this is the case, it is computationally expensive to compute the period matrices of varieties in the family via homotopy continuation. We train neural networks that can quickly and reliably guess the complexity of the Gauss-Manin connection of a pencil of hypersurfaces. As an application, we compute the periods of 96% of smooth quartic surfaces in projective 3-space whose defining equation is a sum of five monomials; from the periods of these quartic surfaces, we extract their Picard numbers and the endomorphism fields of their transcendental lattices.


Introduction
There are two ways to study deformations of algebraic varieties. One is purely algebraic, through the explicit polynomial equations of the family. The other is transcendental, through the variation of Hodge structures. The translation of the algebraic to the transcendental is achieved through the Gauss-Manin connection associated to the family. The differential equations governing the flat sections of the connection are those that trace out the variation of Hodge structures in the corresponding flag variety.
Experimentation suggests that algebraic deformations defined by simple equations can give rise to devilishly difficult differential equations that are well beyond our ability to integrate. In a limited capacity, this article will be concerned with the following two questions: Why are some Gauss-Manin connections so complicated? How can we pick deformations with milder connections?
We approach these questions from a practical point of view. The software presented in [Ser19] allows us to compute the variation of Hodge structures of a pencil of hypersurfaces. In this article, the measure of the complexity of the Gauss-Manin connection for a pencil is the amount of time it takes to integrate its flat sections.
Our investigation led us to exploratory data analysis. We trained deep neural networks that perform well in guessing the complexity of Gauss-Manin connections of pencils of hypersurfaces, see Sections 5.8 and 6.4. Their performance can be further improved by giving a local snapshot of the Gauss-Manin connection at a couple of points as in Section 2.5.2.
Our goal in computing the variation of Hodge structures is to carry the Hodge structure (i.e. periods) of one variety onto another. We determine the periods of 96% of smooth quartic surfaces in P 3 that can be expressed as the sum of five monomials, each with coefficient 1. This allows us to determine their Picard number, with a small chance of error [LS19].
In order to carry out a computation at this scale, we use the software in [Ser19] as an underlying engine and add a parallelization layer, data caching mechanisms, checkpointing, and fault-tolerance. Our neural network training software is also available in this code base. We make our code available for general use 1 .
1.1. Pencils of hypersurfaces. The Hodge structure on a smooth hypersurface Y = Z(g) ⊂ P n+1 C of degree d can be represented by a matrix of periods P g ⊂ C m×m , where m depends only on (n, d), as in Section 2.3. The method given in [Ser19] of computing P g involves deforming g to another smooth hypersurface X = Z(f ) ⊂ P n+1 C of degree d whose periods are already known. To begin, one may take X to be a Fermat type hypersurface whose periods can be expressed by closed formulas.
Given such a pair (f, g) we will consider the pencil of hypersurfaces defined by (1 − t)f + tg, which deforms X to Y . Explicitly representing the variation of Hodge structures from X to Y as in Section 2.5, one can compute (i.e. numerically approximate) a matrix P f,g ∈ C m×m such that if P f is a period matrix of X then P f,g · P f is a period matrix of Y . We call such P f,g a period transition matrix. If Z(h) ⊂ P n+1 C is another smooth hypersurface of degree d. The product of the period transition matrices P f,h and P h,g gives a period transition matrix from X to Y .
There is a large variation on the time to compute P f,g in the inputs f, g. The computation of P f,g is often time consuming, taking hours or days, but for some inputs the computation of P f,g could take only a few seconds. A critical observation is that it is sometimes faster to compute P f,h P h,g than it is to compute P f,g directly. This suggests searching for a sequence of polynomials f =s 0 , s 1 , . . . , s k =g for which P si,si+1 is easy to compute for each i. A random search based on simple heuristics was employed in [Ser19, §3.1] to find such sequences.
Unfortunately, it is difficult to predict whether the computation of each P si,si+1 will terminate within k seconds without actually running the computation for k seconds. Instead, we wish to anticipate the difficulty of such a computation so that we may discard difficult pencils in favor of friendlier ones. In this article, we approach this prediction problem with deep learning.
We cast the problem of discovering a good sequence into one of finding a short path in a weighted graph. Let W be a finite set of homogeneous polynomials, all of the same degree, containing f and g. Consider the complete graph G with vertex set W and some weight function ϕ defined on the edges. For an edge e of G, one may define ϕ(e) to be the number of seconds it takes to compute the transition matrix P e . Our ultimate goal is to identify a path in G that connects f and g and has small total weight. (We elaborate on how to choose W in Section 3.1.3.) If the weight function ϕ was known, then finding such an optimal path could be solved using standard graph traversal strategies such as Dijkstra's algorithm or the A * -algorithm [CGR96,RN16]. The problem we face is that the cost of evaluating ϕ at an edge e is just as expensive as computing P e itself. To address this, we enlist data-driven learning models to help us guess whether ϕ(e) is reasonably small or not.
We train our models on a random subset E of the edge set of the graph G. We collect information on E by attempting a computation on each edge e in E that is a representative fragment of the computation needed to evaluate P e . The models learn to recognize if, given an edge e of G, the computation e → P e will terminate in a reasonable amount of time. Informed by the predictions of these models we then traverse a path in G from f to g. See Section 3 for more details on the general method and Section 5 for our implementation of the models.
Unfortunately, the edges e for which P e can be readily computed are generally rare. Consequently, when k is a reasonably small threshold, the spanning subgraph of G whose edge set is {e ∈ E(G) | ϕ(e) ≤ k} is often sparsely connected or disconnected. One drawback of our method is that there may not exist any good path from f to g inside G, in which case we will have wasted time trying to discover one. Enlarging the polynomial set W or the threshold k may solve this problem, but at the cost of a longer training time. Further improvements may require theoretical advances into the nature of the Gauss-Manin connection. Nevertheless, what we give here significantly improves the computation time spent searching for a good connection when one exists, see Section 6.4. 1.2. Deep learning in algebraic geometry. Deep learning can be used to find elliptic fibrations [HL19], to recognize isomorphism classes of groups and rings [HK19], and for approximating the solutions of high-dimensional partial differential equations [SS18]. Although there are symbolic algorithms for these tasks, they are impractical. Deep learning methods are employed to boost performance. The trade off for this gain in performance is the unreliability of the output.
However, predicting the solution to a problem is not the only way one can apply deep learning methods to mathematical computations. Instead, deep learning methods can be used to assist a more reliable method by providing dynamically generated heuristics, thereby improving performance while preserving reliabilitysee for instance [HEW + 19]. This is our approach here, see Remark 4.1.
1.3. Outline. In Section 2 we give an overview of the period computation strategy for hypersurfaces. In Section 3 we explain the problem from a computational point of view and describe the resource management strategy we employ. Section 4 gives an overview of deep learning methodology for the non-specialist. In Section 5 we describe our implementation of deep learning methods. In Section 6 we apply our code to five-monomial quartics and list their Picard numbers as well as their isomorphism classes. us to use their computational resources: Dartmouth College, Harvard University, Leibniz University Hannover, Max Planck Institute (MPI) MiS Leipzig. This project began while all three authors were at MPI MiS, we thank this institution for providing a stimulating environment. In addition, we thank Pierre Lairez for helpful comments.
Avinash Kulkarni has been partially supported by the Simons Collaboration on Arithmetic Geometry, Number Theory, and Computation (Simons Foundation grant 550033) and by the Forschungsinitiative on Symbolic Tools at TU Kaiserslautern during the course of this project.

Period computation
is a degree d homogeneous polynomial. By the Lefschetz hyperplane theorem all cohomology groups of X are trivial (either 0 or Z), except for the middle (singular) cohomology group H n (X, Z).
We need to represent two kinds of structure on the cohomology groups: the integral structure and the Hodge decomposition on H n (X, C). We recall their definition and summarize their method of computation here. For an exhaustive account, see [Voi07a,Voi07b]. The definitions that are more specific to this article are introduced in Section 2.5.
2.1. The integral structure. As an abstract group, H n (X, Z) is isomorphic to Z m for some m. After choosing an identification ψ : H n (X, Z) ∼ → Z m the intersection product on H n (X, Z) can be represented by an m × m matrix I with integral entries. An integral structure on H n (X, C) = H n (X, Z) ⊗ Z C refers to an identification of the sublattice H n (X, Z) ⊂ H n (X, C).
Although we will suppress this from notation, whenever we refer to a trivialization ψ of the integral cohomology, we also mean a determination of the intersection product I on Z m . The integral structure on H n (X, C) can be represented by a map ψ C := ψ ⊗ C.
2.2. The Hodge structure on cohomology. The Hodge decomposition on H n (X, C) = H n (X, Z) ⊗ Z C C m is a direct sum decomposition: where H p,q (X) is the space of (p, q)-forms.
The Hodge pieces H p,q (X) do not vary holomorphically in X. Therefore, it is more natural for variational problems to consider the Hodge filtration F (X) = n p≥ H p,n−p for = 0, . . . , n. Of course, on an individual hypersurface, one can recover the decomposition from the filtration and vice versa.
Using a generic hyperplane section of X, we can define the hyperplane class h ∈ H 2 (X, Z) in cohomology. If n = dim X is even then h n/2 ∈ H n (X, Z) is called a polarization. In this case, the primitive part of the cohomology H n (X, K) 0 is the orthogonal complement of h n/2 , where K is any ring. If n is odd then we set H n (X, K) = H n (X, K) 0 . Let m 0 := dim C H n (X, C) 0 . The restrictions F (X) ∩ H n (X, C) 0 of the Hodge filtrations to the primitive cohomology will be denoted by F (X) 0 .
2.3. The period matrix. Both the integral structure and the Hodge structure are discrete invariants in isolation. After all, we know H n (X, Z) Z m and we know the dimensions of the pieces of the Hodge decomposition (e.g. [Ara12, §17.3]). The difficulty is putting these pieces together, which requires transcendental invariants -the periods of X.
Definition 2.1. Let us call P ∈ C m×m a period matrix on X if there is an isomorphism ψ : H n (X, Z) Z m such that for each = 0, . . . , n, the first dim C F (X) rows of P span ψ C F (X) ⊂ C m . Similarly, we define a primitive period matrix P 0 ∈ C m0×m0 .
Remark 2.2. It is clear that a period matrix P can be obtained from the primitive period matrix P 0 and vice versa. See [LS19,§7] for an explicit computation.
2.4. The Griffiths basis for cohomology. Different identifications of H n (X, Z) with Z m will change P by the action of the discrete group GL(m, Z). However, there is a continuous family of choices to be made in choosing a basis for F (X). For a hypersurface X, one can resolve the latter indeterminacy by specifying a well-defined basis for the filtration.
Let S = C[x 0 , . . . , x n+1 ], Jac(f X ) = (∂ 0 f X , . . . , ∂ n+1 f X ) be the Jacobian ideal and R = S/ Jac(f X ). Since f X is smooth, R is a finite dimensional algebra over C. Let us write R for the quotient of the homogeneous part S / Jac(f X ) .
For each ≥ 0, Griffiths [Gri69] defines a residue map: where vol P n+1 C is a natural generator of the twisted canonical bundle of P n+1 , namely Ω n+1 P n+1 /C (n + 2), and is given by These residue maps descend to an isomorphism on the quotients: Impose the grevlex ordering on S and consider the ideal of leading terms lt(Jac(f X )) of Jac(f X ). The grevlex ordering gives a well-defined sequence of monomials which descend to a basis of S/ lt(Jac(f X )), and therefore to a basis of R. Their residues ω 1 , . . . , ω m in appropriate degrees yields a basis of the primitive cohomology H n (X, C) 0 . Definition 2.3. The basis ω 1 , . . . , ω m of the primitive cohomology constructed above is a well-defined basis which respects the filtration. We will call this basis the (grevlex) Griffiths basis on X.
Definition 2.4. A primitive period matrix P 0 as in Definition 2.1 will be called a primitive grevlex period matrix of X = Z(f X ) if the i-th row of P 0 equals j=1 is the grevlex Griffiths basis on X. Any period matrix obtained by extending P 0 as in Remark 2.2 will be called a grevlex period matrix of X.
Any two grevlex period matrices of X differ by the action of the discrete group GL(m, Z). Moreover, after fixing a coordinate system on the integral primitive cohomology, the grevlex period matrix is uniquely defined.
2.5. Transition matrices for periods. Consider a pair of smooth hypersur- Definition 2.5. We will call a matrix P f,g ∈ C m0×m0 a period transition matrix if for any grevlex primitive period matrix P f on X = Z(f ), the product P f,g · P f is a primitive grevlex period matrix on Y = Z(g).
A method for computing a period transition matrix P f,g is explained in [Ser19]. A much simpler method is available when g is a linear translate of f and we will do this example in Section 2.5.1. We will briefly outline the general method of [Ser19] here.
(3) Find the matrix B expressing the change of basis from the basis above to the grevlex basis at t = 1. This is done by the computation of normal forms.
(4) For each i = 1, . . . , m find a differential operator D i ∈ C(t)[ ∂ ∂t ] such that D i · res(p i ) = 0. These differential equations annihilate the rows of the grevlex period matrix of X t .
(5) Find a path γ in C \ S where S is the set of values of t for which X t is singular. (6) Compute the transition matrix for the space of solutions of D i at t = 0 and t = 1 obtained by homotopy continuation along γ. (7) Using the indicial equation of D i , find the set of derivatives of res(p i ) whose values at t = 0 would determine res(p i ). (8) Express the derivatives of res(p i ) at t = 0 in terms of the grevlex Griffiths basis at t = 0. Multiply this expression on the right with the transition matrix of (6) and take the first row. (9) Form the m × m matrix whose i-th row is the row obtained at the step above and multiply it by B on the left to get P 0,1 .
The most time expensive step is step (6). However, most attempts to compute the period transition matrix fail on step (4). This is because it requires a Gröbner basis computation for Jac(f t ), numerous normal form computations, and expressions of elements in Jac(f t ) in terms of the generators (∂ 0 f t , . . . , ∂ n+1 f t ).
2.5.1. Period transition matrices of linear translates. If two hypersurfaces are linear translates of one another, then a period transition matrix between them can be easily computed. We will do this here as an instructive example. This is also a computation we will use later in the article in Section 6.
The general linear group GL(n + 1, C) acts on the coordinates on P n+1 linearly. The induced action on the coordinate ring S = C[x 0 , . . . , x n+1 ] is given by u · f (x) = f (x · u t ) where u ∈ GL(n + 1, C), f ∈ S, u t is the transpose of the matrix u, and x = (x 0 , . . . , x n+1 ) is treated as a row vector.
With X = Z(f ) a smooth hypersurface as before, suppose that Y = Z(g) with g = u · f for some u ∈ GL(n + 1, C). If P X is a period matrix of X then it certainly works as a period matrix for Y . However, even if P X is a grevlex period matrix on X it need not be grevlex on Y . We describe the computation of a period transition matrix P f,g below. Let p 1 , . . . , p m ∈ S be the polynomials whose residues give the grevlex Griffiths basis on X. Let φ : Y ∼ → X be the isomorphism induced by u. If P X is a primitive grevlex period matrix on X then P X is a primitive period matrix on Y whose rows represent the residues of u · p 1 , . . . , u · p m .
We can use Griffiths-Dwork reduction (see [Gri69] or a summary [Ser19, §2.4]) on Y to write each res(u · p i ) in terms of the grevlex Griffith basis on Y . If N is a matrix whose rows store the coordinates of res(u · p i ) in the grevlex basis then N −1 · P X will be a primitive grevlex period matrix on Y . In other words, N −1 is a period translation matrix from X to Y . This computation is implemented in PeriodSuite 2 as the function translate_period_matrix.
2.5.2. First order Gauss-Manin connection. We would like to detect which of the differential operators D i appearing in Item (4) of Section 2.5 would be easy to integrate, before we compute D i .
Collectively, these D i define flat sections of the Gauss-Manin connection of the family X t (see for example [Voi07a, §9.3]). Computing the Gauss-Manin connection itself is simpler than computing its flat sections, but still not quite fast enough for rapid testing. Computing the Gauss-Manning connection evaluated at a single point, however, is very fast. Furthermore, these evaluations give an impression of how complicated D i might be (see Section 5.7). Let be a basis for the primitive cohomology on a hypersurface X. Say X is the fiber of a family X t at t = t 0 . Then the residue of the polynomials p i on the family X t will give a basis for the primitive cohomology for all t near t = t 0 . Differentiating these forms with respect to t and evaluating at t = t 0 gives new elements in the cohomology of X. Expressing these new forms in gives an m × m matrix. We call this matrix the first order Gauss-Manin matrix.
Let us note that if X and the polynomials p i are defined over a subfield K ⊂ C then the corresponding first order Gauss-Manin matrix will have entries in K.
See [KO68] for more on this topic. In our applications, we will take K = Q.

Computational scheme
We will now describe our approach to the problem of searching for a good path between polynomials, with the goal of transferring their periods from one to the other. There are two variations of this problem that we are interested in solving. We will state and explain these variations and then generalize them to a common framework in Section 3.1. In Section 3.2 we will explain how we operate in this abstract framework.
3.1. Two types of problems and a general framework. We are interested in solving two types of problems. In the first problem, we are given a pair of polynomials V = {f, g} depicting smooth hypersurfaces of the same degree, and our goal is to compute the period transition matrix P f,g . In the second problem, we are given a (possibly large) set of polynomials V , and our goal is to compute the periods of all elements in V given the periods of any one of them. Let us present some example problems we will solve later on.
3.1.1. First problem: computing periods for one target polynomial. For a general (f, g), a direct computation of P f,g is infeasible even with low precision. The path from f to g must be broken into simpler pieces. The strategy we discussed in Section 1.1 is to find a sequence of polynomials f = h 0 , h 1 , . . . , h s = g such that each intermediate period transition matrix P hi,hi+1 is easily computable. Their product would then give P f,g . This problem was already investigated in [Ser19, §3.1] and a crude heuristic was developed there. We will develop this heuristic further in this section.
3.1.2. Second problem: computing periods for many hypersurfaces. Here we are given a (possibly large) set of polynomials V and we would like to be able to compute the periods of all elements in V given the periods of any one of them. For example, the set V may be the set V n consisting of all smooth quaternary quartics that are expressed as the sum of n distinct monomials with coefficients equal to 1, e.g.
An attempt to compute the periods, and therefore the Picard numbers, of all elements in V 5 was made in [LS19]. Many of the elements in V 5 were out of reach at the time. We apply our methods to compute the periods of most elements in V 5 , see Section 6.
Given a sequence of polynomials f = h 0 , h 1 , . . . , h s = g such that each intermediate period transition matrix P hi,hi+1 is easily computable, and such that the periods of f are known to some precision, we easily obtain the periods of each intermediate quartic by the partial products P hi,hi+1 . . . P h0,h1 P f . Thus, it suffices to determine a set of paths which connect the vertices of V as opposed to computing the periods one by one.
3.1.3. General framework. Both problems can be slightly generalized to fit into the following framework. Suppose we are given a set V containing one element whose periods are known and the periods of all the others are sought.
The set V may not have enough pairs f, g ∈ V such that P f,g is directly computable. In this case, we need to construct a larger set W containing V which introduces many pairs (f, g) ∈ V × W such that P f,g is directly computable. The construction of this W can be based on human heuristics.
Letting K W be the complete graph with vertex set W , we now wish to solve the following problem: Given V ⊂ W , find a tree T ⊂ K W such that the vertex set of T contains V and the computation of the period transition matrix for each edge in T is feasible. Constructing W is a balancing act. If the search space W is too large, it may be impractical to find a good tree T inside K W , even though one may exist. Conversely, if W is too small, we may be able to search the entire space but find that there is no good tree T inside K W .
For our first type of problem (Section 3.1.1), with V = {f, g}, one may take the following set: where supp denotes the monomial support of a polynomial and coefs denotes the set of coefficients of a polynomial. For the second type of problem we will consider V = V 5 . In this case, it seems natural to take something in between Remark 3.1. There are numerous modifications to this problem that can help attain the final goal. We will only point out one: it is often sensible to prune K W by a problem specific heuristic before embarking on a search for T .
3.2. Searching a computation graph for an efficient tree. In this section, we no longer discuss period computations. We are interested only in abstracting the problem laid out in Section 3.1.3 to, essentially, a problem in resource management in a computational exploration.
Let G = (W, E) be a graph with vertex set W and edge set E. Consider a computationally expensive program P which takes as input an edge e ∈ E and returns an output P e . We will consider a cost function ϕ : E → R + that measures the difficulty of performing the computation e → P e . Common examples of ϕ(e) would be the number of arithmetic operations needed to complete the program e → P e or the amount of memory required to execute the program. For our purposes, it is important that a lower bound for ϕ can be determined in finite time. Henceforth, we fix one such ϕ and refer to it as a complexity function.
Remark 3.2. In the context of this paper, W is a set of polynomials. For each edge e, our procedure outputs a period transition matrix P e to a fixed degree of precision. The value of ϕ(e) may then be the number of seconds it took to perform the computation e → P e on a fixed computer. These specifics will not be relevant for the rest of this section, as we continue with the abstraction above.
One could rephrase the problem introduced in the previous section as follows: Given V ⊂ W , find a tree T ⊂ G such that V is contained in the vertex set of T and ϕ(T ) = e∈T ϕ(e) is small, if not minimal.
If we know the value of ϕ on each edge, then this problem becomes a standard minimization problem. However, at the outset, the evaluation of ϕ is just as expensive as the computation e → P e itself. This leads to the following brute-force strategy.
Algorithm 1 Brute force with thresholding Attempt the computation n → P Q(n) .

4:
if successful then 5: return G 8: Increment n 9: return Fail {(One may apply multiple rounds with enlarged W and k to eventually succeed.)} 3.2.1. Brute-force strategy. The brute force method (Algorithm 1) attempts to compute every edge in the graph until the desired tree is constructed. If the computation of an edge takes too long, then the computation is aborted and the computation for the next edge is begun. In practice, we perform the computations in the job queue in parallel, but the size of the edge set (≈ 10 6 ) is so much larger than the number of cores available to us (≈ 10 2 ) that this serial conceptualization is not too far off the mark.
3.2.2. Informed brute force. It is clear that the choice of the ordering Q on the edge set will make a dramatic impact on the total time it takes to find T . In ideal cases, the first few edges in the sequence might be easily computable and sufficient to form a good tree T . We could then stop searching early. In unfortunate cases, all of the edges e for which ϕ(e) > k might be queued up first in the list, in which case we would end up spending all our time trying to compute impossible edges. For this reason, we set finding a good ordering Q as our top priority.
A strategic choice of a statistical model will be used to give us a function φ : E → [0, 1] such that if the computability score φ(e) is close to 1 then ϕ(e) is likely at most k. Moreover, as we perform computations, we can refine our model in order to improve the reliability of this function φ. This suggests the modification of the brute force algorithm given by Algorithm 2. Collect data during this computation in order to retrain the models and improve the reliability of φ(e). 7: if result is not Fail then 8: return result 9: Possibly enlarge W or increase k. 10: Retrain φ. 11: goto 1, using the improved φ.
Along with classical regression models, we use deep neural networks for regression and classification tasks. In the next two sections we will provide a brief introduction to neural networks, in which we define neural networks and outline standard training procedures. We will then describe our implementation of such deep learning models (i.e. choosing hyperparameters) for the computational problems at hand.

Deep learning models
In this section we provide a quick overview of neural networks for the mathematician who does not specialize in applied or computational mathematics. A good reference for this subject is [GBC16].
Section 4.1 establishes an analogy between neural networks and polynomial interpolation. Neural networks are defined in Section 4.2 together with the class of functions associated to them. In Section 4.3, we outline the method of gradient descent on this class of functions. 4.1. Neural networks to approximate functions. A neural network is essentially a framework for approximating functions. Let ϕ : R n → R m be an unknown "target" function that we wish to approximate, and let Γ ϕ ⊂ R n+m be the graph of ϕ. Suppose that we have access to a finite subset T sampled from Γ ϕ . Assuming strong hypotheses on T and ϕ, we can in some cases closely approximate or even recover ϕ.
A classical example of such a scenario is univariate polynomial interpolation. If the target function ϕ is known to be a polynomial with degree bounded above by some known k, then we can closely approximate ϕ from T (provided T is sampled sufficiently well) by solving a reasonably-sized linear system. This example demonstrates a good use-case for polynomial interpolation. However, the linear system could be defined by a Vandermonde matrix that is extremely ill-conditioned, making a solution very difficult and sensitive to noise.
Thus interpolation may not always be the most appropriate way to learn ϕ. In many real-world scenarios (e.g. that encountered in this work), we cannot assume that ϕ is polynomial. A function ϕ that is not sufficiently regular may require interpolating polynomials of very high degree. This can result in a highdimensional linear system requiring a large matrix inversion, which can be a computationally prohibitive task. In this work, we relax our requirement of a perfect fit between model and data, and use a neural network architecture with a variation of gradient descent to solve such a regression problem.
In analysis, it is common practice to construct a class of functions C which exhibits so-called good approximating properties. Once an appropriate class has been chosen, we try to find a sequence of functions {ϕ n } ∞ n=1 ⊂ C which converges in some sense to ϕ. As will be discussed below, a chosen neural network architecture determines such a class 3 C and provides an algorithm for constructing an approximating sequence {ϕ n } ∞ n=1 using T ⊂ Γ ϕ that will hopefully converge to ϕ. This is more generally referred to in computer science as regression.
Remark 4.1. We seek to approximate a complexity function (as in Section 3.2) associated to the period computations, rather than approximating the period function itself. This is because our intended use for the periods require hundreds of digits of accuracy and, ideally, rigorous bounds on error. Approximating functions with neural networks typically capture large scale features of a function. Attaining high precision, let alone bounding the error, is not one of the highest priorities of the method. Furthermore, we can better tolerate error in the complexity function as we will use it only as a heuristic to order our computations (Section 3.2.2).

4.2.
The class of functions associated to a neural network. We define a neural network as a triplet N = (V, E, A) where V = (V 0 , . . . , V k+1 ) is a sequence of real vector spaces with a fixed coordinate system By an affine transformation we mean the composition of a linear map The non-linear endomorphisms E i are typically chosen to be of a very specific form in order to facilitate computations and to enable gradient propagation. They do not all need to be the same function. Unless specified otherwise, we will always take each E i to be a rectified linear unit, also denoted ReLU. When the architecture (V, E) is fixed, we suppress it from notation and associate to each neural network (V, E, A) the function f A : V 0 → V k+1 defined as follows: for some fixed k, referred to as a hyperparameter of the network. The name hyperparameter refers to certain characteristics (e.g. the depth k of the network and the width dim(V i ) of each layer) that define the class C of functions allowed by the neural network. We will call each A i • E i a layer of the neural network N. We refer to layer i = 0 (resp. i = k) as the input (resp. output) layers; the remaining layers are called hidden layers. The naming suggests the layers' use and the asymmetry apparent in the construction of f A . The affine transformations A i are parametrized by variables called weights and variables called biases. In this setting we will equate a neural network simply with a composition of such layers.
It is remarkable that something as simple as the incorporation of ReLU functions vastly expands the space of functions that can be represented by neural networks. Without the non-linearity of the endomorphisms E, the function f A would simply be an affine transformation.
4.3. Gradient descent using neural networks. Fixing the architecture (V, E) of a neural network yields a parametrized family of neural networks (V, E, A), which can be associated to the family of functions C = {f A | A}. For each i > 0, the affine transformation A i has n i n i+1 + n i+1 trainable parameters, so the total number of parameters in f A is That is, we have a parametrization Our immediate goal is to find a sequence A (k) ∈ R N , k ≥ 0, such that the sequence of functions {f A (k) } ∞ k=0 converges in a sense to an approximation of our target function ϕ : V 0 → V k+1 . We will first describe the distance measure that defines convergence for our experiments.
4.3.1. Loss function. Given a function g : V 0 → V k+1 we wish to quantify the goodness of our network's function approximation, i.e. how far g is from being equal to ϕ. The only information we are given about ϕ is the finite subset T ⊂ Γ ϕ ⊂ V 0 × V k+1 of its graph; this subset is called the network's training set. As we have a fixed coordinate system on V k+1 , we will use the induced Euclidean be the loss function defined as follows: As in the classical philosophy of regression, we say that to fit ϕ well given T is to minimize L with respect to A. In the following we discuss two variations of a popular iterative method with the aim of achieving the minimization arg min A∈R N

L(A; T).
The following class of algorithms provides a sequence {f A (k) } ∞ k=0 whose Lvalues will hopefully (and in some limited cases, provably) be decreasing. 4.3.2. Gradient descent. Given T ⊂ Γ ϕ we wish to find A ∈ R N minimizing the error function L(·, T). Our restriction to feedforward networks forbids feedback loops and, therefore, allows for an easy evaluation of the gradient ∇L(·, T) of L(·, T) at any given point A via backpropagation -see [GBC16] for definitions and details.
Choose an initial point A (0) ∈ R N and a sequence of step sizes γ : N → R >0 . This sequence is typically either constant, or converging to zero. Inductively define the following sequence: In the applications we have in mind, the size of T will be too large to make the execution of this method feasible. In examples where #T is large, one might opt for some variant of stochastic gradient descent; one such method is described below. 4.3.3. Stochastic minibatch gradient descent. Stochastic minibatch gradient descent differs from gradient descent in that one trains on a random subset of T at each step, instead of T itself. To do this, fix a batch size b ∈ N, and define A (k) inductively as follows: Choose a random subset T k ⊂ T of size b and let (4.3.8) 4.3.4. Hyperparameter selection. In order to use a neural network architecture (V, E) to approximate the function ϕ, for which a subset T ⊂ Γ ϕ is known, one requires the following: An error function L, step sizes γ, batch sizes b, distributions to randomly choose subsets of T and to choose a starting point A (0) . Collectively, these choices (V, E, L, γ, b) are called hyperparameters. This defines the class C as described in Section 4.1.
Selecting hyperparameters can be more of an art than a science. Although principled selection strategies have been proposed, e.g. [CD15] provides an excellent survey of some such methods, this problem remains largely unsolved for general learning problems. The rate of convergence of an iterative algorithm, and even whether the algorithm converges at all, can depend heavily on parameter choice. For example, a step size chosen too small will cause the algorithm to crawl slowly to a local minimum, whereas a step size that is too large might cause the algorithm to diverge. Experimentation is required to select hyperparameters in a way so that the stochastic gradient converges rapidly and to a reasonable approximation of ϕ.
In the next section, we will discuss our implementation of such a learning algorithm, including how we chose the hyperparameters for the problem of approximating a complexity function as in Section 3.2.

Implementation
Now that we have reviewed the fundamentals of neural networks, we return to the problems explained in Section 3.2. We propose data-driven approaches to learning a proxy to the complexity function ϕ of that section.
The design space for our experiments includes a choice of labelled training dataset T, a statistical model (classical or neural), and hyperparameters that define that model. Sections 5.1-5.3 are concerned with the design of the dataset on which we will learn. Section 5.4 compares useful intermediate representations of this dataset, a subdiscipline called feature extraction. Section 5.8 discusses the experiments with 4-and 5-monomial quartic surfaces that motivated our choice of hyperparameters. 5.1. Computability score. The goal of the learning task in this section is to obtain a function φ : E → [0, 1] which assigns to each edge e ∈ E a probability that the computation e → P e of the period transition matrix will terminate on our hardware. For any such φ we will refer to φ(e) as a computability score of e. The idea is to choose φ so that φ(e) is easy to evaluate but φ(e) is indicative of the computability of our entire algorithm (e → P e ) on e.
Remark 5.1. One could, in principle, set φ = 1 1+ϕ where ϕ is the complexity function introduced in Section 3.2. However, the difficulty in explicitly computing ϕ motivates us to find a more practical, data-driven solution.
Assuming the kind of function φ we seek is moderately well behaved, we could approximate it via statistical learning methods using only finitely many pairs (e, φ(e)). We will choose a random subset E ⊂ E and assign a value φ(e) for each e ∈ E . The resulting set of pairs (5.1.9) T = {(e, φ(e)) | e ∈ E } will be used to train a statistical model to obtain a modest guess for what φ should be. The remaining pairs E\E will be used for testing and validation. We will next define the E and φ that we have chosen for our application to the problems introduced in Section 3.1.

5.2.
Experimental input space E. The abstract discussion in this chapter applies in the full generality of Section 3.1. Recall the notation of Section 3.1.2, where V k ⊂ Q[x, y, z, w] 4 denotes the k-nomial data set, that is, the set of fourvariable homogeneous polynomials that are the sum of k distinct monomials all with coefficient 1. For this section, we will constrain ourselves to fixed sets of polynomials -i.e. to the complete graphs on V 4 and V 5 . We denote by E the edges of the graph V k that is eventually to be traversed. For instance, the 4-nomial data set is defined by the complete graph on V 4 and thus #E = V4 2 . 5.3. Experimental computability score φ. In this section we define our choice of computability score φ. Experimentally we find that one critical subroutine of our larger algorithm e → P e tends to present a computational bottleneck. The global computation of P e requires computing as many ODEs as columns of P e and then numerically integrating them. We will let φ represent how long it takes to compute only the first of these ODEs on a fixed computer. That is, φ : E → R ≥0 is defined so that ε(e) is the computation time for the first ODE required for P e .
Remark 5.2. We experimented with including additional label information in training, such as the degree and order of the first ODE. However, we did not observe an advantage to using this additional data in training our neural networks. This may be attributed to the visible relation between time, degree and order in Figure 1.

5.3.1.
A practical complexity measure. We generate many training samples of φ, evaluated on edges of the graphs based on V 4 and V 5 . Since φ(e) can be arbitrarily large, we choose a threshold and terminate the computation of φ(e) after 30 seconds. In that case, we know φ(e) ≥ 30. To capture this process we define a binary label β(e) ∈ {0, 1}, where β(e) = 0 if φ(e) < 30 (successful computation) and β(e) = 1 if φ(e) ≥ 30 (failed computation).
We will refer to the resulting edge-time correspondences {(e, β(e))} for V k as k-nomial data sets. These pairs will either be used for training or testing, in which case we may refer to e as an input data sample, and to β(e) as that sample's label. Our architecture choice is in part guided by the structure of this data. In our learning task, it ends up being easier to learn the binary classification β rather than the nuanced function φ itself; in fact, in learning we obtain a proxy for φ given by the per-class probability distribution on E. Thus we replace our regression problem (approximate φ) with a binary classification problem (approximate β, and get a proxy for φ en route). We refer to these two functions interchangeably in the development of our strategy. We must be careful with how we represent an edge (f, g) (i.e, a pair of quartic surfaces). Finding a good representation is called feature extraction and can be seen as enlarging the input space V 0 in a meaningful way. As the following section details, we are able to compute extra data associated to (f, g) that empirically relates to the complexity ϕ(f, g).

Feature extraction.
We approximate the computability score function φ : E → [0, 1] via a statistical model that takes as input a vector representation of e ∈ E. We represent each edge e by the concatenation of the coefficient vectors of the two polynomials that are the endpoints of that edge. It must be noted that in using this vector representation, we have discarded some information: the value of φ is linked inextricably to the fact that endpoints of e are polynomials, a characteristic of which the model is no longer aware. This guides the idea that our chosen learning method can better guess the value of φ(e) if it is provided more than just the coordinate representation of e.
One additional piece of domain information that we found useful was the firstorder Gauss-Manin connection, see Section 2.5.2. We can efficiently compute the first-order Gauss-Manin connection at a few points of the pencil corresponding to e. In practice, the complexity of these matrices correlates with the computability score of e, as shown in Sections 5.5.1 and 5.7.
We will use the following notation for these matrices. First, fix t 1 , . . . , t s ∈ C.
For each e let M e = (M e,1 , . . . , M e,s ) be a sequence of such matrices, where M e,i is the first order Gauss-Manin connection evaluated at t = t i on the pencil corresponding to e. We will design our models to work with (e, M e ). In practice, we take s = 1 or s = 2 with t 1 = 0 and t 2 = 1 because we observed little benefit in increasing s further. 5.5. Studying the complexity of the input. When the polynomials used in E have rational coefficients, the connection matrices M e will also be rationalvalued. For example, when we are dealing with quartic surfaces defined over Q, we have that (e, M e ) is an element of Q n , where n = 2 × 35 + s × 21 × 21. In this section, we will consider various measures of complexity of the matrix M e . Any correlation between "easy" statistics (defined in the next section) of M e , and the algorithm-fragment runtime φ(e) is not obvious. This motivates the use of more complex data-driven models such as neural networks. 5.5.1. Complexity of cohomology matrices. Because the standard implementation of neural networks work with floating point arithmetic, the subtleties of computing with a rational number are lost. As a large portion of the computation e → P e is exact, the computation time is affected by the "height" of the rational numbers involved. For this reason, we will modify the entries of M e to better represent the complexity of its entries.
Let us define the following function on rational numbers: (5.5.10) where m 1 , m 2 ∈ Z, m 2 > 0 and lcm(m 1 , m 2 ) = 1. The value ψ(m) of a rational number m is a more faithful representation of the complexity of computing with m then would be a floating point approximation of m. The following variation will also be used: (5.5.11) ψ entropy : m 1 m 2 → log(|m 1 |) 2 + log(m 2 ) 2 .
Various complexity statistics can be extracted from a complexity matrix M . Each statistic is a function Ψ(M ) on a 3-tensor M = (M ijk ) of rational numbers, such as M e . We list three options on Table 2. Each column of Figure 3 corresponds to an entry in Table 2. For each Ψ, we plot Ψ(M e ) against the time it takes to compute the algorithm fragment defining φ, for those edges such that β(e) = 1. The edges for which β(e) = 0 are omitted from the top row of Figure 3 because we had to terminate their computation prematurely, effectively assigning them all the same value of 30. The bottom row of Figure 3 shows that successful edges do tend to have lower matrix statistics than failing edges. However, these distributions are not bimodal enough to make the statistics good classifiers in isolation.
On the other hand, we see no striking patterns between φ(e) and Ψ(M e ) in the first row of Figure 3. From this we conclude that the matrices M e are useful in terms of classifying successes from failures (i.e. approximating β), but their statistics alone are not sufficient to regress (i.e. approximate φ) within the class of successes.
We can gain more insight by viewing each tensor M e as a multi-channel image. We can visualize a matrix (M ij ) as a rectangular image with the ij-th entry colored a shade of blue: the darker the shading, the larger the value of ψ(M ij ). See Figure 4, which has M e,1 in the first row and M e,2 in the second row for four e from the 4-monomial data set. We will use convolutional neural networks to better convey the spatial relationship of entries in a matrix to the neural network. 5.6. Preprocessing the data. 5.6.1. Dimension reduction: principal component analysis. In the edge sets E we consider, the pair of polynomials in e for each e ∈ E are sparse. As a result, we are using an unnecessarily large ambient space to represent e. This in turn, dilutes the capacity of a neural network to learn. Instead we perform a principal component analysis (PCA) on a given E to compress the data stored in e by projecting onto a subspace with minimal loss in information. When considering quartic surfaces, E lies in a 70 = 35+35 dimensional ambient space. However, projecting onto a well chosen smaller subspace cuts down the dimension with minimal loss of information. For example, Figure 5 shows that the first 23 principal components almost entirely recovers the coordinates of edges in the complete graph over V 4 . This number 23 works well also on V 5 . We thus compress each e to a 23 dimensional vector, down from 70. 5.6.2. Balancing the dataset. Unfortunately, even the algorithm fragment of Section 5.3.1 fails to terminate in the vast majority of pairs we considered. As a result, a training set T as in (5.1.9) will consist almost entirely of failed edges of the form (e, 0). This incentivizes the neural network to produce the constant 0 function as an approximation of φ.
As a remedy, we employ a standard method to disincentivize the neural network from converging to a constant. This method is done by over-sampling the under-represented class (here, edges e with β(e) = 1) so that the training set consists of an equal portion of both classes. 5.7. Learning the computability score. A crucial decision is whether to use neural networks at all. Classical statistical methods are simpler to operate and when they are more successful, they can give additional insight about the dataset. However, neural networks have a broader applicability and they can succeed when the classical methods fail. In Section 5.7.1 we show that we are in the latter case, the neural networks consistently outperform the classical methods.
In this subsection we will work with the 4-monomial dataset, i.e. the complete graph on V 4 . We tried several classification techniques as candidates for "learning" φ : E → [0, 1]. It turns out that we get a good estimate for φ if we just try to approximate the classification function β : E → {0, 1}. The methods we consider output functions f : E → [0, 1] and we must couple them with a cut-off value τ so that β(e) = 0 is most likely correlated with f (e) ≤ τ .
In order to rank different methods, we drew in Figure 6 the receiver operating characteristic curves of ten different methods. The curve for each method is obtained by varying the cutoff value τ , with true positive rate and true negative rate as the axes. The closer a curve is to the top left vertex, the better the corresponding method performs. The dotted line is the idealized curve for the method of random guessing.
The best-performing method (deep ensemble) is a composition of the next best two neural network strategies, namely multilayer perceptron and convolutional NN. The fourth place is occupied by a method that is related to both neural networks and classical statistical methods (Gaussian kernel SVM), this is reflected here by its performance. The remaining six are classical methods, all of which underperform in this task. In Section 5.7.1 we describe these methods in greater detail. 5.7.1. Classical statistical methods. We tried seven standard binary classifiers to predict φ from the 4-nomial training set T. Since they were clearly outperformed by the neural networks, we will not discuss them in depth. The methods we used are: logistic regression with an L 2 penalty, a regularized support vector classifier with a linear kernel and regularization parameter 1, a regularized support vector classifier with degree-2 radial basis function kernel and regularization parameter 1, a random forest with 10 trees, per-tree maximum depth 5, an AdaBoost classifier with at most 50 estimators, linear discriminant analysis, and Gaussian Naive Bayes. Figure 6 suggests that the deep classifiers will outperform the classical methods. This may simply be caused by the fact that deep learning methods have more internal parameters and can, therefore, provide better approximations. 5.7.2. Deep neural networks. We will introduce two neural networks and then combine their results into an ensemble method. The first one, a neural network, follows the standard formulation in Section 4; which is a multi-layer perceptron (MLP). In dealing with the quartic fewnomial dataset, we decided on an architecture with five hidden layers, each of width 100. As input, it takes only the edges e ∈ E after compression via principal component analysis as in Section 5.6.1.
Its output is a single value in R. To restrict the codomain to the interval [0, 1] we apply the inverse logit function. After training, the neural network gives an approximation φ MLP : E → [0, 1] of the computability score.  Our second neural network is a two-channel convolutional neural network, a variation of the standard neural network explained in Section 4 to better detect patterns in visual data. This neural network will be trained using the 3-tensors M e encoding the first order Gauss-Manin connections. We will however apply the complexity function ψ from Section 5.5.1 to each entry of M e before giving it as input. The output is adjusted as with the first neural network (MLP) above so We improve on the approximations of the two neural networks by defining the function (5.7.12) This is what we call the "deep ensemble" method (or just ensemble method), as illustrated in Figure 8. Being the product of two functions, which are essentially probability functions, the function φ ensemble is more cautious in returning a value close to 1. We chose this approach because attempting the computation e → P e for edges that do not terminate can be very costly. We prefer a computability score that has a low false positive rate.
5.8. Performance on datasets. We now describe our experimental setup and the performance of our neural networks. Recall that V n is the set of smooth quartics that are the sum of n monomials. = 5778. These numbers are small enough that we could run our test computation (as in Section 5.3.1) on all edges in E. This gives us complete information and allows us to evaluate the performance of our neural network.
We train the ensemble neural network defined in (5.7.12) on E and test it on E := E \ E . We wish to know whether the neural network correctly predicts if the edges in E can be successfully integrated in a short amount of time. We list in Table 9 Table 9. The performance of the ensemble neural network on V 4 with varying proportion α of edges used for training.
Taking α = 0.9, in Figure 7 we plotted the ROC curves of the multilayer perceptron component of the ensemble for different widths of its first three layers. We decided on a width of 500 based on this figure. 5.8.2. Four monomial quartics. The performance of the ensemble method trained solely on V 4 did not perform well on V 5 . One plausible reason is that the first order Gauss-Manin connections on V 4 are all quite simple, whereas on V 5 there is a broader range of complexity displayed by these matrices. In other words, the structure in V 4 does not extrapolate well to V 5 .
To address this, we retrained the ensemble method on 5-monomials. Consider the edges E of the complete graph with vertices V 5 . This is a much larger set, since 3348 = #V 5 #V 4 = 108, so #E = 5, 602, 878. With the ratio α defined as in Section 5.8.1, we display the performance for various α in Table 10 Table 10. The performance of the ensemble neural network on V 5 .

Application
In this section we will give an application of our software package and analyze the performance improvement of using neural networks.
6.1. Five monomials. One application of our software is to the set V 5 of smooth polynomials, each of which are the sum of five monomial terms, with coefficients 1, see Section 3.1.2 for this notation. Our goal is to compute the Picard numbers of these 5-nomial quartics and to find unexpected isomorphisms in this set.
6.2. List of results. Note that the results of this section depend on finding integral linear relations between periods that are only known approximately. Therefore, the results below may contain errors due to insufficient precision. We tried to mitigate this risk by working with 300 digits of precision, when we could attain it.
Of 3348 smooth 5-nomial quartics, there are only 161 S 4 -symmetry classes. Out of this 161 we could reach 154. We computed their Picard numbers and we list their frequency below. Table 11. Picard number frequency for 5-nomial quartics ρ 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 frequency 8 6 0 32 1 9 2 13 0 22 0 9 1 13 4 7 1 20 6 0 In contrast to the 180, 000+ quartics in the database given in [LS19], we are essentially after 161 quartics. However, the database in loc. cit. was built by random exploration to find easy to compute quartics. Here, our quartics V 5 are fixed. For this reason, we were forced to use an order of magnitude more computation time to build our graph (a CPU year versus decade). In fact, some of these quartics are so difficult to reach that our "optimal" tree T spanning 154 S 4 -classes has a diameter of length 21.
Remark 6.1. The Picard number for every 5-nomial quartic in this database was verified using a complementary method that uses the crystalline cohomology of finite reductions [CS20]. This method produces guaranteed upper bounds.
For each of our polynomials, we ran through as many prime reductions as was necessary to have the minimum of the upper bounds attained thrice. In each case, this minimum agreed with the Picard numbers we computed. Remark 6.2. As we compute the Picard lattice and not just the Picard number, we can use [LS19, §3] to find classes of smooth rational curves in the Picard group. We observed that only 7 of the 154 quartics had Picard groups that could not be generated over Q with the polarization and classes of smooth rational curves of degree ≤ 3.
6.2.1. The missing vertices. While our search method succeeded for the vast majority of quartics, we were unable to reach seven of the 5-nomial quartics in V 5 , up to isomorphism. These seven quartics are listed in Table 12. Moreover, we list the Picard number bound obtained by crystalline_obstruction [CS20] as explained in Remark 6.1.
We also tried to brute force every edge from V 4 ∪ V 5 to establish a connection to the three quartics in this table with Picard number 2. Allowing three hours per connection to find only the first ODE succeeded only in connecting these vertices to one another. Naturally we could not try every edge with the three hour time limit. Therefore, it is conceivable but very unlikely that good connections exist.
This prompts us to ask: Which feature of these quartics is responsible for making them inaccessible? Picard number x 3 y + y 3 z + y 3 w + z 3 w + xw 3 ≤ 2 xy 3 + z 4 + x 3 w + y 2 zw + xw 3 ≤ 2 x 4 + y 3 z + xz 3 + x 3 w + yw 3 ≤ 2 y 3 z + xyz 2 + xz 3 + x 3 w + yw 3 ≤ 3 x 3 y + y 3 z + z 3 w + z 2 w 2 + xw 3 ≤ 3 x 2 y 2 + x 3 z + yz 3 + y 3 w + xw 3 ≤ 18 xy 3 + x 3 z + xyzw + z 3 w + yw 3 ≤ 19 6.2.2. Isomorphism classes. Using the Torelli theorem for K3 surfaces [LP80] we can check if the K3 surfaces in our list admit non-trivial isomorphisms. Here, we work only with 154 (of 161) representatives of the S 4 -symmetry classes for which we could compute the periods. We compared their period vectors modulo an integral change of basis for homology. The method of computation is described in Section 6.3.2. We found 9 isomorphism classes of size 2 and 2 isomorphism classes of size 4, all other isomorphism classes appear to be of size 1. We display the non-trivial isomorphism classes in Table 13.
6.3.1. Additional simplifications. The symmetric group S 4 acts on the 5-nomials by permuting the four variables. This is a linear action of the projective space and we can use Section 2.5.1 to compute period translation matrices at essentially no cost. This connects the elements in each S 4 -equivalence class.
In order to translate the periods of one polynomial p to another q, we need the period matrix of p. However, if p is particularly resistant to our computations then we can compute only the first row of the period matrix of p -reducing the work load by a factor of 21. With this first row we are still able to compute the Picard number and isomorphism class of p. However, p becomes a dead-end; we can no longer use p to compute the periods of another polynomial q.
6.3.2. Computing isomorphism classes. The isomorphism class of a K3 surface depends only on its periods [LP80]. In particular, that of the first row of its period matrix.
Suppose w 1 , w 2 ∈ C 22 H 2 (X, C) are periods of two K3s X 1 and X 2 . To detect if X 1 and X 2 are isomorphic, we need to determine if there exists a constant c ∈ C * and an isometry N ∈ Z 22×22 , N : H 2 (X 1 , Z) ∼ → H 2 (X 2 , Z), such that (6.3.13) w 1 · N = cw 2 .
Using approximations of w 1 and w 2 , this can be translated into a problem of finding short lattice vectors as we describe below. The integral relations annihilating w 1 and w 2 cause a difficulty here. So we first compute the Picard groups Pic(X i ) w ⊥ i ⊂ Z 22 . If the rank of the Picard groups are distinct then X 1 and X 2 are not isomorphic.
If ρ := rk Pic(X 1 ) = rk Pic(X 2 ), construct T (X i ) = Pic(X i ) ⊥ ⊂ Z 22 . We can view w i as an element in T (X i ) ⊗ C C 22−ρ . Let v i ∈ C 22−ρ be the new vector corresponding to w i . The surfaces X 1 and X 2 are isomorphic if and only if there exists c ∈ C * and N ∈ Z (22−ρ)×(22−ρ) that satisfy (6.3.14) v 1 · N = cv 2 .
Following [LS19, §2.3], we describe how to find such an N if it exists. That is, if ·, · represents the intersection product on T (X 2 ) then we need to solve for (6.3.15) v 2 , v 1 · N v 2 = v 2 , v 2 v 1 · N , which is linear in N . We now use LLL [LLL82] to find possible integral solutions N . The function isomorphisms_of_k3s in PeriodSuite implements this procedure.
6.4. Performance on applications. In the end, our goal is to improve the computation time for periods. In Section 5.8 we analyzed the predictive power of the neural network in isolation, but we did not consider its effect on period computation. We describe the influence of the neural network on the task of computing the periods here.
We demonstrate the effect of using a neural network on a small example. The problem here was faced on a larger scale, and faced repeatedly, as we sought to complete the calculations for Table 11. Given a 5-nomial f ∈ V 5 , we look for 4-or 6-nomials g ∈ V 4 ∪ V 6 such that the period transition matrix for the edge (f, g) is easy to compute. We used such connections to zig-zag from V 5 to V 4 or V 6 and back to V 5 in order to establish new connections between 5-nomials.
For this example, we chose a random subset S ⊂ V 5 with 100 elements so that each element in S is from one of the 161 distinct S 4 -equivalence classes. For each f ∈ S we consider the edge set (6.4.16) The average size of E f for f ∈ S is 29. We compare two methods of exploring the edges E f , one aided by our neural network and one unaided. We used the neural network (ensemble) that was trained on V 5 but not on V 6 so that there is no extra training time. In particular, the neural network is faced with a data set for which it has not been trained ; nevertheless it performs well.
For the unaided strategy, we picked 10 elements from each E f randomly and tried to compute these edges. For the aided strategy, we sorted E f using our neural network and picked the top 10. The unaided strategy had a 54.4% failure rate as opposed to 33.9% for the aided strategy. Consider the table below that records the frequency of elements in S that had n successful edges for n ∈ {0, . . . , 10}. Table 15. Frequency of successes # of connections 0 1 2 3 4 5 6 7 8 9 10 Unaided method 21 3 5 9 10 16 8 7 3 5 13 Aided method 21 0 0 4 4 5 6 4 7 2 47 The first column of this table shows that, in both cases, 21 vertices had 0 successful edges. This demonstrates the fact that some edges are intrinsically difficult to move away from. Neural networks can sort edges according to difficulty, but they cannot help if every edge is impossible.
On the other hand, we see that the aided method establishes far more connections to V 6 . For instance, 47 vertices in S had all 10 of their chosen edges successful with the aided method as opposed to 13 with the unaided method. In practice, this computation would then be repeated for each successful connection, which means that the advantage grows exponentially.
Our computations for the Picard ranks of V 5 took a CPU decade. The approach presented in this paper allowed us to repeatedly pare down hundreds of thousands of possible edges to a manageable, but likely to succeed, subset. The mini-computation in this section demonstrates the benefit of including a neural ensemble model in the algorithmic pipeline en route to computing period matrices of smooth quartic hypersurfaces.