Laying the Foundations for Riemann Matrices in Sage

I recently uploaded some aforementioned code for computing Riemann theta functions to the Sage Trac server for review. (Trac #6371). It’s the result of about two months of work and even though there can be some improvements here and there it’s in a stable state. I’m hoping that my advisor, Bernard Deconinck, will take a moment to review it since, after all, he’s the one who wanted Riemann theta in Sage.

The next step of my programming project is to start implementing the components necessary for computing Riemann matrices and the miscellaneous algebro-geometric routines included in Maple’s algcurves package. (e.g. Computing the genus of an arbitrary plane algebraic curve.) What’s nice about this stage in the project is that the necessary components that need to be implemented in Sage are, more or less, independent of each other. This is ripe for collaborative work in a future Sage Days.

Since there are several components, I decided to make a “component map”. Over the next few months I’ll be crossing out these entries piece by piece.

(Above image subject to change.)

Riemann Theta Speed Tests

Riemann theta is in a stable position at the moment and with symbolics now enabled (to be discussed in a future post) it’s time to fine tune some of the steps. At some point I’d like to revisit the computation of the major semiaxis / radius, which involves a root-finding algorithm and evaluating the incomplete gamma function; as well as cacheing of the integer points by implementing the uniform approximation formulas. But the most obvious place in the algorithm where we can get some speed improvement is the actual calculation of the finite sum approximation of Riemann theta. It all boils down to computing thousands to tens of thousands of quadratic forms.

I used Cython to convert all of the involved vectors and matrices into C double complex data types and used CBLAS to evaluate the matrix-vector products and inner products.

First of all, note that the “naive” Sage implementation (Blue) is already a factor of 6-10 faster than Maple. (Red) With the Cython optimization and use of CBLAS the speed up from Maple increases to a respectable factor of 8-24 with an average of 12. (Purple) Further improvements can be still made by exploiting the fact that the imaginary part of the Riemann matrix \Omega is positive definite and symmetric.

Once uniform approximation formulas are implemented in Sage’s Riemann theta I would like to then to a time test on evaluating theta across many data points.

NSF Research Proposal

[Every year the National Science Foundation offers a graduate research grant for first and second year graduate students. The grant provides three years of research funding and is one of the most popular and sought-after grants. One of the components of the application is to write a two-page research proposal that outlines the motivation and goals for my Ph.D research project. Since people are sometimes curious about what it is I exactly do in grad school, I though that I'd share the text of my research proposal here on my blog. I removed the references to papers and journal articles. If anyone would like to see those references I'll be happy to post them here. -- Chris]

With the recent advances in computational tools, experimental and computational mathematics are rapidly growing and increasingly useful fields of study. Their importance is demonstrated by the recent NSF award to Brown University to form the Institute for Computational and Experimental Research in Mathematics. My primary goal, with the guidance of my advisor, Bernard Deconinck, is to provide the infrastructure to make Abelian functions as computationally accessible as trigonometric and hyperbolic functions thus allowing experimental advances in non-linear waves, combinatorial optimization, complex analysis, number theory, and algebra.

Proposed Research and Intellectual Merit

Our first goal is to implement algorithms for computing Abelian functions. Abelian functions can be expressed in terms of homogenous rational functions of the Riemann theta function and its derivatives. Our focus is to compute Riemann theta for arbitrary algebraic curves, that is with no restrictions to the hyperelliptic or nonsingular cases.

An implementation of Riemann theta and related algorithms already exists in Maple via the algcurves package with the underlying algorithms developed by Bobenko, Deconinck, Heil, van Hoeij, Patterson and Schmies. My plan is to develop an implementation of this code in Sage, an open-source mathematics software system. There are several reasons for writing this port. Most importantly, it is an opportunity to familiarize myself with the theory and computations of Riemann theta functions. To quote William Thurston on the use of computation as a tool in advancing understanding in mathematics,

“The standard of correctness and completeness necessary to get a computer program to work at all is a couple of orders of magnitude higher than the mathematical community’s standard of valid proofs.”

Additionally, the authors of the Maple implementation of Riemann theta and related functions would like to see a port to an open-source platform since the open-source environment is conducive to collaborative, community-based development. This also allows the developers to better upkeep the code – upkeep that is currently lacking in the Maple implementation. Finally, this project provides an opportunity to update the original algorithms with more efficient ones, such as Couveignes’ topological approach to computing the monodromy group of a complex algebraic curve. The aim is to complete this phase of the project within the first year of funding.

Our second goal is to design and implement algorithms for computing Fay’s prime form, integrating differentials of the second and third kind on Riemann surfaces, and solving integrable equations. The study of the prime form is motivated by a current collaboration with Bernd Sturmfels, Cynthia Vinzant, and Daniel Plaumann on the applications of Riemann theta functions in computing bitangent lines of quartics and determinantal representations of so-called Helton-Vinnikov curves. Improving techniques for integrating differentials will improve existing algorithms such as the computation of the Abel map. As for integrable equations, we will specifically focus on finite-genus solutions to the Korteweg-de Vries, Kadomtsev-Petviashvili, and Non-linear Schrodinger equations and compare the results with experimental data such as that obtained from the Army Corps of Engineering Station at Duck, NC.

A more utopian, longer term goal is to produce a constructive solution to the Schottky problem. That is, given a Riemann matrix can we determine if it originates from a compact Riemann surface and, if it does, compute which algebraic curve in fact generates this Riemann matrix, up to the usual equivalences? Shiota proved that the Kadomtsev-Petviashvili equation characterized all such Riemann matrices. However, no constructive algorithm to recover the Riemann surface is known.

Broader Impacts

There are many different international communities interested in the study of Abelian functions and, in particular, Riemann theta functions. I was provided full funding to present on computing bitangent lines using Riemann theta functions at the ICMS “Higher-genus sigma function and applications” conference in Edinburgh, UK which was organized by mathematicians and scientists interested in the computational aspects of special functions. I was also awarded funding to present in January 2011 at the Mittag-Leffler “Algebraic Geometry with a view towards applications” conference in Oslo which is mainly organized by algebraic geometers. The interest in computing Riemann theta functions by these two distinct communities demonstrates the wide range of fields impacted by the study of this topic.

In research areas such as non-linear waves and combinatorial optimization the work produced as a result of this fellowship will allow scientists and engineers to compare data with models, make predictions, and visualize Abelian functions; the latter of which is key to educating researchers and the public on the theory and applications of mathematics. The aim is that computing with Abelian functions will one day be as accessible as the standard tools scientists and engineers use today to conduct computational and experimental research.

Follow

Get every new post delivered to your Inbox.