\documentclass[a4paper, 11pt]{amsart}

\newfont{\cyr}{wncyr10 scaled 1100}

\usepackage[left=2.7cm,right=2.7cm,top=3.5cm,bottom=3cm]{geometry}
\usepackage{amsthm,amssymb,amsmath,amsfonts,mathrsfs,amscd}
\usepackage[latin1]{inputenc}
\usepackage[all]{xy}
\usepackage{latexsym} 
\usepackage{url}
\usepackage{longtable}
\usepackage{threeparttable}
\usepackage{algorithmic}
\usepackage{array}

\theoremstyle{plain}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{corollary}[theorem]{Corollary}
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{conjecture}[theorem]{Conjecture}

\theoremstyle{definition}
\newtheorem{definition}[theorem]{Definition}
\newtheorem{question}[theorem]{Question}
\newtheorem{examplewr}[theorem]{Example}
\newtheorem{ass}[theorem]{Assumption}

\theoremstyle{remark}
\newtheorem{obswr}[theorem]{Observation}
\newtheorem{remarkwr}[theorem]{Remark}

\newenvironment{obs}{\begin{obswr}\begin{upshape}}{\end{upshape}\end{obswr}}
\newenvironment{remark}{\begin{remarkwr}\begin{upshape}}{\end{upshape}\end{remarkwr}}
\newenvironment{example}{\begin{examplewr}\begin{upshape}}{\end{upshape}\end{examplewr}}


\DeclareMathOperator{\dR}{dR}
\DeclareMathOperator{\sym}{Sym}
\DeclareMathOperator{\Sym}{Sym}
\DeclareMathOperator{\sgn}{sgn}
\DeclareMathOperator{\spec}{spec}
\DeclareMathOperator{\la}{la}
\DeclareMathOperator{\ab}{ab}
\DeclareMathOperator{\cx}{cx}


\DeclareMathOperator{\II}{II}
\DeclareMathOperator{\Betti}{B}





\newcommand{\cF}{\mathcal F}
\newcommand{\cH}{\mathcal H}
\newcommand{\cC}{\mathcal C}
\newcommand{\V}{\mathcal V}
\newcommand{\E}{\mathcal E}
\newcommand{\cE}{\mathcal E}
\newcommand{\g}{\gamma}
\newcommand{\G}{\Gamma}
\newcommand{\Gap}{\Gamma_0^D(p M)}
\newcommand{\Ga}{\Gamma_0^D(M)}
\newcommand{\hGa}{\hat{\Gamma }_0(M)}
\newcommand{\Q}{\mathbb{Q}}
\newcommand{\Z}{\mathbb{Z}}
\newcommand{\F}{\mathbb{F}}
\newcommand{\K}{\mathbb{K}}
\newcommand{\C}{\mathbb{C}}
\newcommand{\I}{\mathbb{I}}
\newcommand{\Ql}{\Q_{\ell}}
\newcommand{\GQ}{G_{\Q}}
\newcommand{\PP}{\mathbb{P}}
\newcommand{\Sel}{\mathrm{Sel}}
\newcommand{\Sha}{\mathrm{Sha}}
\newcommand{\Gal}{\mathrm{Gal\,}}
\newcommand{\GL}{\mathrm{GL}}
\newcommand{\Cl}{\mathrm{Cl}}
\newcommand{\Tr}{\mathrm{Tr}}
\newcommand{\Div}{\mathrm{Div}}
\newcommand{\Fil}{\mathrm{Fil}}
\newcommand{\cond}{\mathrm{cond}}
\newcommand{\Frob}{\mathrm{Frob}}
\newcommand{\End}{\mathrm{End}}
\newcommand{\Jac}{\mathrm{Jac}}
\newcommand{\disc}{\mathrm{disc}}
\newcommand{\Id}{\mathrm{Id}}
\newcommand{\MS}{\mathrm{MS}}
\newcommand{\Aut}{\mathrm{Aut}}
\newcommand{\Fr}{\mathrm{Fr}}
\newcommand{\diag}{\mathrm{diag}}
\newcommand{\Up}{\mathrm{U}_p}
\newcommand{\ord}{{\mathrm{ord}}}
\newfont{\gotip}{eufb10 at 12pt}
\newcommand{\gn}{\mbox{\gotip n}}
\newcommand{\gm}{\mbox{\gotip m}}
\newcommand{\ga}{\mbox{\gotip a}}
\newcommand{\gb}{\mbox{\gotip b}}
\newcommand{\gf}{\mbox{\gotip f}}
\newcommand{\go}{{\mathcal O}_{A_f}}
\newcommand{\OAf}{{\mathcal O}_{A_f}}
\newcommand{\gp}{\mbox{\gotip p}}
\newcommand{\gL}{\mbox{\gotip L}}
\newcommand{\gd}{\mbox{\gotip d}}
\newcommand{\cO}{{\mathcal O}}
\newcommand{\cM}{{\mathcal M}}
\newcommand{\cP}{{\mathcal P}}
\newcommand{\cL}{{\mathcal L}}
\newcommand{\cG}{{\mathcal G}}
\newcommand{\GG}{{\mathcal G}}
\newcommand{\cD}{{\mathcal D}}
\newcommand{\om}{{\omega}}
\newcommand{\ra}{\rightarrow}
\newcommand{\lra}{\longrightarrow}
\newcommand{\hra}{\hookrightarrow}
\newcommand{\qbar}{\bar\Q}
\newcommand{\Norm}{\mathrm{N}}
\newcommand{\NS}{\mathrm{NS}}
\newcommand{\SL}{{\mathrm {SL}}}
\newcommand{\SO}{{\mathrm {SO}}}
\newcommand{\n}{{\mathrm{n}}}
\newcommand{\Pic}{{\mathrm{Pic}}}
\newcommand{\CH}{{\mathrm{CH}}}
\newcommand{\AJ}{{\mathrm{AJ}}}
\newcommand{\tr}{{\mathrm{tr}}}
\newcommand{\si}{\sigma}
\newcommand{\eps}{\epsilon}
\newcommand{\ta}{{\tau}}
\newcommand{\R}{{\mathbb R}}
\newcommand{\M}{{\mathrm{M}}}
\newcommand{\HC}{\mathcal F_{\rm har}}
\newcommand{\PGL}{{\mathrm{PGL}}}
\newcommand{\Supp}{\mathrm{Supp}}
\newcommand{\ep}{\varepsilon}

\def \mint {\times \hskip -1.1em \int}

%\newcommand{\mint}{\times\!\!\!\!\!\hskip-0.3em\int}

\newcommand{\SP}{\vspace{0.5cm}}
\newcommand{\T}{\mathbb T}
\newcommand{\new}{\mathrm{new}}
\DeclareMathOperator{\Hom}{Hom} \DeclareMathOperator{\Meas}{Meas}
\newcommand{\X}{\mathbb X}
\newcommand{\Y}{\mathbb Y}
\newcommand{\U}{\mathbb U}
\newcommand{\Rat}{{\rm Rat}\bigl(\PP^1(\Q_p),\C_p\bigr)}
\newcommand{\Ker}{{\mathrm Ker}}
\newcommand{\res}{\mathrm{res}}

\newcommand{\rom}{\mathrm}
\newcommand{\fr}{\mathfrak}
\newcommand{\cl }{\mathrm{cl}}

\newcommand{\longmono}{\mbox{$\lhook\joinrel\longrightarrow$}}
\newcommand{\longleftmono}{\mbox{$\longleftarrow\joinrel\rhook$}}
\newcommand{\longepi}{\mbox{$\relbar\joinrel\twoheadrightarrow$}}
\newcommand{\mat}[4]{\left(\begin{array}{cc}#1&#2\\#3&#4\end{array}\right)}
\newcommand{\smallmat}[4]{\bigl(\begin{smallmatrix}#1&#2\\#3&#4\end{smallmatrix}\bigr)}
\newcommand{\dirlim}{\mathop{\varinjlim}\limits}
\newcommand{\invlim}{\mathop{\varprojlim}\limits}

\newcommand{\D}{\mathbb D}
\DeclareMathOperator{\Ext}{Ext}
\DeclareMathOperator{\Rep}{Rep}
\DeclareMathOperator{\mhs}{MHS}
\DeclareMathOperator{\et}{et}
\newcommand{\bbH}{\mathbb{H}}
\DeclareMathOperator{\mm}{MM}
\numberwithin{equation}{subsection}
\numberwithin{theorem}{subsection}


\include{thebibliography}

\begin{document}

\title[Computation of triple Chow-Heegner points]{The effective computation of iterated integrals and Chow-Heegner points on triple products}

\author{Henri Darmon, Michael Daub, Sam Lichtenstein and Victor Rotger}
%\contrib[with an appendix by]{William Stein}

\begin{abstract}
FIXME
\end{abstract}

\address{H. D.: Montreal, Canada}
\email{darmon@math.mcgill.ca}
\address{M. D.: Berkeley, California}
\email{mwdaub@math.berkeley.edu}
\address{S. L.: Stanford, California}
\email{saml@math.stanford.edu}
\address{V. R.: Departament de Matem\`{a}tica Aplicada II, Universitat Polit\`{e}cnica de Catalunya, C. Jordi Girona 1-3, 08034 Barcelona, Spain}
\email{victor.rotger@upc.edu}



\maketitle

\section*{Introduction}

The goal of this paper is to describe a complex-analytic algorithm for the computation of \textit{triple Chow-Heegner points}. Fix cuspidal eigenforms $f,g\in S_2(\Gamma_0(N))$ and assume $f$ is a newform with rational Fourier coefficients. To $g$ is associated a Hecke correspondence $T_g\subset X_0(N)\times X_0(N)$, which gives rise to a rational point $P_g\in J_0(N)(\mathbf Q)$. (A more precise definition of $P_g$ is given below.)
The triple Chow-Heegner point $P_{g,f}$ associated to the 3-tuple of modular forms $(g,g,f)$ is the image of $P_g$ in $A_f$, the elliptic curve quotient of $J_0(N)$ associated to the newform $f$. 

It is shown in \cite{DRS} that the rational point $P_{g,f}$ can be computed as an element of the analytic curve $A_f(\mathbf C)$ in terms of a certain iterated path integral (in the sense of \cite{Chen}). This formula is amenable to numerical computation, which we have implemented using the free software package {\tt SAGE}.
%\footnote{Source code is available at <insert URL>.} 
The results of \cite{DRS} together with work of Yuan-Zhang-Zhang give a criterion (cf. \cite[Theorem 1]{DRS}) for when the points $P_{g,f}$ are non-torsion. This criterion implies that triple Chow-Heegner points comprise a collection of non-torsion points on many elliptic curves $A_f$ of rank 1. Our algorithm makes these points readily computable in practice for many  elliptic curves $A_f$ of small conductor.

While the analytic formula is not the only way of computing the points $P_{g,f}$ (see the Appendix) our approach has a theoretical advantage: it requires knowing only the Hodge class $\xi_g$ associated to the cycle $T_g$. In future work, we hope to adapt the algorithm described below to compute Chow-Heegner points \cite{DRS2} associated to Hodge classes on ``modular varieties'' (related to Kuga-Sato varieties), such as classes $\xi$ arising from modular forms with complex multiplication. The rationality of Chow-Heegner points computed in this manner could provide numerical evidence for certain open cases of the Hodge conjecture.

In \S  1 we recall necessary facts about itegrated integrals and related ingredients for our main algorithm.
In \S  2 we specialize to the case of modular curves, define the points $P_{g,f}$ precisely, and write down an explicit analytic formula for them. In \S  3 we describe in detail an algorithm for evaluating this formula numerically. The algorithm is illustrated with numerical examples in \S  4. Some tables of triple Chow-Heegner points on elliptic curves of small conductor are presented in \S  5, along with some discussion of a few phenomena apparent from this data. 

% One interesting question is to determine the index $[A_f(\mathbf Q): \mathbf Z P_{g,f}]$ in terms of arithmetic invariants attached to the elliptic curve $A_f$. By BSD, this index should be related to the ``fudge factor'' connecting the derivative $L(f\otimes g\otimes g, 2)$ with the height of $P_{g,f}$; cf. \cite[Theorem 4.1]{DRS} and Yuan-Zhang-Zhang FIXME:add ref. (What is the fudge factor? Do our data match is?) The numerical data in our tables could be useful for investigating this matter.
\section{Preliminaries}
\subsection{}
Let $F$ be a number field (we take $F = \mathbf Q$ in the sequel)  and fix an embedding $\iota:F\hookrightarrow \mathbf C$.
% Let $\bar F$ denote the algebraic closure of $F$ in $\C$ and set  $G_F:=\Gal(\bar F/F)$. 
Let $X$ be a smooth, complete algebraic curve of genus $g\geq 2$ 
over $F$, and let $Y=X\setminus \{ \infty\}$ be the complement of a single point in $X(F)$. For a smooth variety $V_{/F}$ (such as $X$ or $Y$) we denote by $V^{\mathrm{an}}$ the complex manifold
$(V\otimes_{F,\iota} \mathbf C)(\mathbf C)$ with its analytic topology.

\subsection{}
The {\em de Rham cohomology}  $H^1_{\dR}(X^{\mathrm{an}},\mathbf C)$ 
is the cohomology of the de Rham complex of smooth $\mathbf C$-valued
 differential forms on $X^{\mathrm{an}}$. 

% There is a canonical identification
% $H^1_{\mathrm{Betti}}(Y^{\mathrm{an}},\mathbf Z) \otimes \mathbf C = H^1_{\mathrm{dR}}(Y^{\mathrm{an}}, \mathbf C)$
% arising from integration of smooth differential forms against smooth chains.

% The vector space $H^1_{\mathrm{dR}}(Y^{\mathrm{an}},\mathbf C)$  is equipped with a
% decreasing filtration $\Fil^\bullet$, the so-called {\em  Hodge filtration}, which is
%  given by
% $$ \Fil^j H^1_{\dR}(Y^{\mathrm{an}}) = \left\{
% \begin{array}{ll}
% H^1_{\dR}(Y^{\mathrm{an}}) & \mbox{ if } j \le 0,\\
% \Omega^1_{\log}(Y^{\mathrm{an}}) & \mbox{ if } j = 1,\\
% 0  & \mbox{ if } j \ge 2,
% \end{array} \right.
% $$
% where $\Omega^1_{\log}(Y^{\mathrm{an}})$ denotes the space of holomorphic differentials on $Y^{\mathrm{an}}$
%  with at worst a
% logarithmic pole at the point $\infty$. 

Because the Riemann surface $X^{\mathrm{an}}$ arises from an algebraic curve over $F$, we can 
identify $H^1_{\dR}(X^{\mathrm{an}},\mathbf C)$ with $H^1_{\dR}(X/F)\otimes \mathbf C$, where 
\begin{equation}\label{deRham-comparison}
H^1_{\dR}(X/F) := \mathbf H^1(0\rightarrow \cO_X \rightarrow \Omega^1_X \rightarrow 0)
\end{equation}
is the {\em algebraic de Rham cohomology} of $X/F$, defined as 
the hypercohomology of the de Rham complex of sheaves of
regular differential forms on $X$.

% The Hodge filtration on $H^1_{\dR}(Y/F)$  is described concretely by
% \begin{equation}
% \label{filtrationH1}
% \Fil^j H^1_{\dR}(Y/F) = \left\{ 
% \begin{array}{ll}
% H^1_{\dR}(Y/F) \simeq \Omega^{II}(Y)/dF(Y) & \mbox{ if } j\le 0, \\
% \Omega^1_{\log}(Y/F) & \mbox{ if } j = 1 \\
% 0 & \mbox{ if } j \ge 2,
% \end{array} \right.
% \end{equation}
% where 
% \begin{enumerate}
% \item
% $\Omega^{II}(Y)$ is the space 
% of  rational differential forms on $Y$ over $F$
%  with vanishing residues at all points of $Y$; 
% these are also called {\em differentials
% of the second kind} on $Y$;
% \item 
% $F(Y)$ is the field of rational functions on $Y$;
% \item 
% $\Omega^1_{\log}(Y/F) = \{ \omega \in \Omega^1(Y) 
% \quad \mbox{ with } \quad  \mathrm{ord}_\infty(\omega)\geq -1\}$, 
% where $\mathrm{ord}_P(\omega)$ 
%  denotes  the valuation of a differential
%  $\omega$ at a closed point $P$. 
% \end{enumerate}

As is well known, the fact that $X$ is a curve means that $H^1_{\dR}(X/F)$ has a particularly simple description in terms of $\Omega^{II}(X)$, the space of \textit{differentials of the second kind} on $X$. By definition, these are rational 1-forms on $X$ with vanishing residues at all points of $X$. 
By the residue formula  we may identify $\Omega^{II}(X)$ with $\Omega^{II}(Y)$, the differentials of the second kind on $Y$. Thus we have a canonical isomorphism
\[H^1_{\dR}(X/F) = \Omega^{II}(Y)/\mathrm{d}F(X),\]
where $F(X)$ is the field of rational functions on $X$.
By applying Riemann-Roch, this description can be simplified: it is not difficult to show that $\Omega^{II}(Y)/\mathrm{d}F(X) \cong \Omega^1(Y)/\mathrm{d}\Gamma(Y,\mathcal O_Y)$. So $H^1_{\dR}(X/F)$ can also be computed as the space of
\textit{regular} 1-forms on $Y$, modulo exact forms. For computational purposes, the latter description is the most useful: we will compute with classes in $H^1(Y)$ using rational 1-forms on $X$, regular away from the point $\infty$.
These are amenable to computation via their Laurent expansions about $\infty$.
\subsection{}

Fix a base point $o\in Y^{\mathrm{an}}$; let $\Gamma := \pi_1(Y^{\mathrm{an}};o)$ denote the fundamental group
of the Riemann surface $Y^{\mathrm{an}}$.
Let $\mathbf Z[\Gamma]$ be the integral group ring on $\Gamma$ and write $I\subset \mathbf Z[\Gamma]$  for its augmentation ideal.
Note that $H_1(X^\mathrm{an},\mathbf Z) = H_1(Y^{\mathrm{an}},\mathbf Z) \cong \Gamma^{\mathrm{ab}}$
(as follows from the well-known presentation for the fundamental group of a Riemann surface), and that this abelian group 
is naturally identified with $I/I^2$.

The {\em path space} on $Y$ based at $o$,
denoted $\mathbf P(Y;o)$,  is  the set of
piecewise-smooth paths 
$$
p:[0,1]\lra Y^{\mathrm{an}}, 
\quad  \mbox{ with  } p(0)=o. $$ 
Let $\pi:{\tilde Y}\to Y^{\mathrm{an}}$ denote the universal covering space of  $Y^{\mathrm{an}}$ corresponding to the choice of basepoint $o$, which
 can be regarded  as the space of homotopy classes in 
$\mathbf P(Y;o)$. Likewise, denote by $\tilde X$ the universal cover of $X$
corresonding to the same basepoint $o$.
%This space  is conformally equivalent to the open unit disc in $\mathbf C$.
 The group $\Gamma $
 acts on $\tilde Y$ 
transitively and without fixed points, 
and the map $p\mapsto p(1)$  identifies
the quotient ${\tilde Y}/\Gamma$ with $Y^{\mathrm{an}}$. 
Recall that if $\eta$ is a closed $C^\infty$ 1-form
(resp. a meromorphic 1-form of the second kind)
 on $X^{\mathrm{an}}$, then it admits a smooth (resp. meromorphic) \textit{primitive function}
 $F_\eta:{\tilde X}\to \mathbf C,$
 defined by the rule
$$F_\eta(p) := \int_0^1 p^\ast \eta.$$


The {\em  basic iterated integral }
attached to a tuple of smooth 1-forms  $\omega_1,\ldots,\omega_n$ on $Y^{\mathrm{an}}$, evaluated along a path 
$p\in \mathbf P(Y;o)$,
is defined to be
\begin{equation*}
\label{eqn:def-ii}
 \int_{p} \omega_1\cdot\omega_2\cdot\ldots \cdot\omega_n := \int_{0\le t_n \le t_{n-1} \le \cdots \le t_1 \le 1}
p^*(\omega_1)(t_1) \cdots p^*(\omega_n)(t_n).
\end{equation*}
The integer $n$ is called the {\em length} of this basic iterated integral.
Note that when $n=2$,
the basic iterated integral attached to
 $\omega$ and $\eta$ can be computed  by the formula
$$ \int_\gamma \omega \cdot \eta =  
\int_\gamma \omega F_{\eta} = 
\int_0^1  \gamma^\ast(\omega F_{\eta}).$$

An {\em iterated integral} 
 is  a linear combination  of basic iterated integrals, perhaps of different lengths, viewed as a function on $\mathbf P(Y;o)$.
 The  length of an iterated
integral is then defined to be the maximum of the lengths
of its constituent basic iterated integrals.

 An iterated integral is said to be 
{\em homotopy invariant} if its value on any path $p$  depends only on the homotopy class of
$p$. The space $\II(Y)$ of homotopy invariant iterated integrals will
 be viewed as a subspace
of the space of $\mathbf C$-valued functions on $\Gamma$.
Extending   $J\in \II(Y)$ to 
the group ring $\mathbf C[\Gamma]$ by
$\mathbf C$-linearity,  we regard $\II(Y)$ as a space of complex functionals on $\mathbf C[\Gamma]$ via the inclusion $  \II(Y) \subset \operatorname{Hom}_{\mathbf C}(\mathbf C[\Gamma],\mathbf C)$.

For each $n$, let
 $\II^{\le n}(Y)$ 
denote the subspace of  homotopy invariant
iterated integrals of length $\le n$. 
Observe that any element $J\in \II^{\le n}(Y)\subset \operatorname{Hom}(\mathbf Z[\Gamma],\mathbf C)$ 
 vanishes on $I^{n+1}$, and hence gives rise to a well-defined  element of
$\operatorname{Hom}(I/I^{n+1},\mathbf C)$. 
The natural map
$ \II^{\le n}\lra \operatorname{Hom}(I/I^{n+1},\mathbf C)$
is an isomorphism.\footnote{FIXME: add reference? This is stated without proof in DRS.}


% The following lemma is helpful in describing $\II^{\le 2}(Y)$ concretely.
% \begin{lemma}
% Let $\eta_1$ and $\eta_2$ be  closed  $C^\infty$ 1-forms  on $Y^{\mathrm{an}}$, and let 
%  $\alpha$
%  be any smooth one-form
%  satisfying the differential equation 
% \begin{equation}
% \label{eqn:diff-hard}
%  d \alpha = \eta_1 \wedge \eta_2.
% \end{equation}
% Then we have the following.
% \begin{itemize}
% \item[(i)] The  iterated integral
% $$ J_{\eta_1\otimes\eta_2,\alpha}(\gamma) = \int_{\gamma} \eta_1 \cdot \eta_2 - \int_{\gamma} \alpha$$
% is homotopy invariant. 

% \item[(ii)] The space $\II^{\le 2}(Y)$ is spanned by
%  iterated integrals of the form $J_{\eta_1\otimes\eta_2,\alpha}$. 
% \end{itemize}
% \end{lemma}
% \begin{proof}
% % The first claim follows directly from the fact that the differential 
% % on ${\tilde Y}$ 
% % appearing in the expression
% % $$ J_{\eta_1\otimes \eta_2,\alpha}(\gamma) = \int_\gamma \eta_1 F_{\eta_2} - 
% % \alpha$$
% % is closed.

% % As for the second, given any $J\in \II^{\le 2}(Y)$ , after  writing $p(J) = \sum \omega_i \otimes \eta_i$,
% % it can be checked that
% % $$ J = \left(\sum J_{\omega_i\otimes \eta_i,\alpha_i}\right)  - J_{0,0,\alpHa}$$
% % for suitable smooth differentials $\alpha_i$ on $Y(\C)$ 
% %   and $\alpha \in H^1_{\dR}(Y(\C))$. 
% See \cite[Prop. 3.1]{Hain2}, for example.
% \end{proof}


We will be interested in numerically evaluating certain iterated integrals $J \in \II^{\leq 2}(Y)$.
 Specificaly, suppose we are given $\omega,\eta\in H^1_{\dR}(X/F)$,
represented as differentials of the second kind, regular on $Y$.  Recall that a differential on a Riemann surface is said to have a
 {\em logarithmic pole} at a point if its expansion in terms of a local parameter $q$
  at this point is of the form $\sum_{n=0}^\infty a_n q^n \frac{dq}{q}$.
Let $\alpha_{\omega,\eta}$ be a meromorphic 1-form on $X$ which is regular on $Y$ 
and is such that the induced differential $\omega F_{\eta}-\alpha$ on $\tilde X$ has at worst a logarithmic pole at (any point lying over) $\infty$. This condition is well-posed because the principal part of $\omega F_\eta$ at at $\tilde x\in \tilde X$
depends only on the image $x$ of $\tilde x$; see \cite[\S 1]{DRS}.
The form $\alpha_{\omega,\eta}$ exists -- and in fact can even be taken to be algebraic and defined over $F$ --  by Riemann-Roch.

\begin{lemma}\label{homotopyinv}The iterated integral $J_{\omega,\eta}:= \int \omega\cdot \eta - \alpha_{\omega,\eta}$, viewed as a function on $\mathbf P(Y,o)$,  is homotopy-invariant.

Moreover, suppose that $\omega$ and $\eta$ represent integral cohomology classes.
Then when $\II^{\leq 2}(Y)$ is identified with 
 $\operatorname{Hom}(I/I^3, \mathbf C)$, the restriction of $J_{\omega,\eta}$
to $I^2/I^3$ is $\mathbf Z$-valued and can be identified with $\omega\otimes \eta$,
viewed as an element of
\[H^1(X,\mathbf Z)\otimes H^1(X,\mathbf Z) \cong (H_1(X,\mathbf Z)\otimes H_1(X,\mathbf Z))^\vee = (I/I^2 \otimes I/I^2)^\vee = (I^2/I^3)^\vee,\]
where $A^\vee$ denotes the $\mathbf Z$-dual of an abelian group $A$.
 \end{lemma}

\begin{proof}
The homotopy invariance of $J_{\omega,\eta}$ follows from the fact that $J_{\omega,\eta}(\gamma) = \int_\gamma \omega F_\eta - \alpha_{\omega,\eta}$,
and the one form on $\tilde X$ in the integrand is holomorphic when restricted to $\tilde Y$. For the second claim, see the discussion at the beginning of \S 1 of \cite{DRS}, and  \textit{loc. cit.}, Lemma 1.1(2). 
\end{proof}


% [FIXME: Is this paragraph necessary?]
% Given  $\omega\in \Omega^1_{\log}(Y^{\mathrm{an}})$ and a
% smooth closed one-form $\eta$ on $Y^{\mathrm{an}}$, consider the expressions 
% $\omega\otimes\eta + \eta\otimes\omega$ and
% $J_{\omega\otimes\eta+\eta\otimes\omega,\alpha}\in \Fil^1 H^1_{\dR}(Y^{\mathrm{an}})$.
% [FIXME: huh?]
% Since
% the differential
% $\omega\wedge\eta+ \eta\wedge\omega$ vanishes identically, the 
% form $\alpha$ of type $(1,0)$ necessarily belongs to 
% $\Omega^1_{\log}(Y(\C))$, and can be chosen to be $0$.
% With this choice, we then observe that 
% \begin{eqnarray*}\label{v}
% \label{square}
% J_{\omega\otimes \eta+\eta\otimes\omega,0}(\gamma) &=& 
% \int_\gamma (\omega \cdot \eta+\eta\cdot \omega) \\ 
% &=& \int_{0\leq s\leq t\leq 1}\gamma^*\omega(s) \gamma^* \eta(t) 
% + \int_{0\leq s\leq t\leq 1}\gamma^*\eta(s) \gamma^* \omega(t)   \\
% &=& \int_{[0,1]\times [0,1]}\gamma^*\omega(s) \gamma^* \eta(t) =  \left(\int_\gamma \omega \right) 
% \left(\int_\gamma\eta\right).
% \end{eqnarray*}
% This element of $\II^{\le 2}(Y)$ can be readily calculated in practice.

%[FIXME: We need to clarify where we are working on $\tilde Y$
%and where we are working on $\tilde X$ = universal cover of $X^{\mathrm{an}}$.]

% Let us represent the cohomology class of a smooth closed 1-form $\eta$ by a 
%  meromorphic differential form of the second kind on 
% $Y^{\mathrm{an}}$,
% making use of the algebraic description of de Rham cohomology.
% Let $\omega$ be a holomorphic 1-form on $X^{\mathrm{an}}$.
% The principal part of the meromorphic $1$-form
%  $\omega  F_{\eta}$ on $\tilde Y$
% at a point $P$ only depends on the image $\pi(P)\in Y^{\mathrm{an}}$; cf. \cite[\S 1]{DRS}.
%  % Indeed, let
% %  $P'=\gamma P\in \tilde Y$, for $\gamma\in \Gamma$, be any other point in 
% % $\pi^{-1}\pi (P)$. Then, comparing   principal parts, we find
% % \begin{equation}\label{principal}
% % \mathrm{pral}_{\gamma P}(\omega F_\eta) = \mathrm{pral}_P(\gamma^*(\omega F_\eta)) 
% % = \mathrm{pral}_P(\omega F_\eta + c\omega)
% % \end{equation}
% % for some constant $c$, and the claim follows because 
% % $\omega$ is regular. 
% In particular, the principal part
% $\underline{\mathrm{pp}}_P(\omega F_\eta)$ at $P\in Y^{\mathrm{an}}$ is well-defined, even though
% the one-form $\omega F_\eta$ is only defined on ${\tilde Y}$. 
% By the  Riemann-Roch theorem,  there exists a meromorphic 
% differential form $\alpha=\alpha_{\omega,\eta}$ on $Y^{\mathrm{an}}$ with at
%  worst a logarithmic pole at $\infty$
%  such that 
% $\underline{\mathrm{pp}}_P(\alpha) = \underline{\mathrm{pp}}_P(\omega\cdot F_{\eta})$
%  for all $P\in X^{\mathrm{an}}$. 
% The one-form $ \omega\cdot F_{\eta}-\alpha$ 
% is then holomorphic on $\tilde Y$. 
% \begin{lemma}
% The iterated integral
% \begin{equation}\label{iv}
% J_{\omega\otimes \eta,\alpha} = \int \omega\cdot \eta - \int \alpha
% \end{equation} 
% belongs to $\II^{\le 2}(Y)$.
% \end{lemma}
% \begin{proof}
% FIXME
% \end{proof}

% , and this
% seems easier  in practice than
% solving the equation
% \eqref{eqn:diff-hard} associated to  a  smooth  representative for the class
% of $\eta$. 

Now consider an integral class $\xi = \sum \omega_i\otimes \eta_i \in H^1(X,\mathbf Z)\otimes H^1(X,\mathbf Z)$.
By the previous lemma, the iterated integral $J_\xi = \sum J_{\omega_i,\eta_i}$
is homotopy invariant and induces a homomorphism
\[J_\xi:H_1(X,\mathbf Z) = I/I^2 \to \mathbf C/\mathbf Z.\]
Fix an auxiliary holomorphic 1-form $\rho \in  H^{1,0}(X_{\mathbf C})\subset H^1(X^\mathrm{an},\mathbf C)$.
Denote by $\Lambda$ the period lattice $\langle \int_\gamma \rho:\gamma\in H_1(X^\mathrm{an},\mathbf Z)\rangle$. The class $\gamma_\rho \in H_1(X^\mathrm{an},\mathbf C)$ which is Poincar\'e dual
to $\rho$ actually belongs to $H_1(X^\mathrm{an},\mathbf Z)\otimes \Lambda$.
Consequently $J_\xi(\gamma_\rho)$ is a well-defined element of $\mathbf C/\Lambda$.


\subsection{}
Let $X_1,X_2$ denote copies of $X$, and $X_{12}$ the diagonal copy of $X$
in $X_1\times X_2$. To a divisor
$Z\subset X\times X = X_1\times X_2$ (defined over $F$) we associate the point
\[P_Z = D_Z - \deg(D_Z)o\in \operatorname{Pic}^0(X),\]
where (recall) $o\in X(F)$ is a fixed base point and we set
$D_Z= (Z\cap X_{12})-(Z\cap X_1) - (Z\cap X_2)$.

We now state the iterated integral formula from \cite{DRS}
for the image of $P_Z$ under the Abel-Jacobi map
\[\mathrm{AJ}_X: \operatorname{Pic}^0(X) \to  \Omega^1(X^{\mathrm{an}})^\vee/H_1(X^{\mathrm{an}},\mathbf Z).\]

Let $\epsilon_o$ be the projector on $\operatorname{Pic}(X\times X)$
defined by
\[\epsilon_o(Z) = Z - i_{1*}\pi_{1*}-i_{2*}\pi_{2*}\]
where $\pi_{1},\pi_2:X\times X\rightrightarrows X$ are the projections
and $i_1,i_2:X\rightrightarrows  X\times X$ are the inclusions of ``vertical and horizontal'' copies of $X$ over the basepoint $o$.

Let $$\mathrm{cl}(\epsilon_o-):\operatorname{Pic}(X\times X)\to  H^1_{\mathrm{dR}}(X^{\mathrm{an}},\mathbf Z)\otimes H^1_{\mathrm{dR}}(X^{\mathrm{an}}, \mathbf Z)$$
denote the composition of the cycle class map and the projector $\epsilon_o$.
(The effect of $\epsilon_0$ is to annihilate $H^2\otimes H^0$ and $H^0\otimes H^2$
factors in the K\"unneth decomposition of $H^2(X\times X)$.)

\begin{theorem}[\cite{DRS}, Corollary 3.6]\label{iterformula}
Suppose $\mathrm{cl}(\epsilon_o Z)$ is represented by $\sum \omega_i\otimes \eta_i$,
where $\omega_i,\eta_i\in \Omega^1(Y)$ (one of each pair being regular at $\infty$,
since $\mathrm{cl}(\epsilon_o Z)$ is a Hodge class). 
Then the image $\mathrm{AJ}_X(P_Z)\in \Omega^1(X^{\mathrm{an}})^\vee/H_1(X^{\mathrm{an}},\mathbf Z)$
is represented by the linear functional
which maps  $\rho\in \Omega^1(X^{\mathrm{an}})$
to 
\[\sum J_{\omega_i,\eta_i}(\gamma_\rho) = \sum \int_{\gamma_\rho} (\omega_i\cdot\eta_i - \alpha_{\omega_i,\eta_i})+\deg(D_Z) \int_o^\infty \rho  \in \mathbf C,\]
where
$\gamma_\rho\in H_1(X^{\mathrm{an}},\mathbf C)$ is Poincar\'e dual
to $\rho\in H^{1,0}(X^{\mathrm{an}})\subset H^1_{\mathrm{dR}}(X^{\mathrm{an}},\mathbf C)$.\qed
\end{theorem}
% \begin{remarkwr}\label{}
%   Some remarks upon the well-definedness of the righthand side of this formula for $P_{Z,f}$ are in order. Note that while the iterated integrals $J_{\xi\otimes\xi',\alpha}$ occuring in the formula are homotopy invariant functionals on the path space of $Y(\mathbf C)$, they are being evaluated on the path $\gamma_f$. Since $\gamma_f$ is defined to be Poincar\'e dual to a certain differential $\omega_f$, it is well-defined only as an element of the quotient
% $\pi_1(X(\mathbf C),o)^{\mathrm{ab}}$ of $\pi_1(Y(\mathbf C),o)$. On the other hand, the theorem states only that the righthand side of the formula gives a well-defined value in the quotient  $\mathbf C/\Lambda_f$ of $\mathbf C$, not in $\mathbf C$ itself. 

% FIXME: Explain why the two ambiguities cancel one another to make the formula in the theorem valid on the nose.
% \end{remarkwr}

% We now describe an alternate construction of $P_Z$ in order to explain how $P_Z$
% gives rise to $F$-rational points on abelian variety quotients of $\operatorname{Jac}(X)$. 

%   one can associate a nulhomologous codimension 2 algebraic cycle $\epsilon\Delta_T \in \mathrm{CH}^2(X\times X\times X)_0$. Over the complex numbers, there is a generalized Abel-Jacobi map
%  \[\mathrm{AJ}:\mathrm{CH}^2(X^{\mathrm{an}}\times X^{\mathrm{an}}\times X^{\mathrm{an}})_0 \to J^2(X^3)\]
%  taking values in the \textit{intermediate Jacobian} $J^2(X^3)$. (FIXME: add refernece, perhaps to Voisin's book. Probably best to avoid defining $J^2(X^3)$ in this paper.) The cycle $\Pi:=\mathrm{diag}(X)\times \mathrm{diag}(X)\subset (X\times X)\times (X\times X)$, viewed as an element of $\mathrm{CH}^2(X^4)$,
%  gives rise to a map $\Phi:\mathrm{CH}^2(X^3)_0\to \mathrm{CH}^1(X)_0=\operatorname{Pic}^0(X) = (\operatorname{Jac} X)(F)$, as well as an analytic avatar of the same, $\Phi^{\mathrm{an}}:J^2(X^3)\to J^1(X)=\Omega^1(X)^\vee/ H_1(X,\mathbf Z)$, which are compatible with the generalized Abel-Jacobi map above and the classical Abel-Jacobi map $\operatorname{Pic}^0(X)\to \Omega^1(X)^\vee/H_1(X,\mathbf Z)$, in the sense that the following diagram commutes.
%  \[\xymatrix{\mathrm{CH}^2(X^3)\ar[rr]^{\mathrm{AJ}}\ar[d]_{\Phi}&&J^2(X^3)\ar[d]^{\Phi^{\mathrm{an}}}\\(\operatorname{Jac} X)(F)\ar[rr]^{\mathrm{AJ}}_{\gamma\mapsto \int_\gamma}&&\Omega^1(X)^\vee/H_1(X,\mathbf Z)}\]
%  For details (as well as generalizations) of this construction we refer to \cite{DRS1}. We mention only the result of that paper which is most relevant to our purposes. Suppose $\omega\in \Omega^1(X)$ is a holomorphic 1-form. This gives rise to an evaluation map $$\mathrm{ev}_\omega: (\operatorname{Jac}X)(\mathbf C)\stackrel{\mathrm{AJ}}{\cong}\Omega^1(X)^\vee/ H_1(X,\mathbf Z) \to \mathbf C/\Lambda_\omega$$
%  where $\Lambda_\omega=\langle \int_\gamma \omega:\gamma\in H_1(X,\mathbf Z)\rangle \subset \mathbf Z$ is the period lattice associated to $\omega$.
%  We define the point $P_{Z,\omega}^{\mathrm{an}}\in \mathbf C/\Lambda_\omega$
%  to be $\mathrm{ev}_\omega \Phi^{\mathrm{an}} \mathrm{AJ}(\epsilon\Delta_Z)$.


% Now let $\pi:\operatorname{Jac}(X)\to E$ be a homomorphism from the Jacobian of $X$
% to an elliptic curve $E_{/F}$, and assume $\ker \pi$ is connected.
% Taking $\rho = \pi^*\omega_E$ to be the pullback of the N\'eron differential on $E$



\section{Triple Chow-Heegner points on modular curves}\label{Numerics}
We now specialize the discussion of the preceding section to the case of classical modular curves $X$. We shall define certain rational points on an arbitrary elliptic curve $E_{/\mathbf Q}$ called \textit{triple Chow-Heegner points}, such that the corresponding points of $E(\mathbf C)\cong \mathbf C/\Lambda_E$ can be computed using iterated path integrals via Theorem \ref{iterformula}.
\subsection{}
Let $N\geq 1$ be an integer and let $X=X_0(N)$ denote the canonical model over $\mathbf Q$ of the classical modular curve of level $N$; write $J_0(N)$ for the Jacobian of $X_0(N)$.
With this choice of $X$ we place ourselves in the setup of \S  1,
taking the ground field $F$ to be $\mathbf Q$
and the point $\infty\in X(\mathbf Q)$ to be the usual cusp at infinity.
Thus $Y:= X_0(N)-\{\infty\}$. (Note that $Y\supsetneq Y_0(N) = X_0(N)-\{\text{cusps}\}$.) We will be deliberately vague concerning our basepoint $o\in Y^{\mathrm{an}}$ for topological constructions, but see \S 3.1 for the relevance of the choice of $o$ when performing explicit computations.

We shall make use of the \textit{Poincar\'e pairing} on $H^1(X)$, which is a symplectic form
\[\langle,\rangle:H^1_{\dR}(X/\mathbf Q)\times H^1_{\dR}(X/\mathbf Q)\to \mathbf Q.\]
If $\omega$ and $\eta$ are smooth 1-forms on $X$, then  $\langle\omega,\eta\rangle:=\frac{1}{2\pi i}\int_X \omega\wedge \eta$.
If $\omega$ and $\eta$ are differentials of the second kind on $X$, holomorphic away from the cusp $\infty$, then the induced pairing on $H^1(Y)$
can also be computed as
\[\langle \omega,\eta\rangle = \res_\infty (F_\omega \cdot \eta) =
 -\res_\infty (\omega \cdot F_\eta),\]
where as above $F_\nu$ denotes the primitive function $\tilde Y\to \mathbf C$
of the differential $\nu$. Given 1-forms $\omega,\eta$ of the second kind on $X$, regular on $Y$, the Poincar\'e pairing of their cohomology classes is thus computable from Laurent expansions of $\omega,\eta$ about $\infty$ by integrating formally.


% Let $\{ \omega_{f_1,1}, ..., \omega_{f_1,n_1},..., \omega_{f_r,1}, ..., \omega_{f_r,n_r}\}$ a basis of the space $H^0(X,\Omega^1)$ of regular differential forms on $X$ over $\Q$, which we label in such a way that, for any given normalized (primitive or not) eigenform $f\in S_2(\Gamma_0(N))$, the $f$-isotypical subspace is $H^0(X,\Omega^1)[f]=\langle \omega_{f,1}, ..., \omega_{f,n} \rangle$. With these notations, $n$ is equal to the degree of the number field $\Q_f$ generated by the fourier coefficients of $f$. And the genus of $X$ is equal to $n_1+...+n_r$. 

% For each such $f$, there exists a collection $\{ \eta_{f,1}, ..., \eta_{f,n}\}$ of differential forms of the second kind on $X$ with a single pole at the cusp at infinity such that $$H^1_{\dR}(X/\Q)[f]=\langle\omega_{f,1}, ..., \omega_{f,n}, \eta_{f,1}, ..., \eta_{f,n}\rangle$$ and 
% \begin{equation*}\label{poincare}
% \langle \omega_{f,i},\eta_{f,j}\rangle =\delta_{ij}.
% \end{equation*}

%  Here, for any differential form of the second kind $\eta$ on $X$ we write $\tilde\eta$ for the local primitive of $\eta$ about $\infty$, normalized so that $\tilde\eta(\infty)=0$. 

% Note that the choice of each $\eta_{f,i}$ is unique up to exact forms $dh$, where $h$ is a rational function on $X$ with a single pole at $\infty$. (FIXME: Is this really true? For example, if $f$ has rational fourier coefficients, and $\omega_f$,$\eta_f$ form a symplectic basis, then so do $\omega_f$, $\eta_f+c\omega_f$ for any constant $c$. I think all we are saying here is that given the cohomology class of $\eta_{f,i}$, the differential form representing it is well-defined up to exact forms.)
\subsection{}
Now let $E_{/\mathbf Q}$ be an elliptic curve of conductor $N$ whose isogeny class corresponds to a newform $f\in S_2(\Gamma_0(N))$ with rational Fourier coefficients. In particular there is a \textit{modular parametrization}
$\pi_E: J_0(N) \to E$, a homomorphism of abelian varieties defined over $\mathbf Q$.
We will assume for the remainder of this paper that $E$ is optimal, so $\ker \pi_E$ is connected. In this case the N\'eron lattice of $E$ coincides with the period lattice $\Lambda_f$ of the differential $\omega_f = 2\pi if(z)\mathrm{d}z\in \Omega^1(X^{\mathrm{an}})$ corresponding to $f$.
% We define the \textit{triple Chow-Heegner point} $P_{g,f}$
% associated to the 3-tuple of modular forms $(g,g,f)$
% to be $\pi_f\circ\Phi(\epsilon\Delta_{T_g})\in E(\mathbf Q)$.
% The main goal of this paper is to describe an algorithm for the effective computation of the points $P_{g,f}$ on an arbitrary elliptic curve $E$.
% The strategy is to compute (approximations to) the images of the points $P_{g,f}$
% in $E(\mathbf C)$ by using the iterated integral formula of Theorem \ref{iterformula}.
The map $\pi_E$ can be computed on complex points explicitly, using the Abel-Jacobi  isomorphism $\mathrm{AJ}_X:J_0(N)(\mathbf C)\cong \Omega^1(X^{\mathrm{an}})^\vee/H_1(X^{\mathrm{an}},\mathbf Z)$, the Weierstrass uniformization $W:\mathbf C/\Lambda_f\cong E(\mathbf C)$, and the analytic parametrization $$\pi_E^\mathrm{an}:\Omega^1(X^\mathrm{an})^\vee/H_1(X^\mathrm{an},\mathbf Z) \stackrel{\text{evaluate at $\omega_f$}}{\longrightarrow} \mathbf C/\Lambda_f.$$
Namely, for $P_{\mathbf C}\in J_0(N)(\mathbf C)$ we have
$\pi_{E,\mathbf C}(P_{\mathbf C}) = W(\pi_E^\mathrm{an}(\mathrm{AJ}_X(P_{\mathbf C}))) = W(\mathrm{AJ}_X(P_{\mathbf C})(\omega_f))$.

Let $\mathbf T_0 = \mathbf Z[\{T_n\}_{n\nmid N}]$ be the ``anemic'' Hecke algebra. Then $\mathbf T_0 \otimes \mathbf Q$ factors as a product $\prod_{g'} K_{g'}$ where $g'$ runs over newforms of all levels $M$ dividing $N$ and $K_{g'} = \mathbf Q(\{a_n(g')\}_{n\geq 1})$ is the number field generated by the Hecke eigenvalues of $g'$. For any divisor $M$ of $N$ and a newform $g \in S_2(\Gamma_0(M))$, denote by $T_g \in \mathbf T_0 \otimes \mathbf Q \cong \prod_{g'} K_{g'}$ the idempotent with 1 in the $K_g$ component and 0 elsewhere. This gives rise to a correspondence in $\operatorname{Pic}(X\times X) \otimes\mathbf Q$, which by abuse of notation will also be denoted by $T_g$.

\begin{definition}\label{}
  The \textit{triple Chow-Heegner point} $P_{g,f}$ corresponding to the 3-tuple $(g,g,f)$ of modular forms is the element $\pi_E(P_{T_g})\in E(\mathbf Q)\otimes \mathbf Q$.
\end{definition}
For generalizations of this definition, see for example \cite{BDP2}.
\begin{remarkwr}\label{}
  In this definition, $P_{T_g}$ is defined as in \S 1 taking $Z=T_g$. However note that $T_g$ might not literally be a divisor on $X\times X$; the correspondence $T_g$ is merely a $\mathbf Q$-linear combination of such divisors. Thus $P_{T_g}$ belongs to $\operatorname{Pic}^0(X)\otimes \mathbf Q = J_0(N)(\mathbf Q)\otimes \mathbf Q$. If we define the ``denominator'' $d_g$ of $T_g\in \mathbf T_0\otimes \mathbf Q$ to be the smallest positive integer such that $d_gT_g$ lies in the image of $\mathbf T_0$ under the inclusion $\mathbf T_0 \hookrightarrow \mathbf T_0 \otimes \m, in my opinion. athbf Q$, then $d_gP_{g,f} \in E(\mathbf Q)$ is rational. See \S 3 for more details.
\end{remarkwr}

\subsection{}
To obtain from Theorem \ref{iterformula} an explicit formula for  a triple Chow-Heegner point $P_{g,f}$ in terms of iterated integrals, we must know the components of $\mathrm{cl}(\epsilon_o T_g)\in H^1_{\dR}(X/\mathbf Q)^{\otimes2}$ when this class is decomposed as a sum of pure tensors.

The action of the Hecke algebra $\mathbf T_0$ on modular forms extends to an action on the de Rham cohomology of $X$. Under this action, we have
$$H^1_{\mathrm{dR}}(X/\mathbf Q) \cong H^1_{\mathrm{dR}}(X/\mathbf Q)[g_1] \oplus \cdots \oplus H^1_{\mathrm{dR}}(X/\mathbf Q)[g_n],$$
indexed by Galois conjugacy classes of newforms of all levels $M$ dividing $N$.

\begin{lemma}\label{}
Let $M|N$ and let $g\in S_2(\Gamma_0(M))$ be a newform.
  Let $\{\omega_{g,1},\ldots, \omega_{g, k},\eta_{g,1},\ldots, \eta_{g,k}\}$
be a collection of differentials of the second kind on $X$ representing a basis for $H^1_{\dR}(X/\mathbf Q)[g]$ that is symplectic with respect to the Poincar\'e pairing; i.e. assume
$\langle \omega_{g,i},\eta_{g,j}\rangle = \delta_{i,j}$ and $\langle \omega_{g,i},\omega_{g,j}\rangle = \langle \eta_{g,i},\eta_{g,j}\rangle = 0.$ % (Kronecker $\delta$).\footnote{Here $\langle,\rangle$ denotes the symplectic pairing $H^1(X)\times H^1(X)\stackrel{\wedge}{\to} H^2(X) \stackrel{\int_X}{\to} \mathbb C$.
% The basis $\{\omega_{g,i},\eta_{g,j}\}_{1\leq i,j\leq k}$ could thus in particular be taken to be a standard ``symplectic'' or ``Darboux'' basis for the symplectic vector space $(H^1(X),\langle,\rangle)$.}
Then  $\mathrm{cl}(\epsilon T_g) = \sum_{i=1}^k \omega_{g,i}\otimes \eta_{g,i} - \eta_{g,i}\otimes \omega_{g,i}.$
\end{lemma}
\begin{proof}
% FIXME: Fill in the details of the proof sketched in Victor's email of May 7.

% Let $\xi_g = \mathrm{cl}(T_g)$. By definition this is an element of $H^2(X\times X)$ and it is characterized by the property that for all $\xi'\in H^2(X\times X)$ we have
% \[\int_{X\times X} \xi_g\wedge \xi' = \int_{T_g} \xi.\]
% This is because the cycle class of $T_g$ is Poincar\'e dual to the homology class of $T_g$, viewed as an analytic subset of $X\times X$.
% (Here, according to Voisin, vol. 1, cor. 11.15, we do \textit{not} pick up a $2\pi i$ factor, so the above should be true ``on the nose''.)

% (Remark: If $T_g$ were smooth, $\xi_g$ can be taken to be the class of a smooth $2$-form with proper support contained in a tubular neighborhood $V$ of $T_g\subset X\times X$ (i.e. a neighborhood isomorphic to a neighborhood of the zero section in the normal bundle of $T_g$), satisfying the condition that
% $\int_{V_z}\xi_g = 1$ for all $z\in T_g$, where $V_z$ is the preimage in $V$
% of the point $z$ under the structural map $N_{T_g/X\times X}\to T_g$.)

% By Voisin, vol 1, prop. 11.20, $\xi_g$ belongs to $H^{1,1}(X)$.

% By the K\"unneth formula, $H^2(X\times X) = H^2(X)\otimes H^0(X)\oplus H^1(X)\otimes H^1(X) \oplus H^0(X)\otimes H^2(X)$. We may assume wlog 
% that $\xi_g = \mathrm{cl}_{1,1}(T_g)$ is concentrated in the (1,1) component.
% (FIXME: why? because this should amount to removing vertical and horizontal components from the correspondence $T_g$, and that should be harmless either because there are none, or because those components have no effect on cohomology or, ...?)

% We can decompose $H^1(X) = \bigoplus_{h} H^1(X)[h]$
% where $h$ runs through Galois conjugacy classes of eigenforms in $S_2(X)$.
% Since $\mathrm{cl}:\mathrm{CH}^1(X\times X)\to H^2(X\times X)$
% is Hecke-equivariant, $\xi_g$ belongs to $H^1(X)[g]\otimes H^1(X)[g]$.
% (JUSTIFY)

% Suppose $g$ has rational Fourier coefficients,
% so $H^1(X)[g]$ is 2-dimensional, spanned by $\omega_g\in \Omega^1(X)[g]$
% (which has rational Fourier coefficients)
% and its complex conjugate $\eta_g = \bar\omega_g$. 
% FIXME: $\eta_g$ should be cohomologous to a differential of the second kind with rational fourier coefficients, satisfying $\langle\omega_g,\eta_g\rangle = 1$, since $\int_X \omega_g \wedge \bar\omega_g = \int_X |\omega_g|^2 = 1$
% if $\omega_g$ is suitably normalized. ???

% In general, suppose $H^1(X,\mathbf Q)[g]$ has $\mathbf Q$-basis 
% given by the classes of holomorphic differentialss $\omega_{g,1},\ldots,\omega_{g,k}$ and differentials of the second kind $\eta_{g,1},\ldots, \eta_{g,k}$,
% and assume this basis is symplectic. 
% As noted above, by Hodge theory, the image of $\xi_g\in H^2(X\times X,\mathbf Z)$ inside $H^2(X\times X, \mathbf C)=\bigoplus H^{p,q}(X\times X)$
% lands in $H^{1,1}(X\times X) \cap H^1(X)\otimes H^1(X)$.
% Thus we can write
% \[\xi_g = \sum a_{i,j} \omega_{g,i}\otimes \eta_{g,j} + b_{i,j}\eta_{g,i}\otimes \omega_{g,j}.\]
% Let $\xi' = \omega_{g,i}\otimes \eta_{g,j}$.
% Then
% \[\int_{T_g} \xi' = \int_{X\times X}\xi_g\wedge \xi'
% = \sum_{i',j'}\int_{X\times X}(a_{i',j'}\omega_{g,i'}\otimes\eta_{g,j'}
% + b_{i',j'}\eta_{g,i'}\otimes \omega_{g,j'})\wedge(\omega_{g,i}\otimes \eta_{g,j})\]
% \[=\sum_{i',j'}a_{i',j'} \left(\int_X \omega_{g,i'}\wedge \omega_{g,i}\right)\left(\int_X \eta_{g,j'}\wedge \eta_{g,j}\right) + b_{i',j'}\left(\int_X \eta_{g,i'}\wedge \omega_{g,i}\right)\left(\int_X \omega_{g,j'}\wedge \eta_{g,j}\right)\]
% by Fubini.
% This last expression is precisely
% \[\sum_{i',j'} a_{i',j'}\langle \omega_{g,i'},\omega_{g,i}\rangle\langle \eta_{g,j'},\eta_{g,j}\rangle + b_{i',j'}\langle \eta_{g,i'},\omega_{g,i}\rangle\langle \omega_{g,j'},\eta_{g,j}\rangle
% = - b_{i,j}.\]
% (Because we chose a symplectic basis.)
% By analogous reasoning one finds that for $\xi'' = \eta_{g,i}\otimes \eta_{g,i}$, one has
% \[\int_{T_g} \xi'' = a_{i,j}.\]
% But because of the definition of $T_n = \{(E,E')\in X\times X: E, E'\text{ are related by a cyclic isogeny of degree }n\}$,
% $T_g$ is a linear combination of symmetric correspondences, so it is itself symmetric. Thus we find $a_{i,j}=-b_{i,j}$.
By definition $T_g$ is a correspondence which acts on $H^1(X)$
as the idempotent projector onto $H^1(X)[g]$. 
The $H^2\otimes H^0$ and $H^0\otimes H^2$ K\"unneth components of $\mathrm{cl}(T_g)$ act trivially on $H^1(X)$. (See \cite[11.5.1]{BL}, for example.)
Thus $\mathrm{cl}(\epsilon_o T_g)\in H^1(X)\otimes H^1(X)$ also acts by projecting onto $H^1(X)[g]$.

The action on $H^1(X)$ of a correspondence $Z\subset X\times X$ whose cycle class is in $H^1(X)\otimes H^1(X)$ can be written in terms of the Poincar\'e pairing.
Using that $\mathrm{cl}(\epsilon_o T_g)$ is a projector on $H^1(X)[g]$, one finds
that $\mathrm{cl}(\epsilon_o T_g) = \sum_{i,j} \langle b_i, b_j\rangle b_i\otimes b_j$ for \textit{any} basis $\{b_1,\ldots, b_{2k}\}$ of $H^1(X)[g]$. From the the claim follows immediately.
\end{proof}
Combining the previous results, we obtain the following formula for $P_{g,f}$. Let $\gamma_f$ be the Poincar\'e dual of $\omega_f$
and let $\omega_{g,1},\ldots, \omega_{g,k},\eta_{g,1},\ldots, \eta_{g,k}$
be (differentials of the second kind which give rise to) a symplectic basis for the $g$-isotypic $\mathbf Q$-subspace 
$H^1(X/\mathbf Q)[g]\subset H^1(X/\mathbf Q)$. Then, recalling that $W$ denotes the Weierstrass uniformization of $E(\mathbf C)$, we find 
that $P_{g,f}\in E(\mathbf Q)\otimes \mathbf Q\subset E(\mathbf C)\otimes \mathbf Q$ can be computed as
\begin{equation}\label{pgf}
  P_{g,f} = W\left(\sum_{i=1}^k \left(\int_{\gamma_f} \omega_{g,i}\cdot \eta_{g,i}-\eta_{g,i}\cdot\omega_{g,i} - 2\alpha_{\omega_{g,i},\eta_{g,i}}\right)\right).
\end{equation}
We emphasize that by Lemma \ref{homotopyinv} (and the discussion immediately following it) the righthand side of \eqref{pgf} depends only on the \textit{homology} class
$\gamma_f\in  H_1(Y^\mathrm{an}, \mathbf Z) = H_1(X^\mathrm{an},\mathbf Z)$ Poincar\'e dual to $\omega_f$.
It can therefore be evaluated by lifting $\gamma_f$ arbitrarily to an element $\tilde \gamma_f\in \pi_1(Y^\mathrm{an};o)$
and evaluating $\sum \omega_{g,i}\cdot\eta_{g,i}-\eta_{g,i}\cdot\omega_{g,i}-2\alpha_{\omega_{g,i},\eta_{g,i}}\in \II^{\leq 2}(Y)$ on any loop in the homotopy class $\tilde\gamma_f$.
%  \begin{example} 
% When $g$ has \textit{rational} Fourier coefficients the $g$-isotypic component is 2-dimensional, and the above formula reads
% $$
% P_{g,f} := \int_{\gamma_f} (\omega_g\cdot \eta_g - \eta_g \cdot \omega_g) -2\int_{\gamma_f} \alpha_{\omega_g,\eta_g} =$$ $$=2(\int_{\gamma_f} \omega_g\cdot \eta_g -\int_{\gamma_f}\alpha_{\omega_g,\eta_g})-\int_{\gamma_f} (\omega_g\cdot \eta_g + \eta_g \cdot \omega_g) \in \C/\Lambda_f \sim E(\C).
% $$
% \end{example}
\subsection{}
Recall that in the above definition of iterated integral, everything depends on the choice of a base point $o$. Likewise, the projector $\epsilon_o$ depends on $o$, and hence \textit{a priori} so does the point $P_{g,f}$. However we have the following.
\begin{lemma}\label{basepointindependence}
The point $P_{g,f}$ is independent of $o$.
\end{lemma}
% FIXME: Should this lemma be in the preceeding section? Is the theorem from DRS valid (i.e. is the RHS of the formula for $P^{\mathrm{an}}_{Z,\omega}$ even well-defined) when $\mathrm{cl}(Z)$ is not orthogonal to $\omega$?
\begin{proof}
Changing the basepoint from $o$ to $o'$ amounts to conjugating the representative path $\gamma_f$ for the homology class Poincar\'e dual to $\omega_f$
by a path $\beta$ from $o$ to $o'$. This manifestly does not affect
the value of the  integral of the meromorphic 1-form $2\alpha_{\omega_{g,i}, \eta_{g,i}}$. Thus the issue is whether we have an identity
\begin{equation}\label{basepoint}
  \int_{\gamma_f} \omega_{g,i}\cdot\eta_{g,i}-\eta_{g,i}\cdot\omega_{g,i} \stackrel{?}{=}\int_{\beta\gamma_f\beta^{-1}}\omega_{g,i}\cdot \eta_{g,i}-
\eta_{g,i}\cdot\omega_{g,i}.
\end{equation}

But by \cite[Exer. 8]{Hain}, for any 1-forms $\omega,\eta$,
loop $\gamma$, and path $\beta$,  we have
\begin{equation}\label{cobp}
\int_{\beta\gamma\beta^{-1}} \omega\cdot\eta = \int_\gamma \omega\cdot \eta + \left|\begin{smallmatrix}\int_\gamma \omega&\int_\gamma \eta\\
\int_\beta \omega&\int_\beta \eta\end{smallmatrix}\right|.  
\end{equation}
In our situation,
 the determinants expressing the difference between the two sides
of \eqref{basepoint} vanish.
Indeed, $\int_{\gamma_f} \omega_{g, i} = \langle \omega_{g, i},\omega_f\rangle = 0 = \langle \eta_{g, i},\omega_f\rangle = \int_{\gamma_f} \eta_{g, i}$,
since the decomposition into isotypic components for the action of the Hecke algebra is orthogonal with respect to the Poincar\'e pairing.
\end{proof}

% However, the point $P_{g,f}$ {\em does depend} on the choice of the cusp $\infty$ we chose in order to define the open curve $Y\subset X= X_0(N)$. (Note that in general $X_0(N)$ has several cusps.)
% We should keep this in mind, and if necessary we may actually write it as $P_{\infty, g;f}$. 
\subsection{}
We record a fundamental property of the points $P_{g,f}$.
\begin{theorem}[\cite{DRS}, Theorem 1]\label{}
Assume that the local signs of Garrett's triple product L-function $L(g,g,f,s)$ at the primes $p\mid \mathrm{gcd}(M,N)$ are all $\varepsilon_p(g,g,f)=+1$.
Then the point $P_{g,f}\in E(\mathbf Q)\otimes \mathbf Q$ is nonzero
(or equivalently, the point $\pi_E(d_g P_{T_g})\in E(\mathbf Q)$ is non-torsion) if and only if the following three conditions hold:
\begin{itemize}
\item[i.] $L(f,1)=0$,
\item[ii.] $L'(f,1)\neq 0$, and
\item[iii.] $L(f\otimes \operatorname{Sym}^2(g),2)\neq 0$.
\end{itemize}
\end{theorem}


\section{Algorithm for effective computation of triple Chow-Heegner points}
We now turn to the question of numerically evaluating formula \eqref{pgf} for a triple Chow-Heegner point $P_{g,f}\in E(\mathbf Q)\otimes\mathbf Q$ for an optimal elliptic curve $E = E_f$. We retain all the notation from \S \S 1-2.

% \begin{remarkwr}\label{}
% Some of the previous assumptions could be relaxed: there are rational points $P_{g_1;g_2;f}\in J_0(N)(\mathbf Q)$ associated to any 3-tuple of modular forms
% (whose respective weights satisfy a certain combinatorial condition).
% Relaxing the hypothesis that $f$ has {\em rational} fourier coefficients 
% amounts to investigating the images of these rational points in an abelian variety quotient $A_f$ of $J_0(N)$ which has dimension $>1$.
%  Relaxing the  assumption that $E$ has rank 1 is only of moderate interest: 
% by \cite{DRS1}, when $E$ has (analytic?) rank $\neq 1$, the resulting points $P_{g,f}$ are all torsion.  
%  The hypothesis on the sign of the functional equation of $g$ is made so that the sign of the functional equation of $L(\mathrm{Sym}^2(g)\otimes f,s)$ is $-1$ and thus $L'(\mathrm{Sym}^2(g)\otimes f,2)$ can be (and is expected often to be) non-zero.
% \end{remarkwr}

The following ingredients occur in the formula \eqref{pgf} for $P_{g,f}$.
\begin{itemize}
\item[1.] The Poincar\'e dual $\gamma_f\in H_1(X,\mathbf C)$
of $\omega_f\in H^1_{\dR}(X^{\mathrm{an}},\mathbf  C)$. 
\item[2.] A collection of rational differentials of the second kind $\omega_{g,1},\ldots, \omega_{g,k},\eta_{g,1},\ldots, \eta_{g,k}$ on $X$, regular away from $\infty$, whose images in $H^1_{\dR}(X/\mathbf Q)$ are a symplectic basis for the $g$-isotypic component $H^1_{\dR}(X/\mathbf Q)[g]$.
\item[3.] Meromorphic differentials $\alpha_{\omega_{g,i},\eta_{g,i}}$ on $X$, regular on $Y$, such that $\omega_{g,i} F_{\eta_{g,i}}-\alpha_{\omega_{g,i},\eta_{g,i}}$ has at worst a logarithmic pole at (any point lying over) $\infty$.
\end{itemize}
These data must be ``known'' in a sufficiently concrete form
to evaluate the iterated integrals occuring in \eqref{pgf}.
It is also desirable to know 
\begin{itemize}
\item[4.] the denominator $d_g$ of the projector
onto the $g$-isotypic component of the chomology of $X$.
\end{itemize}
This last item will allow for the computation of a point in $E(\mathbf Q)$, as opposed to one in $E(\mathbf Q) \otimes \mathbf Q$. This section is devoted to methods of computing these four ingredients.

\subsection{Evaluating iterated integrals}\label{integration}
Let $J = \sum \omega_i\cdot \eta_i -\alpha_i\in \II^{\leq 2}(Y)$ be a homotopy-invariant iterated integral of length $\leq 2$ on $Y$, expressed in terms of differentials of the second kind on $X$, regular on $Y$. We seek to compute the righthand side of formula \eqref{pgf}, which is $J(\gamma)$ for a particular choice of $J$ and (homotopy class of) path $\gamma\in \pi_1(Y;o)$. As remarked earlier, said formula actually depends only on the homology class $\gamma_0$ of of $\gamma$.
This homology class belongs to  $H_1(Y^\mathrm{an},\mathbf Z) = H_1(X^{\mathrm{an}},\mathbf Z)$, which
is the abelianization of the quotient $\pi_1(X^{\mathrm{an}},o)=\bar \Gamma_0(N)$
of $\Gamma_0(N)$ by the smallest normal subgroup containing the elliptic and parabolic elements. 
To evaluate $J(\gamma_0)$ for $\gamma_0\in H_1(Y^{\mathrm{an}},\mathbf Z)$,
we lift $\gamma_0$ arbitrarily to a path $\tilde \gamma$ in $\tilde Y$ based at $o$.
If we choose the basepoint $o$ away from the divisor of cusps on $X$,
then $o$ can be lifted to an element $\tau_0$ in the upper half-plane
$\mathfrak H$, regarded as a cover\footnote{FIXME: this cover is ramified at the elliptic points, so some additional care in required to make this step rigorous.} of $Y_0(N)$. 


The path $\tilde \gamma$ can then be viewed as a path in $\mathfrak H$ from $\tau_0$ to $\gamma\tau_0$, where $\gamma\in \Gamma_0(N)$ is a lift
of $\gamma_0$. 

\begin{lemma}\label{}
Suppose $\gamma_0$ is Poincar\'e-dual to $\rho$.
As an element of $\mathbf C/\Lambda_\rho$, we have
  $$J(\gamma_0) = \sum \int_{\tau_0}^{\gamma \tau_0} \omega_i F_{\eta_i} - \alpha_{\omega_i,\eta_i}$$
where we conflate 1-forms on $X$ with their pullbacks to $\mathfrak H^* = \mathfrak H\cup\{\infty\}$.
Moreover, $F_{\eta_i}$ has Laurent expansion about $\infty\in \mathfrak h^*$
given by formally integrating the Laurent expansion of $\eta_i$
about the cusp $\infty\in X$.
\end{lemma}
\begin{proof}
Clear.
\end{proof}

Given any differential form $\lambda$ of the second kind on $X$, and any $\gamma\in\Gamma_0(N)$, let
$$ I(\lambda;\gamma) := \int_{\tau_0}^{\gamma\tau_0} \lambda.$$
(As above, in the righthand side of this expression $\lambda$ is conflated with its pullback to $\mathfrak H^*$.)
By the residue formula,
this expression is independent of the choice of path on the upper half-plane $\mathfrak H$  from $\tau_0$ to $\gamma\tau_0$. The $\Gamma_0(N)$-invariance of
$\lambda$ also shows that this expression is independent of the choice of  base point 
$\tau_0\in \mathfrak H$, which justifies suppressing $\tau_0$ from the notation. 

If $\lambda$ instead denotes a differential of the second kind on  $\tilde X$ then the integral above still makes sense
but depends on both the basepoint $o$ and the chosen lift of $o$ to $\tau_0\in \mathfrak H$. We will primarily be interested in evaluating such integrals in the context of \eqref{pgf}, for which the choice of basepoint is ultimately irrelevant. (This is because the Poincar\'e dual of the homology class of $\gamma$, is orthogonal
to the 1-forms in the iterated integral giving rise to the path integral we seek to evaluate; cf. Lemma \ref{basepointindependence}.)
However, as we are about to see, for the purposes of algorithmic efficiency it is necessary to break up the path of integration into pieces which can be computed relatively quickly. The integrals over these pieces may no longer be basepoint-independent: when we express $\gamma$ as a product of computationally-amenable elements $\gamma^{(j)}\in \Gamma_0(N)$, the corresponding 
homology classes may no longer lie in (the Poincar\'e dual of) the orthogonal complement of $H^1_{\mathrm{dR}}(X/\mathbf Q)[g]$.
Thus for a general meromorphic 1-form $\lambda$ on $\tilde X$ and a general $\gamma\in\Gamma_0(N)$, we adopt the notation
\[I_{\tau_0}(\eta;\gamma) = \int_{\tau_0}^{\gamma\tau_0}\eta\]
to emphase the dependence on the choice of basepoint.

By meromorphicity, for $\lambda$ as above (defined on either $X$ or $\tilde X$) the integral $I(\lambda;\gamma)$ or $I_{\tau_0}(\lambda;\gamma)$ can be computed by integrating term-wise a Laurent expansion for $\lambda$
using the fundamental theorem of calculus.
Thus, in practice, one computes the Laurent expansion for the primitive 
$F_\lambda$ about $\infty \in X$ (or a choice of $\tilde \infty\in \tilde X$ lying over $\infty$), regarded as function given by a convergent power series
in $q = e^{2\pi i \tau}$ on $\mathfrak h$, and evaluates it at $\tau_0$ and $\tau_0'=\gamma \tau_0$. The larger the imaginary parts of $\tau_0$ and $\tau'_0$ are, the faster this series converges and the fewer coefficients of the Laurent series of $\lambda$ are necessary to approximate $I(\lambda;\gamma)$ or $I_{\tau_0}(\lambda;\gamma)$ to a give degree of accuracy.
Writing $\gamma = \begin{smallmat} abcd\end{smallmat}$, it is well-known that the best compromise between $\mathrm{Im}(\tau_0)$ and $\mathrm{Im}(\tau'_0)$ is achieved when we choose $\tau_0=-\frac{d}{c} + \frac{1}{|c|} i$ (cf., for example, \cite[p. 35]{cremona}).
This {\em optimal basepoint for $\gamma$} will be denoted $\tau^*_\gamma$.

With this remark in mind, we take the following approach to computing $J(\gamma_0)$ as in the lemma above. First compute Laurent expansions for the differentials $\omega_i,\eta_i, \alpha_{\omega_i,\eta_i}$. Then find a ``good'' expression for the homology class $\gamma_0\in H_1(X^{\mathrm{an}},\mathbf C)$, writing it as a $\mathbf C$-linear combination of classes $\gamma_0^{(j)}\in H_1(X^{\mathrm{an}},\mathbf Z)$ which lift to elements $\gamma^{(j)}\in \Gamma_0(N)$ with small lower-left entries $cN$. Finally, calculate approximations to the integrals $I_{\tau_0}(\omega_i F_{\eta_i}; \gamma^{(j)})$ and $I(\alpha_{\omega_i,\eta_i};\gamma^{(j)})$. The appropriate linear combination of these integrals is an (approximate) representative for the coset $J(\gamma_0)\in \mathbf C/\Lambda_\rho$.)

To calculate $I(\alpha_{\omega_i,\eta_i};\gamma^{(j)})$, one is free to change the basepoint from $\tau_0$
to the optimal basepoint $\tau_j^*:=\tau_{\gamma^{(j)}}^*$ for $\gamma^{(j)}$, since $\alpha_{\omega_i,\eta_i}$ is defined on $X$
and not only on $\tilde X$. The same is \textit{not} true for $\omega_i F_{\eta_i}$.
To evaluate $I_{\tau_0}(\omega_i F_{\eta_i};\gamma^{(j)})$ we appeal to the following lemma.
\begin{lemma}\label{}
$I_{\tau_0}(\omega_i F_{\eta_i};\gamma^{(j)}) = I_{\tau_j^*}(\omega_i F_{\eta_i}; \gamma^{(j)}) -
I(\eta_i;\gamma^{(j)}) \int_{\tau_0}^{\tau_j^*}\omega_i$.
\end{lemma}
Observe that every term on the righthand side can be computed using the fundamental theorem of calculus, evaluating powerseries only at the points $\tau_0$ and $\tau_j^*$. In particular, taking $\tau_0 = i/N$,
each such evaluation converges at least as fast as an evaluation at $\tau_j^*$, so this formula for the integral is ``optimally efficient''. 
\begin{proof}[Proof of the lemma]
Since $\lambda = \omega_i F_{\eta_i}$ is a holomorphic 1-form on $\mathfrak H$, its integral along a closed contour vanishes. Thus 
\[I_{\tau_0}(\lambda;\gamma^{(j)}) = I_{\tau_j^*}(\lambda;\gamma^{(j)}) + \int_{\tau_0}^{\tau_j^*}\lambda 
- \int_{\gamma^{(j)} \tau_0}^{\gamma^{(j)} \tau_j^*}\lambda.\]
To evaluate the second term on the righthand side, we observe
that $\omega_i$ comes from a 1-form on $X$, so it is $\Gamma_0(N)$-invariant; it thus pulls back to itself along the fractional linear transformation defined by
$\gamma^{(j)}$.
On the other hand,
\[I(\eta_i;\gamma^{(j)})=\int_{\gamma^{(j)}}\eta_i = F_{\eta_i}(\gamma^{(j)}\tau) - F_{\eta_i}(\tau),\qquad \text{for all }\tau\in \mathfrak H.\]
Hence $(\gamma^{(j)})^* F_{\eta_i} = F_{\eta_i} + I(\eta_i;\gamma^{(j)})$.
So $(\gamma^{(j)})^* \lambda = \lambda + I(\eta_i;\gamma^{(j)})\omega_i$, and we find
\[\int_{\gamma^{(j)} \tau_0}^{\gamma^{(j)} \tau_j^*}\lambda = \int_{\tau_0}^{\tau_j^*} (\gamma^{(j)})^*\lambda
= \int_{\tau_0}^{\tau_j^*}\lambda + I(\eta_i;\gamma^{(j)}) \int_{\tau_0}^{\tau_j^*}\omega_i,\]
which yields the lemma.
\end{proof}

\begin{remarkwr}\label{}
We warn the reader  that possibly $I_{\tau_0}(\omega_i F_{\eta_i};\gamma^{(j)}) \neq \int_{\gamma^{(j)}}\omega_i\cdot \eta_i$ (regarding $\gamma^{(j)}$ as an element of $\pi_1(Y^\mathrm{an};o)$.
Indeed, the iterated integral $\omega_i\cdot \eta_i$ need not even be homotopy invariant (!) so $\int_{\gamma^{(j)}}\omega_i\cdot \eta_i$ is ill-defined. In particular, one \textit{cannot}
relate $I_{\tau_0}(\omega_i F_{\eta_i};\gamma^{(j)})$ to $I_{\tau_j^*}(\omega_i F_{\eta_i};\gamma^{(j)})$
using the change-of-basepoint formula \eqref{cobp} for iterated integrals.
\end{remarkwr}


% This is indeed possible thanks to the following

% \begin{lemma}
% Let $\bar{\Gamma}_0(N)$ denote the quotient of $\Gamma_0(N)$ by the normal closure of the subgroup generated by elliptic and parabolic elements. Then $\bar{\Gamma}_0(N)/\bar{\Gamma}_0(N)^{ab}$ is generated by elements $\gamma = \begin{smallmat}abcd\end{smallmat}\in \Gamma_0(N)$ with $c=N$.
% \end{lemma}

% \begin{proof}
% To be given. Perhaps only works under some hypothesis on $N$ (e.g. $N$ is prime). Modify the statement accordingly if necessary.
% \end{proof}
% FIXME: I think this lemma is just false. It fails for $N=40,42,88,105$ for example. I guess it is possible it works when $N$ is prime, but that's really too restrictive for our purposes anyway. (Also I have no idea how to prove it even for prime $N$.) Regardless, we should be able to work around it, although currently the code seems to give nonsense when we modify hypgens(N) in such a way that it first tries to produce generators as in the lemma, then allows matrices with $c=2N$, then $c=3N$, etc., until a generating set is found.
% I think the issue is basepoints, but I am currently quite confused. -- Sam
% \begin{remarkwr}\label{}
%   Even when the previous lemma\footnote{Assuming we can prove it!} does not apply, there is a simple algorithm to produce a set of generators $\gamma$ as in the statement of the lemma (when such exists). Indeed, one simply iterates through hyperbolic elements of the desired form until the corresponding homology classes span the integral homology of $X_0(N)$, as can be checked via standard modular symbols algorithms. 
% \end{remarkwr}

To efficiently evaluate the integrals in \eqref{pgf} using the method just explained, it is therefore necessary to know:
\begin{itemize}{}
\item[a.]  the homology class $\gamma_f$ as a $\mathbf C$-linear combinationation
of class $\gamma_0^{(j)}$ whose lifts to $\Gamma_0(N)$ have small lower-left entries $cN$; and,
\item[b.] Laurent expansions about $\infty$ for a symplectic basis $\omega_{g,i},\eta_{g,j}$
of $H^1_{\mathrm{dR}}(X/\mathbf Q)[g]$ and the forms $\alpha_{\omega_{g,i},\eta_{g,i}}$.
\end{itemize}
In the rest of this section we turn to the task of computing these data.




\subsection{Calculating a symplectic basis for $H^1_{\mathrm{dR}}(X/\mathbf Q)[g]$}\label{symplecticbasis}
The calculation of a basis for the deRham cohomology can be carried out
by first writing down a modular function $u$ -- that is, a rational function on $X=X_0(N)$ -- which is regular away from $\infty$. Such a function exists by Riemann-Roch and a $q$-expansion for one such function can be computed explicitly using the Dedekind eta-function, as explained in the next subsection.

Using a modular symbol algorithm, one can compute $q$-expansions for a basis
of $S_2(\Gamma_0(N), \mathbf Q)$; cf. \cite{S2}, for example. Write $\omega_1,\ldots, \omega_{t}$ for the corresponding holomorphic 1-forms on $X$,
where for convenience we set $t = p_a(X) = \dim S_2(\Gamma_0(N))$.


Define $\eta_1=u\omega_i$, which is a differential of the second kind by the residue theorem, and let $\mathcal B = \{ \omega_1,...,\omega_{t}, \eta_1,...,\eta_{t}\}\subset H^1_{\mathrm{dR}}(X/\mathbf Q)$ be the corresponding set of cohomology classes.  A simple application of Riemann-Roch shows the following.
\begin{lemma}\label{h1basis}
  The set $\mathcal B$ is basis for $H^1_{\mathrm{dR}}(X/\mathbf Q)$ whenever $\infty$ is not a Weierstrass point on $X$ and $u$ has a pole of order $t+1$ (i.e., the smallest possible) at $\infty$. \qed
\end{lemma}
\begin{proof}
Since $\infty$ is not a Weierstrass point on $X$, we may assume that $\ord_\infty (\omega_i) = i-1$, and thus $\ord_\infty (\eta_i) = i-t-2$. For any differential of the second kind $\omega'$, we can find a linear combination of $\eta_1,\ldots,\eta_t$ and $dh$ for an appropriate rational function $h$ having the same principal part as $\omega'$. Thus the difference is holomorphic, and lies in the span of $\{\omega_1,\ldots,\omega_t\}$.
\end{proof}
\begin{remarkwr}\label{oggremark}
By a result of Ogg \cite{Ogg}, the cusp $\infty$ is not a Weierstrass point when the level $N$ is prime, or more generally when $N = pM$ for prime $p$ and
an integer $M\geq 1$ such that $X_0(M)$ has genus zero and $p \nmid M$. Even if $u$ has a pole of order $> g(X)+1$, the set $\mathcal B$ may still be a basis of $\mathrm H^1_{\mathrm{dR}}(X/\mathbf Q)$. This can be checked by computing the matrix for the Poincare pairing, and in every example we have computed this is the case.

When $\infty$ is a Weierstrass point, there is a rational function with a single pole at $\infty$ of order $\leq g(X)$. When $u$ is taken to be such a function, then the set $\mathcal B$ will never be a basis. Indeed, since $\infty$ is a Weierstrass point, there exists a holomorphic differential form $\omega$ with order of vanishing $\geq g(X)$ at $\infty$. Then $u\omega$ is still holomorphic, and thus lies in the span of $\{\omega_1,\ldots,\omega_t\}$. But $u\omega$ also is in the span of $\{\eta_1,\ldots,\eta_t\}$ by definition of the $\eta_i$, giving rise to a linear dependence relation. Hence, in order for $\mathcal B$ to be a basis, it is necessary for $u$ to have a pole at $\infty$ of order greater than the order of vanishing at $\infty$ of any holomorphic differential.
% For example, of the 60 squarefree values of $N$ in $[2,100]$, 26 have the property
% that the optimal $u$ has a pole of order $g(X_0(N))+1$ at infinity.
\end{remarkwr}

Given one basis $\mathcal B$ for $ H^1_{\mathrm{dR}}(X/\mathbf Q)$ -- for example, one computed as above -- it is then a matter of linear algebra to produce a better basis which is adapted to the action of the Hecke algebra. Note that the usual formula for the action of $\mathbf T$ on holomorphic modular forms in terms of $q$-expansions extends to weakly holomorphic modular forms, such as 1-forms of the second kind on $X$, as well.\footnote{FIXME: add reference?} Thus, using $q$-series
for the elements of the basis $\mathcal B$, we can write down the
matrix $[T_p]\in\mathrm{Mat}_{2t\times 2t}(\mathbf Q)$ which describes the action of any $T_p\in \mathbf T$
on de Rham cohomology. By finding the eigenspaces of finitely many such matrices\footnote{FIXME: Add a reference to the bound on the number of generators of $\mathbf T$ acting on modular forms of level $N$; must explain why a similar \textit{effective} bound holds also for all of $H^1$. Remark that in practice only very small Hecke operators are required when $N$ is small, but that using our strategy to write down an essentially ``random'' basis in terms of $\omega$s and $u\omega$s, the rational numbers showing up even in the matrix of $T_2$ seem to grow complicated exponentially fast (as a function of $N$).} we can write down $\mathbf Q$-bases for each isotypic component of $H^1$. 
Using these it is elementary to produce the desired symplectic bases $\{\omega_{g,i},\eta_{g,j}\}$ for each isotypic component $H^1_{\mathrm{dR}}(X/\mathbf Q)[g]$.

\subsection{Modular units and $\eta$-products}\label{eta} 
The preceding discussion raises the question of how to compute the rational function $u$ used to write down an initial choice of basis $\mathcal B$ for $H^1_{\mathrm{dR}}(X/\mathbf Q)$. To ``compute $u$'' means to compute its Laurent expansion about $\infty$.

%  When $\infty$ is not a Weierstrass point

% Having a $q$-expansion for this function $u$ -- and also knowing that $u$ has the smallest possible pole order at $\infty$ ($t+1$ in the case that $\infty$ is non-Weierstrass) -- is also relevant for the computation of the differentials $\alpha_{\omega_{g,i},\eta_{g,i}}$ explained below.

Recall that the \textit{modular units} $U$ (for $\Gamma_0(N)$) are the multiplicative group 
of modular functions $u\in \mathbf C(X)^\times$ with divisor supported on the cusps of $X=X_0(N)$. % A certain subgroup of $U$ is amenable to explicit computations.
\begin{definition}\label{etadef}
The {\em eta group} $U_\eta$ is the group of rational functions $u\in \mathbf Q(X)$ of the form $$u(q) =
\prod_{0<d\mid N} \eta(q^d)^{r_d},$$ where $\eta(q)=q^{1/24} \prod_{n>0} (1-q^n)$ is the classical eta function,
 and $\{ r_d\}_{d\mid N}$ is a collection of integers satisfying the following conditions.
\begin{itemize}
\item[i.]  $\sum_{d\mid N} r_d =0$,
\item[ii.] $\prod_{d\mid N} d^{r_d}\in \mathbf Q^\times$ is a square,
\item[iii.] $(n_d):=A_N\cdot (r_d)$ is a vector (indexed by divisors $d$ of $N$) of integers divisible by 24, where $A_N$ is the $\sigma(N)\times \sigma(N)$-matrix whose entry indexed by $(d,d')$ is $\frac{N\cdot (d,d')^2}{d d' (d',N/d')}$.
\end{itemize}
\end{definition}
Work of Newman and Ligozat shows that such functions are indeed modular units on $X$; that is, $U_\eta\subset U$. In fact more is true:

\begin{proposition}\label{QUeta}
$\mathbf{Q}\otimes U_\eta = \mathbf{Q}\otimes U$.
\end{proposition}
\begin{proof}
It easy to see that the set $\{ \frac{a}{d}:\, d\mid N, a\in (\mathbf Z/(d,N/d)\mathbf Z)^\times\}\subset \mathbf P^1(\mathbf Q)$ is a complete set of representatives of the cusps of $X$. The subspace $\mathbf Q \otimes U_\eta \subset \mathbf Q\otimes U$ coincides with $\mathbf Q\otimes U'$,
where $U'\subset U$ consists of modular units which have the same valuation at any two cusps $a/d$, $a'/d$ with the same denominator; cf. \cite[Prop. 2]{Gonzalez}. %  Indeed, an easy computation shows that, for each $d\mid N$:
% \begin{equation}\label{orderu}
% \mbox{the order of }\, u(q) = \prod_{d\mid N} \eta(q^d)^{r_d} \, \mbox{ at cusps  }\, a/d \, \mbox{ is }\, n_d/24,
% \end{equation} 
% from which one inclusion follows. The reverse inclusion also holds because $\mathrm{det}(A_N)\in \Q^\times$.
% For example, if $N$ is square-free it follows that $\mathbf Q \otimes U_\eta=\mathbf Q\otimes U$ because in this case $(d,N/d)=1$ for all $d\mid N$ and thus the set of cusps of $X$ is $\{ \frac{1}{d}:\,d\mid N\}$.%  But in general this restriction might impose non-trivial linear conditions on $\Q\otimes U$ and $\Q\otimes U_\eta $ might be properly contained in $\Q\otimes U$.
This implies the proposition in light of the next lemma, since an element  $u\in U\subset \mathbf Q(X)$ has the same valuation at any two Galois-conjugate cusps. 
\end{proof}
\begin{lemma}\label{}
  Let $d|N$. Then the cusp $1/d$ is rational if and only if $(d,N/d)=1$.
More generally, the Galois orbit of the cusp $1/d$ is the set of cusps
$a/d$ with $a$ relatively prime to $(d, N/d)$.
\end{lemma}
\begin{proof}[Proof of the lemma]
We prove the first statement using the results of \cite[\S 1.3]{Stevens}.
Namely, it is known that the cusps of $X$ are rational over $\mathbf Q(\zeta_N)$, and the Galois action can be described explicitly as follows \cite[Thm. 1.3.1]{Stevens}: if $\tau$ is the automorphism of $\mathbf Q(\zeta_N)$ which sends $\zeta_N\mapsto \zeta_N^n$ and $n'\in \mathbf Z$ is chosen so that $nn'\equiv 1\pmod N$ then $\tau$ sends the cusp $x/y$ to $x/n'y$. In particular, it follows straightforwardly that $[1/d]^\tau = [n'/d]$. 
Thus it suffices to prove that the cusps $n'/d$ and $1/d$ coincide for all $n'$ relatively prime to $N$, if and only if $(d, N/d)=1$. This fact can be shown by an elementary argument using the conditions for the integer  matrix sending $1/d$ to $n'/d$ to be in $\Gamma_0(N)$.
The second statement is proved similarly\footnote{FIXME: double check this!}.
\end{proof}

% The group $U_\eta$ is the sub-lattice in $\mathbf Q\otimes U_\eta$-vector space given by those functions whose $q$-expansion is of the form
% \begin{equation}\label{latticeUeta}
% u = \frac{1}{q^N}+\sum_{n>N} c_n q^n, \, c_n\in \Z.
% \end{equation}


% In any case, cusp $\infty$ is equivalent to $1/N$, and cusp $0$ is equivalent to $1$; since both $1$ (resp. $1/N$) is the single cusp with denominator $1$ (resp. $N$), a
By the Riemann-Roch theorem, there exist nonconstant rational functions on $X$ which are regular away from $\infty$. The proposition implies that an integer power of such a function belongs to the subgroup $U_\eta\subset U$, which yields the following.
\begin{corollary}\label{uexists}
There exists an eta product $u\in U_\eta$
which is regular away from $\infty$.\qed
\end{corollary}
It is thus possible to compute the rational function $u$ required in the compation of a basis for $H^1_{\mathrm{dR}}(X)$ as an eta product.
% to write down a basis for de Rham cohomology, a natural approach is to find a vector $(r_d)_{d|N}$ such that the corresponding eta product $u$  is of this form. By Lemma \ref{h1basis}, we do obtain a basis $\mathcal B=\{\omega_1,\ldots, \omega_{t}, u\omega_1,\ldots, u\omega_{t}\}$ in this manner if $\infty$ is non-Weierstrass and  $u$ has pole order $t+1$ at $\infty$. Generally this is not the case. (While the fact that $\infty$ is non-Weierstrass implies that there exists $w\in U$ regular away from $\infty$ and with pole order exactly $t+1$ at $\infty$, the inclusion $U_\eta\subset U$ may be strict, so it may be necessary to pass to an integer power $u = w^n$ in order to obtain a function that can be written down explicitly as an eta product.)
%  But by Remark \ref{oggremark} we can check whether $\mathcal B$ is a basis in any given case. 
%   Unfortunately, even for prime level, it not always the case that the optimal $\eta$-quotient $u$ -- that is, the $\eta$-quotient which is holomorphic away from $\infty$ with a pole of minimal order at $\infty$ -- has the smallest pole order $g(X)+1$ at $\infty$  consistent with Riemann-Roch. (While $U_\eta \otimes \mathbf Q = U\otimes \mathbf Q$ for $N$ squarefree, the inclusion $U_\eta\subset U$ thus seems to be strict most of the time, even for squarefree level.)
 A practical approach to finding the vector $(r_d)_{d|N}$ giving rise to the $u$ we seek is to apply a mixed-integer linear programming algorithm: one minimizes the pole order $-n_N$ of $u$ at $\infty$ subject to the criteria of Newman-Ligozat in Definition \ref{etadef} and the condition that the orders $n_d$ of $u$ at other cusps are non-negative. 
\begin{remarkwr}\label{}
By minimizing the pole order of $u$ at $\infty$, we may compute using the method of the previous subsection  a basis $\mathcal B$ for the de Rham cohomology of $X$ to a desired degree of precision using as few Fourier coefficients as possible for the cusp forms $\omega_1,\ldots,\omega_t$.
It is desirable that this minimal pole order $-n_N$ equal $t+1$.
This condition is relevant for the computation of
$\alpha_{\omega_{g,i},\eta_{g,i}}$ (cf. \S \ref{compute}), as well as to apply Lemma \ref{h1basis} from the previous subsection.
Unfortunately it does not always hold; see the discussion in \S \ref{compute}.
\end{remarkwr}
\begin{remarkwr}\label{ueff}
To determine the complexity of the algorithm described in this paper (see \S \ref{complexity}), it is necessary to bound effectively (as a function of $N$) the order of the pole at $\infty$ of the eta quotient $u$ in Corollary \ref{uexists}. This can be done by examining the proof of Corollary \ref{uexists}.
By the Riemann-Roch theorem, there is a nonconstant rational function $w$ on $X$ which is regular one $Y$
and has a pole of order $\leq t = \text{genus}(X)$ at $\infty$.
From the formula for $\text{genus}(X)$ as a function on $N$, one can thus extract the bound $-\ord_\infty(w)= O(N\log\log N)$; cf. \cite{genera}. (We adopt the convention that unless decorated with a subscript, an expression $O(-)$ denotes
a bound with an absolute implied constant.)
The proof of \cite[Prop. 2]{Gonzalez}, which was invoked to show Proposition \ref{QUeta}, 
shows that $w^\mu$ belongs to $U_\eta$ for an integer $\mu = O(\det A_N)$.
Combining this with the explicit formula \cite[Prop. 1]{Gonzalez} for $\det A_N$
readily yields the estimate $-\ord_\infty(u) = O(e^{C(\log N)^2})$
for an absolute constant $C$.
\end{remarkwr}


 \subsection{Computing the Poincar\'e dual $\gamma_f$ of $\omega_f$}
\label{computepoincaredual}

Assume that $\{\gamma^{(j)}\}$ is a collection of elements of $\Gamma_0(N)$
with small lower-left entries $cN$, whose homology classes $\gamma_0^{(j)}$ generate $H_1(X^{\mathrm{an}},\mathbf Z)$. By a brute-force search it is straightforward to find such elements $\gamma^{(j)}$ in practice. (For small $N$, often one need take $c$ no greater than 2 or 3.)
% Expressing each $m_i,m_i'$ above as a $\mathbf Z$-linear combination of the $s_j$'s, we can thus efficiently compute the periods $I(\gamma_i,\omega_j)$
% and $I(\gamma_i',\omega_j)$ as linear combinations of periods $I(\sigma_i,\omega_j)$ which can be calculated using relatively few Fourier coefficients for $\omega_j$.




For any $m\in H_1(X^{\mathrm{an}},\mathbf C)$, write $\eta_m$ for the Poincar\'e dual of $c$; conversely, for any differential $\eta$ of the second kind on $X$, let $m_\eta\in H_1(X^{\mathrm{an}},\mathbf C)$ denote the Poincar\'e dual
of its cohomology class. We normalize the Poincar\'e duality isomorphism so that it is characterized by the property
\begin{equation}\label{pd}
\langle \eta_m,\eta\rangle  = \int_m \eta.
\end{equation} 

% (In \cite[Ch. III, (1.1.1)]{FK} they rather agree to say that "Poincar\'e dual" means $\int_m \eta = \int_{X_0(N)} \eta_m\wedge \eta$; we use a different convention, which differs from \cite{FK} by the factor $2\pi i$. If we stick to this, we should at least be consistent all the time. Check throughout that we didn't put a wrong factor $2\pi i$ anywhere.)

The vector space $H_1(X^{\mathrm{an}},\mathbf C)$ is also equipped an intersection product, which is related to the Poincar\'e pairing by the formula %(cf.\, \cite[Ch. III, (1.1.2)]{FK})
\begin{equation*}
m\cdot m_\eta = \frac{1}{2 \pi i}\langle \eta_m,\eta\rangle.
\end{equation*}
The homology of $X$ also admits a natural action of the Hecke algebra, compatible with the action on cohomology via Poincar\'e duality. For any normalized eigenform $f\in S_2(\Gamma_0(N))$ and any $m\in H_1(X^{\mathrm{an}},\mathbf C)$, write $m^f\in H_1(X^{\mathrm{an}},\mathbf C)[f]$ for the projection of $f$ onto the $f$-isotypic component of homology.
Simiarly, for $\eta\in H^1_{\mathrm{dR}}(X/\mathbf Q)$ write
$\eta^f$ for its projection onto the $f$-isotypic component.

We can assume that via the method described above a symplectic basis
$$\mathcal S = \{\omega_{f,1},\ldots, \omega_{f, n}, \eta_{f,1},\ldots, \eta_{f, n}\}$$ for $H^1_{\mathrm{dR}}(X/\mathbf Q)[f]$ has already been computed.

\begin{lemma}\label{fact} Fix $\gamma_1$, $\gamma_2\in \Gamma_0(N)$ and let $m_1$, $m_2\in H_1(X^{\mathrm{an}},\mathbf Z)$ denote the corresponding homology classes on $X$.%, induced by any path in $\mathfrak h$ joining any base point $\tau_0\in \mathcal H$ with $\gamma_i \tau_0$.
For any normalized eigenform $f\in S_2(\Gamma_0(N))$, we have
$$m_1^f\cdot m^f_2=\frac{1}{2\pi i}\sum_{i=1}^{n} I(\omega_{f,i};m_1)I(\eta_{f,i};m_2) - I(\omega_{f,i};m_2)I(\eta_{f,i};m_1),$$
where $\omega_f  = 2\pi i f(z)\mathrm{d}z$ is the 1-form corresponding to $f$.
\end{lemma}

\begin{proof}
Let $\eta_k=\eta_{m_k}$ and write $\eta^{f}_k=\sum c^{(k)}_i\omega_{f,i} + \sum d^{(k)}_i\eta_{f,i}$. Then we compute $m^f_1\cdot m^f_2 =  \frac{1}{2 \pi i} \langle \eta^f_1,\eta^f_2\rangle = \sum_i \frac{1}{2 \pi i} (c^{(1)}_i d^{(2)}_i-c^{(2)}_i d^{(1)}_i )= \frac{1}{2 \pi i} \sum_i (I(\omega_{f,i};\gamma_1)I(\eta_{f,i},\gamma_2)-I(\eta_{f,i};\gamma_1)I(\omega_{f,i};\gamma_2) )$. 
\end{proof}


% Assume from now on that $f=\sum_{n=1}^\infty a_n q^n$ has {\em integral} fourier coefficients and let $E/\Q$ be the elliptic curve corresponding to $f$ by the construction of Eichler-Shimura. 
% Let us also assume that the sign in the functional equation of $L(f,s)$ is $-1$, so that $E(\Q)$ is expected to have odd rank.


Now assume $f$ is the newform with rational Fourier coefficients which parametrizes the elliptic curve $E$, and as above denote by $\omega_f$ the corresponding holomorphic 1-form. 
Then the $f$-isotypic components of the homology and cohomology of $X$ are two-dimensional.
Write $\eta_f$ for the ``complementary'' form of the second kind such that $\{\omega_f,\eta_f\}$ is a symplectic basis for $H^1_{\mathrm{dR}}(X/\mathbf Q)$; in particular, $\langle \omega_f, \eta_f\rangle = 1$.
 Let $\gamma_f^+$ (resp. $\gamma_f^-$) be a generator of the plus (resp. minus) eigenspace of
the $f$-isotypic component of $H_1(X^{\mathrm{an}},\mathbf Z)$ under the action of complex conjugation; these are unique up to sign. Since the splitting into plus and minus subspaces only takes place over $\mathbf Z[\frac{1}{2}]$, we have $\gamma_f^+\cdot \gamma_f^-\in \pm 2^{\mathbf{N}}$; after adjusting the signs if necessary, we can assume that $\gamma_f^+\cdot \gamma_f^-=2^n$ for some $n\geq 0$. This determines the pair $(\gamma_f^+, \gamma_f^-)$ up to a sign.

Let $ \Omega_f^\pm = I(\omega_f,\gamma_f^\pm).$
Note that $\Omega_f^+$ and $\Omega_f^-$ generate a lattice $\mathbf Z \Omega^+  + \mathbf Z \Omega^-\subset \mathbf C$ which is independent of the choice of the sign above,
and is commensurable with the N\'eron lattice $\Lambda_f$ of $E$. 
Finally, let 
$$ \gamma_f  := \frac{1}{2^{n+1}\pi i} (\Omega_f^+\gamma_f^- - \Omega_f^- \gamma_f^+)\quad \in H_1(X^{\mathrm{an}},\mathbf C)[f].$$

\begin{proposition}\label{dual}
The homology class $\gamma_f$ is Poincar\'e dual to the cohomology class of $\omega_f$.
\end{proposition}

\begin{proof}
We only need to check \eqref{pd} for $\eta=\eta_f$, and $I(\eta_f; \gamma_f) =  1 = \langle \omega_f,\eta_f\rangle$
by \ref{fact}.
\end{proof}
The preceding discussion reduces the computation of $\gamma_f$ 
to finding the homology classes  $\gamma_f^\pm$,
from which the periods $\Omega_f^\pm$ are readily obtained by integrating.
The classes $\gamma_f^\pm$ can be calculated using modular symbols, for which we refer to \cite{S2},
and then expressed in terms of the ``good'' basis $\{\gamma^{(j)}_0\}$.
% The calculation of a basis of the first homology group $H_1(X^{\mathrm{an}},\mathbf Z)$ and the space $H^1(X,\Omega^1)$ of regular differential forms on $X$,  and the action
% of Hecke operators on these two spaces can be carried out using the standard modular
% symbol algorithms. (FIXME: add a reference to Cremona or Stein)
%  Hence, we can assume we have at our disposal 
% \begin{enumerate}
% \item[(i)] a basis $\omega_1,\ldots,\omega_{g(X)}$ of regular differentials of $X$ and 

% \item[(ii)] a collection $\{ \gamma_1,\gamma_1',...,\gamma_{g(X)},\gamma_{g(X)}' \}\subset \Gamma_0(N)$ of hyperbolic elements such that, if we choose a base point $\tau_0\in \mathcal H$ and write $m_i$ (resp. $m'_i$) for the homology class of the projection to $X$ of any continuous path in $\mathcal H$ from $\tau$ to $\gamma_i \tau_0$ (resp. $\gamma'_i\tau_0$), then 
% \begin{equation}
% H_1(X,\Z) = \langle m_1,m'_1,...,m_{g(X)},m'_{g(X)}\rangle_\Z.
% \end{equation}
% \end{enumerate}

% In particular we are able to compute the periods 
% \begin{equation}
% I(\gamma_i,\omega_j)=I(m_i,\omega_j), \quad I(\gamma'_i,\omega_j)=I(m'_i,\omega_j),
% \end{equation}
% and the above two basis can be found such that $I(\gamma_i,\omega_j)=I(\gamma'_i,\omega_j) = 0$ for $i\ne j$




% \subsection{Computation of iterated integrals}
%  We now turn to analyze the algorithm to compute the element 
% $$
% z_{g,f} := \int_{\gamma_f} (\omega_g\cdot \eta_g -\eta_g\cdot \omega_g).
% $$ 

% Notice that the difference between the points $P_{g,f}$ and $2 z_{g,f}$ is just the line integral $\int_{\gamma_f} \alpha = \langle \omega_f,\alpha\rangle$, which is {\em rational number} because of the definition \eqref{poincare} we gave of the Poincar\'e pairing. In the next subsection we will discuss two methods for computing this rational adjustment factor.



% In other words, 
% $$
% [z_{g,f}]=[P_{g,f}] \in \C/\langle \Lambda_f, \,\Q\rangle. 
% $$

% Hence in practice we can often exploit the following strategy: 

% \begin{itemize}


% \item The element $\gamma_f\in \C[\Gamma_0(N)]$ and the differentials of the second kind $\omega_g$, $\eta_g$ have been already computed in points (I) and (II), respectively. In (I) it was also explained how to compute $\gamma_f$ and how to choose an optimal base point for integrating. Be careful though when computing in practice these iterated integrals and choosing the base points... The way we do it will have to be explained with care in the final version of this article.

% \end{itemize}

%  Write $\tilde{\eta}_g$ for the primitive of $\eta_g$ on $\mathcal H$.
% FIXME: Explain what we mean by the integral
% $\int_{\gamma_f} \omega_g \tilde{\eta}_g -\eta_g \tilde{\omega}_g$
% and how to compute it relatively efficiently using the Fourier expansions of $\omega_g$ and $\eta_g$
% and
% the coefficients of $\gamma_f\in (\bigoplus \mathbf Z[\gamma_i]\oplus \mathbf Z[\gamma_i])\otimes \mathbf C$, where
% the $\gamma_i$ and $\gamma_i'$ comprise our ``nice'' collection
% of hyperbolic elements of $\Gamma_0(N)$ whose images generate the integral homology of $X$. 



\subsection{Computing the adjustments $\int_{\gamma_f}\alpha$}\label{computealpha}
At this point we are already able to compute
\[z_{g,f}:= \sum_{i=1}^k \int_{\gamma_f}(\omega_{g,i}\cdot\eta_{g,i}-\eta_{g,i}\cdot \omega_{g,i}),\]
part of the righthand side of formula \eqref{pgf}.

We describe two approaches for computing the difference $\Delta_{g,f}=P_{g,f}-z_{g,f}$. 
\subsubsection{}\label{compute}
The direct approach is to
compute $q$-expansions for each 1-form $\alpha = \alpha_{\omega_{g,i},\eta_{g,i}}$ explicitly. We do not know how to do this in general, but the method we now describe works under the following hypothesis:
\begin{equation*}\label{hyp}\tag{$\dagger$}
\begin{array}{l}\text{The point $\infty$ is not a Weierstrass point, and the optimized eta product $u$ of \S \ref{eta}} \\
\text{has a pole of order 
$t+1 = p_a(X)+1$  at $\infty$.} \end{array}
\end{equation*}
Assume \eqref{hyp} holds. Then the meromorphic  differentials $u\omega$ as $\omega\in H^0(X,\Omega^1_X)$ varies are all of the second kind, regular on $Y$, and have poles of all orders $2,3,\ldots, t+1$ at $\infty$. 

Recall that the defining property of $\alpha$ is that its principal part at $\infty$ agrees with that of  $\xi = \omega_{g,i} F_{\eta_{g,i}}$ on $\tilde X$, modulo $\mathrm{d}q/q$; i.e., $\xi-\alpha$ has at worst logarithmic poles. % So provided
% we can compute a 1-form $\lambda$ of the second kind on $X$, regular on $Y$, which 
% has the same \textit{leading term} as a given differential such as $\xi$,
% it is possible to compute $\alpha$. Indeed, simply build a suitable linear combination of these 1-forms $\lambda$ which matches the principal part of $\xi$ except for the $\mathrm{d}q/q$ term.
The symplectic basis $\{\omega_{g,i},\eta_{g,i}\}$ for $H^1_{\mathrm{dR}}(X/\mathbf Q)[g]$ is obtained by applying linear operations to a basis $\mathcal B$ consisting of 1-forms which are either holomorphic or of the type $u\omega$ for holomorphic $\omega$. In particular, $\xi$ has poles of order $\leq t$
at the points lying over $\infty$.
Thus there exists holomorphic $\omega$ such that  $u\omega$ has exactly the same principal part as $\xi$, modulo $\mathrm{d}q/q$. So assuming \eqref{hyp}, $\alpha$ can be computed explicitly knowing only Fourier expansions of $u$ and of a basis for $S_2(\Gamma_0(N))$. 
The resulting explicit Laurent expansions for $\alpha$ can then be integrated over $\gamma_f$ using the same approach discussed above for computing $z_{g,f}$, to find an approximation to the rational number $\Delta_{g,f} = P_{g,f}-z_{g,f}$. See \S \ref{37} for an example computation of this sort.

Unfortunately, \eqref{hyp} is not always satisfied, but it does hold in some situations. For example, using the results of Newman and Lizogat combined with the formula for the genus of $X_0(N)$, it is an easy exercise to show that \eqref{hyp} holds if $N=p$ for a prime $p \equiv 1$ (mod 12), $N=2q$ for a prime $q \equiv 1$ (mod 4), or $N=3r$ for a prime $r \equiv 1$ (mod 3). On the other hand, it never holds if $N=p$ for a prime $p\not\equiv 1$ (mod 12).



\subsubsection{}\label{guess}
When \eqref{hyp} fails, an alternative approach is required. One such method is to make an educated guess as to the value of $\Delta_{g,f}$. This also has the advantage of avoiding computationally expensive integral evaluations.
Note that $\Delta_{g,f}=\sum_i\langle \omega_f, \alpha_{\omega_{g, i},\eta_{g, i}}\rangle$ is a rational number because $\omega_f$ and each $\alpha$  are algebraic and defined over $\mathbf Q$. The method of lattice reduction is well-suited to guessing the value of an unknown rational number.

 We now sketch how this might be done. First compute the elliptic logarithm $Q\in \mathbf C/\Lambda_f$ of a generator of the Mordell-Weil group $E(\mathbf Q)$. (This can be done in various ways. We remark that the interest of computing the Chow-Heegner points $P_{g,f}$ is not as a tool for computing $E(\mathbf Q)$, so appealing to an independent algorithm to obtain $Q$ is not circular reasoning.)
Next find a basis $\{b_1,b_2\}$ for $\Lambda_f$.
Since $d_g P_{g,f}$ corresponds a rational point of $E$, some integer multiple of it
must be a $\mathbf Z$-linear combination of $b_1,b_2$, and $Q$.
Using LLL or another lattice reduction algorithm, find an approximate  dependence relation
\[Dd_gz_{g,f} + A_1 b_1 + A_2 b_2 + NQ + M = 0,  \qquad A_1, A_2, D, M, N \in \mathbf Z.\]
There will be such a relation with $M/Dd_g = \Delta_{g,f}$, indicating that up to $D$-torsion in $E$, $d_gP_{g,f}$ maps to $N$ times the chosen generator of $E(\mathbf Q)$. 

Unfortunately it is not easy in practice to compute $\Delta_{g,f}$ in this manner. The problem is that a prohibitively large degree of accuracy is usually necessary to identify the ``correct'' dependence relation as above, since in general the rational number $\Delta_{g,f}$ may have fairly large height. See \S \ref{43} for an example.



\subsection{Computing the denominator $d_g$}
The final ingredient to be computed is the denominator $d_g$, or the smallest positive integer such that $d_gT_g \in \mathbf T_0$. This can be accomplished by computing the Hecke algebra in {\tt SAGE}. Under the inclusion $\mathbf T_0 \hookrightarrow \mathbf T_0 \otimes \mathbf Q \cong \prod_g K_g$, the operator $T_n$ is sent to the vector $(a_n(g))_g$. Since each eigenvalue is an algebraic integer, then the image lies in $\prod_g \mathcal{O}_{K_g}$. Therefore, $\mathbf T_0$ can be embedded (as a $\mathbf Z$-module) in $\mathbf Z^t$, where $t$ is the genus of $X_0(N)$. The image of any Hecke operator $T_n$ in $\mathbf Z^t$ can be computed by finding the image of $a_n(g)$ in $\mathcal{O}_{K_g}$ with respect to an integral basis.

It is well known that $\mathbf T_0$ is generated as a $\mathbf Z$-module by $\{T_n\}_{1 \leq n \leq r, (n,N)=1}$, where $r = \lceil \frac{N}{6}\prod_{p \mid N} (1 + \frac{1}{p})\rceil$; see for instance \cite{S1}. Hence, the image of $\mathbf T_0$ as a submodule of $\mathbf Z^t$ can be computed by taking the submodule generated by the images of a finite number of Hecke operators. It is then a simple matter to find $d_g$ for each newform $g$.

\subsection{Remark on complexity}\label{complexity}
The complexity of the computations we have described
is primarily determined by
 the {\em number} $n_D$ of Fourier coefficients required
to compute $z_{g,f}$ (and also the correction $\Delta_{g,f}$, if using the method of \S \ref{compute}) to a given number $D$ of digits of accuracy.
In this subsection we sketch a method for obtaining a bound on $n_D$ in terms of $N$.

\subsubsection{}
Write the Fourier expansion of $u$ as 
$$u(\tau) = \sum_{n\geq -n_0} c_n q^n, \qquad q = e^{2\pi i \tau}, \tau \in \mathfrak H.$$ 
Let the principal part of $u$ at $\infty$ be
\[\mathrm{pp}_\infty (u) = \sum_{-n_0 \leq -m \leq 0} \frac{d_m}{m} q^{-m},\qquad d_m = mc_{-m}.\]
In \cite{ono}, Bringmann and Ono prove an exact formula for the Fourier coefficients of harmonic Maass forms, of which weakly holomorphic modular functions such as $u$ are examples. To avoid introducing unnecessary notation, we state only the very special case of their result applicable to our situation. We remark that long ago Rademacher used the circle method to prove a similar exact formula for the coefficients of the $j$-function \cite{rad}, and a modification his argument would probably yield a simpler and more direct proof of the special case we require.
Using \cite[Thm. 1.1]{ono}, one can express the coefficients $c_n, n>0$ in terms of the coefficents $d_m$,
 the order-1 $I$-Bessel function $I_1(z)$,  and the Kloosterman sum
\[K(-m,n,c):=\sum_{\substack{0<v< c\\ (v,c)=1}} \exp\left(\frac{2\pi i}{c}(nv + m\bar v)\right),\]
where $\bar v$ denotes the multiplicative inverse of $v$ modulo $c$.
Namely, \textit{loc. cit.} yields the formula
\begin{equation}\label{coeff}
  c_n = 2\pi  \sum_{-n_0\leq -m \leq 0} d_m \sum_{\substack{c>0\\ c\equiv 0\pmod N}} (\tfrac{m}{n})^{1/2} \frac{K_0(-m,n,c)}{c}I_1\left(\frac{4\pi\sqrt{|mn|}}{c}\right),\qquad n>0.
\end{equation}
By Remark \ref{ueff}, we have 
\begin{equation}\label{n0est}
  n_0 = -\ord_\infty(u) = O(e^{C(\log N)^2}).
\end{equation}
 for an absolute constant $C$.
We can trivially bound the numbers $d_m$ as follows.
Let $\xi_r(x) = re^{2\pi i x}$ for $0<r<1$
and set $y=-\frac{1}{2\pi}\log r>0$.
The Cauchy integral formula applied to the meromorphic function $U(q) = \sum c_n q^n$
of $q$ in the unit disk gives
\[\frac{d_m}{m} = \frac{1}{2\pi i} \int_{\xi_r} \frac{U(q)}{q^{m+1}}\operatorname{d}q = e^{2\pi m y} \int_0^1 
u(x+iy)e^{-2\pi inx}\operatorname{d}x.\]
Taking $y=1$, say, we thus have
\[|d_m| \leq m e^{2\pi m}\int_0^1 |u(x+i)|\operatorname{d}x
= m e^{2\pi i m} \int_0^1 \prod_{d|N} |\eta(dx+id)|^{r_d}\operatorname{d}x,\]
where $\eta$ is the Dedekind eta function and $(r_d)_{d|N}$ is the vector giving rise to $u$.
If $B$ is any absolute bound for the holomorphic function $\eta$
on $\mathfrak H^* \cap \{\mathrm{Im}(z) \geq 1-\epsilon\}$,
then we obtain the estimate (recall Definition \ref{etadef})
\begin{equation}\label{dest}
  |d_m|\leq me^{2\pi m} B^{\sum r_d} = m e^{2\pi m} \leq n_0 e^{2\pi n_0}.
\end{equation}
From \eqref{coeff} and \eqref{dest}, standard estimates for Kloosterman sums, and asymptotics for the $I$-Bessel function, one obtains by the method of \cite[\S \S 5.1-2]{BrPh} the estimate
\[c_n = O(N^{5/4}n_0^{7/4} n^{-3/4} \exp(2\pi n_0 + N^{-1}4\pi \sqrt{nn_0})).\]
In light of \eqref{n0est} this yields
\begin{equation}\label{cnest}
  c_n \leq A n^{-3/4} e^{B N^{C\log N}\sqrt n},
\end{equation}
for some absolute constants $A,B,C$.


\subsubsection{} 
The coefficients $c_n$ determine the Fourier coefficients of the 1-forms occuring in the formula \eqref{pgf}
for $P_{g,f}$. Unfortunately the relationship is indirect, as the construction of the 1-forms $\omega_{g,i},\eta_{g,i},
\alpha_{\omega_{g,i},\eta_{g,i}}$ involves multiplying $u$ against a basis of cusp forms for $\Gamma_0(N)$
and then performing a lot of linear algebra. By Deligne's proof of the Ramanujan-Petersson conjecture,
the cusp forms have $n$th coefficient of size $O_\epsilon(n^{\frac{1}{2}+\epsilon})$.
It follows that $n$th Fourier coefficient of an element of the basis $\mathcal B$
for $H^1_{\mathrm{dR}}(X/\mathbf Q)$ computed in \S \ref{symplecticbasis} has size
\[O(P(n) e^{BN^{C\log N}\sqrt n}),\]
for absolute constants $B,C$ and a universal polynomial $P(n)$.
To compute the 1-forms $\omega_{g,i},\eta_{g,i},\alpha$, linear algebra operations
are performed on this basis, which spans a vector space of dimension $\text{genus}(X) = O(N\log\log N)$.
It thus seems likely that a careful analysis of the linear algebra operations performed
would yield a bound 
\begin{equation}\label{final}
  O(Q(n,N \log N) e^{B N^{C\log N}\sqrt n})
\end{equation}
for the $n$th Fourier coefficient of \textit{any} 1-form integrated in the course of computing \eqref{pgf}.
Here $B,C$ are absolute constants and $Q(X,Y)$ is a universal polynomial independent of $N$.

Suppose $\lambda$ is such a 1-form (on $X$ or $\tilde X$),
and consider the problem of integrating the pullback of $\lambda$ to $\mathfrak H^*$
along a path from $\tau_1$ to $\tau_2$. 
By the method explained in \S \ref{integration},
we can assume that $\text{im} \tau_1,\text{im}\tau_2 \geq (c_{\mathrm{max}, N}\cdot N)^{-1}$,
where $c_{\mathrm{max}, N}\cdot N$ is the largest of the lower-left entries
of the elements $\gamma^{(j)}\in \Gamma_0(N)$ introduced at the beginning of \S \ref{computepoincaredual}.
Recall that these consisted of a collection of elements which span $H_1(X^\mathrm{an},\mathbf Z)$
and have lower-left entries as small as possible. We do not know how to bound $c_{\mathrm{max},N}$ in terms of $N$,
although in practice it seems to be very small ($c_{\mathrm{max},N}\leq 3$ for $N\leq 500$)\footnote{FIXME: I just put this in as a rough guess; go back and add the right numbers.}

If the Laurent expansion for $\lambda$ about $\infty$ (or a lift of $\infty$ to $\tilde X$)
is $\lambda = \sum a_\lambda(n) \frac{\operatorname{d}q}{q}$,
then setting $\tau_j = x_j + iy_j$ for $j=1,2$ (where $y_j \geq (c_{\max, N} N)^{-1}$), we have
\[\int_{\tau_1}^{\tau_2} \lambda = \sum_{n\gg -\infty} \frac{a_\lambda(n)}{n} (e^{2\pi i n x_2}e^{-2\pi n y_2}
- e^{2\pi i n x_1}e^{-2\pi n y_1}).\]
Our problem is to determine $n_D$
such that the tails of these sums are bounded by the requisite precision, say $10^{-D}$.
It clearly suffices to ensure that
\[S(n_D):=\sum_{n\geq n_D} n^{-1} |a_\lambda(n)| e^{-2\pi n y} \leq 10^{-D}.\]
Granting \eqref{final}, we have
\[S(n_D) \ll \sum_{n\geq n_D} n^{-1} Q(n, N\log N) e^{B N^{C\log N} \sqrt n - 2\pi n / c_{\max, N} N}
\ll \sum_{n\geq n_D} e^{-2\pi n/ Nc_{\max, N}} =  \frac{e^{-2\pi n_D/ N c_{\max, N}}}{1-e^{-2\pi/Nc_{\max, N}}}.\]
This shows that provided \eqref{final} holds, we have the following estimate for $n_D$ in terms of $D$ and $N$:
\begin{equation}\label{nDest}
  n_D = O(DNc_{\max, N}),
\end{equation}
where the implied constant is absolute.
As remarked earlier, for practical purposes it seems that  $c_{\max,N}$  can be treated as a constant.
\section{Numerical examples}
\subsection{Example: {\tt 37a1}}\label{37}

 Take $N=37$ in the setup of our algorithm. In this setting, the space of regular differentials on $X=X_0(37)$ is spanned by $\omega_f$ and $\omega_g$,
 which are associated to elliptic curves over $\mathbf Q$ (labeled {\tt37a1} and {\tt 37b1} in Cremona's database) of ranks $1$ and $0$, respectively.

 By computing the periods attached to $\omega_f$ and $\omega_g$, it can be checked that the
  classes of the matrices
%  $$
%  \gamma_1 = \left(\begin{array}{cc} 12 & -1 \\ 37 & -3 \end{array}\right), \quad
%  \gamma_2 = \left(\begin{array}{cc} -3 & -1 \\ 37 & 12 \end{array}\right), \quad
%  \gamma_3 = \left(\begin{array}{cc} 15 & 2 \\ 37 & 5 \end{array}\right), \quad
%  \gamma_4 = \left(\begin{array}{cc} 5 & -1 \\ -74 & 15 \end{array}\right)$$


%  The same can be said about the matrices
 $$
 \gamma_1=\begin{pmatrix}2&-1\\37 &-18\end{pmatrix}, \,\gamma_2= \begin{pmatrix}3&-1\\37 &-12\end{pmatrix}, \,\gamma_3= \begin{pmatrix}5&2\\37&15\end{pmatrix} ,\, \gamma_4=
 \begin{pmatrix}14 & 3\\ 37 &  8\end{pmatrix} 
 $$
  generate the rational homology of $X$.
These are a ``nice'' basis for the homology in the sense that the lower left entries are exactly 37 (rather than $37c$ for $|c|>1$), so the integral $\int_\tau^{\gamma_i\tau}\lambda$ can be evaluated efficiently for any meromorphic differential 1-form on $\mathbf H$ or $X_0(37)$ which is regular away form $\infty$,
using relatively few Fourier coefficients of $\lambda$.

 If we denote by $[\gamma]$ the homology class attached to the group element $\gamma$,
 then
 \begin{eqnarray*}
 \gamma_f^{+}&=&\quad\quad\frac{-1}{2}[\gamma_2] +\frac{1}{2}[\gamma_3] -\frac{1}{2}[\gamma_4] 
 \\
 \gamma_f^{-}&=&[\gamma_1] -2[\gamma_2]
 \end{eqnarray*} 
 generate the $f$-isotypic part of the integral homology of $X$. The superscripts 
 indicate the eigenvalue of complex conjugation acting on the homology class.

 To obtain differentials of the second kind representing classes in the deRham cohomology, we consider the elements of the
 form 
 $$ \eta_1 = u\cdot \omega_f, \qquad \eta_2 = u\cdot \omega_g,\quad \qquad \mbox{ where } 
 u = q^{-3} \prod_{n=1}^\infty (1-q^n)^2 (1-q^{37n})^{-2}.$$

 It is not hard to check (by calculating the periods along $\gamma_f^\pm$ and $\gamma_g^\pm$) 
 that the classes of $\omega_f, \omega_g, \eta_1$ and $\eta_2$ generate the deRham cohomology
 of $X$. Furthermore, by finding the matrix $M$ of the Hecke operator $T_2$
acting on $H^1_{\dR}(X_0(37))$ with respect to the basis $\omega_f,\omega_g,\eta_1,\eta_2$, and then determining the eigenspaces of $M$, one finds that
 \begin{eqnarray*}  
 \eta_f &=&  \frac{1}{4}(-37\omega_g+4\eta_1 -8  \eta_2), \\
 \eta_g &=&  \frac{1}{4}(37\omega_f-6\eta_1 + 10 \eta_2)
 \end{eqnarray*}
 are in  the $f$ and $g$ isotypic parts of the deRham cohomology respectively, and  $$ \langle \omega_f,\eta_f\rangle = \langle \omega_g,\eta_g\rangle = 1.$$

When one computes the Poincar\'e dual $\gamma_f$ of $\omega_f$,
one finds (with our normalization) that it is 
\[\gamma_f = \frac{1}{2\pi i}\left(\Omega_E^-\left([\gamma_2]-[\gamma_3]+[\gamma_4]\right)-\Omega_E^+\left(-[\gamma_1]+2[\gamma_2]\right)\right).\]
Here 
\[\Omega_E^+ \approx 2.993458646...,\qquad \Omega_E^-\approx (2.45138938...)i\]
are the real and imaginary periods of the elliptic curve $E$
corresponding to $f$ (labeled '37a1' in Cremona's database).

By computing with principal parts, one finds that
$\mathrm{pp}_\infty(2\omega_gF_{\eta_g})
= \mathrm{pp}_\infty(\frac{1}{2}(\eta_1-\eta_2))\mod \frac{\mathrm{d}q}{q}$.
Thus $\frac{1}{2}(\eta_1-\eta_2) = \alpha_{\omega_g;\eta_g}$
and integrating this over $\gamma_f$ yields $-1$ (to many digits of precision).

Computing the iterated integral $\int_{\gamma_f}(\omega_g\cdot\eta_g-\eta_g\cdot\omega_g)$ via the algorithm we have described yields the ``raw'' point
\[z_{g,f} = -1.40936100075... + (1.22569469099...)i.\]
Thus $P_{g,f} = z_{g,f}-\int_\gamma\alpha_{\omega_g;\eta_g}
= 0.4093610075... + (1.22569469099...)i$.

Now $E(\mathbf Q)$ is generated by the point $p = (0:-1:1)\in \mathbf P^2(\mathbf Q)$. The elliptic logarithm of this point in $\mathbf C/\Lambda_E$
is $P \approx 2.06386593094656...+(1.22569469099340...)i$.

By a lattice reduction algorithm one easily finds the linear dependence relation
\[2P_{g,f} - 8\Omega_E^+ - 7\Omega_E^- + 12 P \approx 0\]
holds to at least 15 digits of accuracy. (In this example, all of the iterated integrals have been computed using 350 Fourier coefficients, and on a laptop comuter the entire computation finished in a matter of seconds.)
This says that the image of $P_{g,f}$ in $E(\mathbf C)$
is equal to $-6p = (6:14:1)$ modulo an irrational 2-torsion point.
To explain the denominator ``2'' that has occured, 
let $\mathbf T_0$ denote the ``anemic'' Hecke algebra generated over $\mathbf Z$
by the Hecke operators $T_p$ for $p\neq 37$.
Then one can compute
using the first few Fourier coefficients of $f$ and $g$
that  the idempotent $e=(0,1)\in \mathbf Q\times \mathbf Q \stackrel{(\star)}{\cong} \mathbf T\otimes\mathbf Q$ 
does not belong to $\mathbf T\subset \mathbf T\otimes\mathbf Q$
but $2e$ does so. (The identification $(\star)$
associates $T_p\otimes 1\in \mathbf T\otimes \mathbf Q$
to $(a_p(f),a_p(g))\in \mathbf Q\times \mathbf Q$.)
Thus in fact the point we have denoted $P_{g,f}$ is the wrong thing;
the true triple Chow-Heegner point on $E$ associated to $(g,g,f)$
is $2P_{g,f}$ which is the rational point $-12P = (1357/841:28888/24389:1)$.
\subsection{Example: {\tt 43a1}}\label{43}
Let $N=43$ and let $E$ be the elliptic curve labeled {\tt 43a1} in Cremona's database. The modular curve $X=X_0(43)$ has genus 3. There are two isotypic components of $H^1_{\dR}(X)$, one of dimension 2 corresponding to the modular form $f$ which parametrized $E$,
and another of dimension 4 corresponding to a newform $g$
with Fourier coefficients in $\mathbf Q(\sqrt2)$, associated to an abelian surface quotient of $J_0(43)$.

In this case, a linear programming algorithm identifies the eta-quotient $u$
which is modular for $\Gamma_0(47)$ of weight 0, holomorphic away from the cusp $\infty$, and with minimal pole order at $\infty$, as
\[u = \frac{\eta(q)^{4}}{\eta(q^{43})^{4}} =  q^{-7}-4q^{-6}+2q^{-5}+8q^{-4}-5q^{-3}-4q^{-2}-10q^{-1}+8+9q+14q^3+O(q^4).\]
Since this has a pole of order $7 > 3+1 = 4$ at $\infty$, $u$ is not optimal for the purpose of our computations.

Nonetheless, computing the residue pairing shows that for a basis of cuspforms
 with rational Fourier coefficients, corresponding to holomorphic 1-forms
$\omega_f, \omega_{g,1},\omega_{g,2}$ on $X$, the collection
\[\omega_f,\omega_{g,1},\omega_{g,2},u\omega_f,u\omega_{g,1},u\omega_{g,2}\]
forms a basis for $H^1_{\dR}(X/\mathbf Q)$.
By finding the matrices of a few Hecke operators with respect to this basis,
one can as in the case $N=37$ produce symplectic bases
\[\omega_f,\eta_f,\qquad \text{and}\qquad \omega_{g,1},\omega_{g,2},\eta_{g,1},\eta_{g,2}\]
for $H^1_{\dR}(X)[f]$ and $H^1_{\dR}(X)[g]$ respectively.

While we can compute the Poincar\'e dual $\gamma_f$ and the ``raw'' point
\[z_{g,f} = \sum_{i=1}^2 \int_{\gamma_f}(\omega_{g,i}\cdot\eta_{g,i}-\eta_{g,i}\cdot\omega_{g,i})\approx -1.1460154...+(2.726364836...)i,\]
the fact that $u$ has such a large pole at $\infty$
prevents us from being able to find the 1-forms $\alpha_{\omega_{g,i};\eta_{g,i}}$ on $\tilde Y$.
Nonetheless, we can use lattice reduction to try to see whether $z_{g,f}$ differs from (the elliptic logarithm of) a rational point of $E$
by an adjustment factor in $\mathbf Q$.
Actually, we should first scale $z_{g,f}$ by the ``denominator''
of the idempotent $e\in \mathbf T\otimes\mathbf Q$
which projects onto the $g$-isotypic component when viewed as an operator on $H^1(X)$; as in the case $N=37$ a computation with Hecke eigenvalues of $g$ and $f$
reveals that the denominator of $e$ is $2$.

The elliptic curve $E$ again has rank $1$, generated by the point $p=(0:-1:1)$,
with elliptic logarithm $P \approx 1.53155105...$.
Lattice reduction reveals the linear dependence (to at least 58 digits of accuracy, the previous computations having been done with 1200 Fourier coefficients):
\[(2z_{g,f}) + 5 \Omega_E^+ - 4 \Omega_E^- - 8 P + \frac{1847467}{1984785} \approx 0.\]

The fact that so many Fourier coefficients were necessary to obtain this relation reflects the fairly large prime factor occuring in the denominator
\[1984785 = 3 \times 5 \times 11 \times 23 \times 523.\]
When the differentials
$\eta_{g,i}$ are expanded with respect to a basis for $H^1_{\mathrm{dR}}(X/\mathbf Q)$ consisting of differentials holomorphic away from $\infty$ and with integral Laurent expansions about $\infty$ -- which thus have integral period over $\gamma_f$ -- the prime factors above arise in the denominators of the coefficients. One thus expects these primes to occur in the denominators of the forms $\alpha_{\omega_{g_i},\eta_{g_i}}$.
% \subsection{Example: {\tt 57a1}}
% This example illustrates that the algorithm works fine even when there are oldforms of level $N$.

% Set $N=57$.  In this case there are four isotypic components of $H^1_{\dR}(X_0(57))$, and $X_0(57)$ has genus 5. The first corresponds to the elliptic curve of rank 1 labeled {\tt 57a1} in Cremona's database;
% the second and third correspond to the other isogeny classes of elliptic curves of conductor 57,
% and the fourth is the 4-dimensional space of oldforms arising from level 19.
% The optimal eta-quotient $u$ has a pole of order 6, so we can compute all the $\alpha$-adjustment factors explicitly.
% Hecke algebra computations show
% that the projectors onto the second, third, and fourth isotypic components
% of $\mathbf T_{\mathbf Q}$ have denominators 24, 3, and 4 respectively.
% The three points $P_i = P_{g_i;f}$ corresponding to the $g_i$-isotypic components for $g_i$ distinct from $f$ (when rescaled appropriately by the aforementioned denominators) are easily seen via lattice reduction to be 
% \[P_1 \equiv 32 P \mod \Lambda_E,\]
% \[P_2 \equiv P_3 \equiv -16P \mod \Lambda_E,\]
% where $P \approx 4.71803361226...$ is the elliptic logarithm of the generator
% $p = (2:1:1)$ of $E(\mathbf Q)$.
% Note that the y-coordinate of $32p$ is a rational number whose denominator has 24 decimal digits; so $P_1$ corrsponds to a rational point of very large height.
% \subsection{$N=99$}

\subsection{Table}
In table 1 we report the triple Chow heegner points which lie on several strong Weil curves of conductor $<200$. With the exception of the curves labeled {\tt 43a1} and {\tt 65a1} in Cremona's database, we report only in those cases where the differnetials $\alpha_{\omega_{g,i},\eta_{g,i}}$ can be computed explicity and integrated using the method explained in \S \ref{compute}.
The exceptional cases were computed using the method of \S \ref{guess}.

The format of the table is as follows. Let $N$ be the conductor of the curve $E$ in question. We report the points as elements of $E(\mathbf Q)\otimes \mathbf Q$. They are ordered according to the ordering of the isotypic components of the space $\mathbf S_2(\Gamma_0(N))$ of cuspidal modular symbols for $\Gamma_0(N)$, as listed via the command {\tt ModularSymbols(N).cuspidal\_subspace().decomposition()} in {\tt SAGE}, skipping the nonsense point 
that corresponds to taking $g=f$.


% . In this table, the newform $f$ which parametrizes $E$
% usually corresponds to the 0th isotypic component in this list;
% \footnote{Exceptions:
% for $E=$ {\tt 153b1}, $E$ corresponds to the 1st isotypic component, and the listed CH points correspond to components labels 0,2,3,4,5,6,7; similarly,  for $E=$ {\tt 171b1}, $E$ corresponds to the 1st isotypic component,
% and the listed CH points corresond to components labeled 0,2,3,4,5,6,7,8. FIXME: What about {\tt 92b1}, {\tt 106b1}, {\tt 176c1}, {\tt 184b1}?} the Chow-Heegner points $P_{g,f}$ correspond to the remaining components in the order listed below (via $g$ is an eigenform in the corresponding component).
The $-\otimes\mathbf Q$ factor accounts for the denominator $d_g$ of the projector onto the corresponding isotypic component of the anemic Hecke algebra. All the curves $E$ in the list are optimal and have rank 1; they are listed by their Cremona label.  The points $P_{g,f}$ are listed as multiples of a generator $P$ for $E(\mathbf Q)$ computed using {\tt SAGE}.

The last two columns indicate the number of Fourier coefficients used in the computation and a lower bound for number of decimal digits of accuracy to which each of the points we have computed (after adjusting by the periods of the $\alpha_{\omega_{g,i},\eta_{g,i}}$) agree (modulo $\Lambda_E$) with the indicated elements of $E(\mathbf Q)\otimes \mathbf Q$.

\renewcommand{\arraystretch}{1.2}

\begin{table}
\caption{Some elliptic curves with non-torsion triple Chow Heegner points}
\begin{tabular}[h]{|m{1.8cm}|m{1.8cm}||r|r|r||r|r|}

\hline
Curve $E = E_f$ &generator $P\in E(\mathbf Q)$ & $P_{g,f}$ (mod torsion) & $d_g$ & $[E(\mathbf Q):\mathbf Z d_gP_{g,f}]$& \multicolumn{2}{m{2.1cm}|}{accuracy and \# coeff.}\\[.3cm]
\hline\hline
% \end{tabular}

% \begin{tabular}[h]{|c||m{1.8cm}||rl||r|r|}
% \hline

{\tt 37a1}& (0,-1)& $-6P$&2 & 12&54 & 800\\[.1cm]
\hline
{\tt 43a1}& (0,-1)& $4P$& 2 & 8 & 58 & 1200\\[.1cm]
\hline
{\tt 57a1}& (2,1) & $\frac43P$& 12 & 16&31 & 800\\
          &       & $-\frac{16}{3}P$ &3& 16&    &    \\
          &       & $-4P $ & 2&8          &    & \\[ .1cm]
\hline
{\tt 58a1}& (0,-1)& $4P$ & 4  & 16&  32& 800 \\
          &       & 0     & 2 & $\infty$ & &     \\[.1cm]
 \hline
{\tt 61a1}& (1,-1)& $-2P$ & 2 & 4&32 & 800\\[.1cm]
\hline
{\tt 65a1}& (-1,1)& $P$ & 2 & 2& 55 & 3200\\
          &       & $P$ & 2 & 2&   &     \\[.1cm]
\hline
{\tt 82a1}& (0,0) & $0$ & 4 & $\infty$& 33 & 1200\\
          &       & $2P$ & 2 & 4&  &      \\[.1cm]

\hline
{\tt 99a1}& (2,0) & $-\frac23 P$ & 12 & 8& 37 & 1600\\
          &       & $0$ & 12 &$\infty$& & \\
          &       & $\frac23P$ & 6 & 4& & \\
          &       & $\frac23P$ & 12 & 8& & \\
          &       & $-\frac23P$ & 6 & 4& & \\[.1cm]
\hline
{\tt 106b1}& (2,1) & $\frac{12}{5}P$ & 10& 24&38&1800\\
           &       & $-\frac43 P$ & 6 & 8 && \\
           &       & $-\frac{11}{3}P$&48& 176& & \\
           &       & $P$ & 16 & 16& & \\
           &       & $\frac{28}{5}P$ & 10 & 56 && \\[.1cm]
\hline
{\tt 122a1}& (1,-3)& $-\frac{16}{13}P$ &  26 & 32&  28& 1600\\
           &       & $-P$ & 16 & 16 & & \\
           &       & $P$  & 16 & 16 & & \\
           &       & $ -\frac{36}{13}P$ & 26 & 72 & & \\[.1cm]

\hline
{\tt 129a1}& (1,-5)& $-\frac{16}{15}P$& 15& 16&  28& 1600 \\
 &&$-\frac43P$&12 &16&&\\
&&$-\frac{20}{7}P$ & 14&40&&\\
&&$ -\frac85P$& 40&64&&\\
&&$-\frac87P$& 14&16&&\\[.1cm]
\hline
{\tt 153a1}& (0,1)& $P$  & 48 & 48& 24 & 1800  \\
           &      & $2P$ & 24 &   48&  &        \\
           &      & 0    & 24 &    $\infty$&&        \\
           &      & $-P$ & 48 &    48&&       \\
           &      & $-P$ & 16 &    16&&       \\
           &      & $-2P$& 24 &    48&&       \\
           &      & $P$  & 16 &    16&&       \\[.1cm]
\hline
\end{tabular}
 \end{table}
% \begin{table}
% \begin{threeparttable}
% \caption{Some elliptic curves with only torsion Chow-Heegner points}
% \begin{tabular}[h]{|c||m{1.8cm}||r|l|c||r|r||c|}
% \hline
% Curve $E=E_f$ &generator $P\in E(\mathbf Q)$ & $P_{g,f}$& $d_g$ & $E(\mathbf Q)_{\mathrm{tor}}$ & \multicolumn{2}{m{2.1cm}|}{accuracy and \# coeff.}&reason\\[.3cm]
% \hline
% \hline
% {\tt 33a1}& rank 0& 0 & 3 & $\mathbf Z/2\oplus \mathbf Z/2$ &53& 800& OK\\[.1cm]
% \hline
% {\tt 34a1}& rank 0 &$(2,-3)$\tnote{a}& 2& $\mathbf Z/6$& 56 & 800 & OK\\[.1cm]
% \hline
% {\tt 35a1}& rank 0 &$(1,-4)$\tnote{a} & 2 & $\mathbf Z/3$ &54 & 800 & OK\\[.1cm]
% \hline
% {\tt 88a1}& (2,-2)& 3 points, all zero & -\tnote{b} & 0& 19&1600& ()\\[.1cm]
% %  $0$ & 16& 19& 1600\\
% %           &       & $0$ & 8 &   &     \\
% %           &       & $0$ & 2 &   &     \\[.1cm]
% \hline
% {\tt 92b1}& (1,1) & 3 points, all zero & -\tnote{b} & 0 & 30 &1200& ()\\[.1cm]
% % $0$ & 2& 30 & 1200\\
% %           &       & $0$ & 15 &   &     \\
% %           &       & $0$ & 20 &   &     \\[.1cm]
% \hline
% {\tt 112a1}& (0,-2)& 4 points, all zero & -\tnote{b} & $\mathbf Z/2$ &  23& 2400&()\\
%  && $\frac{1}{2}\cdot$(irrational 2-torsion point)\tnote{c}& 16& & &&\\[.1cm]
% % $0$ & 8&23&2400 \\
% %   & & $0$ & 16 & & \\
% %   & & $0$ & 8 & & \\
% %    & & $0$ & 16 & & \\  
% % & & $0$ & 16 & & \\[.1cm]
% \hline

% {\tt 117a1}& $(-\frac14,\frac{15}8)$ & 4 points, all zero & -\tnote{b} & $\mathbf Z/4$ & 23 & 1200&()\\[.1cm]
% % & $0$ & 4 & 23&1200 \\
% %            &                         & $0$ & 8 &   &     \\
% %            &                         & $0$ & 8 &   &     \\
% %            &                         & $0$ & 8 &   &     \\[.1cm]
% \hline
% {\tt 124a1}& (-2,1)& 4 points, all zero & -\tnote{b} & 0 & 30 & 1600 & ()\\[.1cm] 
% % $0$& 6 & 30 &1600 \\
% % &&0&6&&\\
% % &&0&66&&\\
% % &&0&44&&\\[.1cm]
% \hline
% {\tt 136a1}& (2,-2)& 6 points, all zero & -\tnote{b} & $\mathbf Z /2$ & 29 & 3600 & () \\[.1cm]
% % $0$& 8& 29& 3600\\
% % &&0&8&&\\
% % &&0&8&&\\
% % &&0&24&&\\
% % &&0&8&&\\
% % &&0&12&&\\[.1cm]
% \hline
% {\tt 148a1}& (-1,2) & 5 points, all zero & -\tnote{b} & 0 & 25 & 1600&()\\[.1cm]
% % &$0$& 12 & 25 & 1600 \\
% %  && 0 & 20 && \\
% %  && 0 & 12 && \\
% %  && 0 & 6 && \\
% %  && 0 & 30 && \\[.1cm]

% \hline
% {\tt 152a1}& (-1,-2)& 6 points, all zero & -\tnote{b} & 0 & 10 & 1200\\[.1cm]
% % 0 & 8 &10&1200  \\
% % && 0 & 24 &&\\
% % && 0 & 24 &&\\
% % && 0 & 16 &&\\
% % && 0 & 24 &&\\
% % && 0 & 6 &&\\
%  \hline

% {\tt 153b1}& (5,-14)& 7 points, all zero & -\tnote{b} & 0 & 25& 1800\\[.1cm]
% % $0$ & 24 & 25 & 1800 \\
% %            &        & $0$ & 24 &    &       \\
% %            &        & $0$ & 24 &    &       \\
% %            &        & $0$ & 48 &    &       \\
% %            &        & $0$ & 16 &    &       \\
% %            &        & $0$ & 24 &    &       \\
% %            &        & $0$ & 16 &    &       \\[.1cm]
% \hline
% {\tt 171b1}& (2, -5)& 8 points, all zero & -\tnote{b} & 0 & 15 & 1200\\[.1cm]
% % 0 & 24& 15 & 1200\\
% %  && 0 & 96 && \\
% % && 0 & 32 && \\
% % && 0 & 32 && \\
% % && 0 & 96 && \\
% % && 0 & 12 && \\
% % && 0 & 24 && \\
% % && 0 & 12 && \\[.1cm]

% \hline
% {\tt 172a1}&(-3,6)& 5 points, all zero &-\tnote{b}& $\mathbf Z/3$& 18 & 1400 & ()\\[.1cm]
% %$0^{12},\quad 0^{20},\quad 0^{42}$, $\quad 0^{30},\quad 0^{28}$& 16&1200\\
% \hline
% {\tt 176c1}&(1,-2)& 7 points, all zero & -\tnote{b} & 0 & 9 & 1600 \\[.1cm]
% %$0^{2^{13}},\quad 0^{2^{11}\cdot 3}, \quad 0^{2^{11}\cdot 3},\quad 0^{2^{14}}$, $\quad 0^{2^{10}\cdot 3},\quad 0^{2^{14}},\quad 0^{2^{13}}$&9&1600\\
% \hline
% {\tt 184a1}&(0,1)& 8 points, all zero & -\tnote{b} & 0 & 30 & 2400 &\\[.1cm]
% %&$0^{2^{11}},\quad0^{96},\quad0^{2^{11}\cdot 3},\quad 0^{2^8\cdot 3^2},\quad 0^{2^{11}}$, $\quad 0^{2^8},\quad 0^{2^5\cdot 3^2\cdot 5},\quad 0^{2^{12}\cdot 5}$&30&2400\\
% \hline


% \end{tabular}
% \begin{tablenotes}
% \item [a] 3-torsion point
% \item [b] various Hecke indices $d_g$ depending on $g$
% \item [c] Namely, ``$(1,1)$'' if we identify $\mathbf E[2](\mathbf C)\cong \mathbf Z/2\oplus\mathbf Z/2$.
% \end{tablenotes}
% \end{threeparttable}
% \end{table}

\subsection{Discussion}\footnote{FIXME:
Check  when possible the conditions in the theorem saying whe $P_{g,f}$ is non-torsion,

Explain why quadratric twists switch the sign of the points.}

% \section{Appendix: Summary of the algorithm}



% Maybe we should include pseudocode for some of the more interesting parts of the algorithm, e.g. eta(N), hypgens(N), finding matrices of hecke operators using the Poincar\'e pairing, finding primes such that the corresponding hecke operators distinguish isotypic components of $H^1$, finding the denominator
% of the projector in $\mathbf T_{\mathbf Q}$ onto the $g$-isotypic component. 

\section{Extensions of the method}
\subsection{Allowing $f$ to be an oldform}\footnote{FIXME: to be added!}



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\newpage  % tempory for my convenience -- william
\section{Appendix by William Stein: An Alternative approach to
  Numerically Approximating Chow-Heegner Points in Some 
 Special Cases}

In this appendix we give an alternative approach to computing points
$P_{E,F} \in E(\Q)$ associated to a pair of elliptic curves $E, F$,
following a construction of Zhang.  We make no claims about the
complexity or rigor of our method, and explain it mainly because it is
simple to understand compared to the method described
elsewhere in this paper, and provides a double check on the more
sophisticated iterated integral computations described above.











\newpage  % tempory for my convenience -- william
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{thebibliography}{99}

\bibitem[BDP1]{BDP1} M. Bertolini, H. Darmon, K. Prasanna, Generalised Heegner cycles and p-adic
Rankin L-series, submitted for publication.

\bibitem[BDP2]{BDP2} M. Bertolini, H. Darmon, K. Prasanna, Chow-Heegner points on CM elliptic curves
and values of $p$-adic L-series, in progress.

\bibitem[BL]{BL} Ch. Birkenhake, H. Lange, \textit{Complex abelian varieties}. Grundlehren der mathematischen Wissenschaften \textbf{302}. Springer-Verlag (1992).

\bibitem[BO]{ono} K. Bringmann, K. Ono, Coefficients of harmonic weak Maass forms, \textit{Proceedings of the 2008 University of Florida conference on partitions, $q$-series, and modular forms}.

\bibitem[BrPh]{BrPh} N. Brisebarre, G. Philibert, Lower and upper bounds for the Fourier coefficients of powers of the modular invariant $j$, \textit{J. Ramanujan Math. Soc.} \textbf{20}, No. 4 (2005), 255-282. 

\bibitem[Ch]{Chen} K.-T. Chen, Iterated path integrals, {\em Bull. Amer. Math. Soc.} {\bf 83} (1977), 831--879.

\bibitem[Cr]{cremona}J. Cremona, Chapter 2: Modular Symbol Algorithms, in \textit{Algorithms for Modular Elliptic curves}. \url{http://www.warwick.ac.uk/~masgaj/book/fulltext/index.html}.

\bibitem[CWZ]{genera} J. A. Csirik, J. L. Wetherell, M. E. Zieve, On the genera of $X_0(N)$, \textit{J. Number Theory}

\bibitem[D]{Dar}
H. Darmon, Cycles on modular varieties and 
rational points on elliptic curves, in {\em Explicit Methods in Number Theory,
Oberwolfach reports} {\bf 6:3} (2009), 1843--1920. 

\bibitem[DRS]{DRS} H. Darmon, V. Rotger, I. Sols, Iterated integrals of modular forms\\ and rational points on elliptic curves, preprint March 2011.

\bibitem[DRS2]{DRS2} H. Darmon, V. Rotger, I. Sols, Diagonal cycles on triple products of Kuga-Sato varieties, in progress.

%\bibitem[FK]{FK} Farkas, Kra, Riemann surfaces.

\bibitem[G]{Gonzalez} J. Gonz\'alez, Equations of hyperelliptic modular curves, {\em Ann. inst. Fourier} {\bf } (), --.

\bibitem[H1]{Hain} R. Hain, Lectures on the Hodge-de Rham theory  of the fundamental group of $\mathbf P^1 - \{0,1,\infty\}$. \url{http://math.arizona.edu/~swc/notes/files/05HainNotes.pdf}.

% \bibitem[H2]{Hain2} R. Hain, The Geometry of the Mixed Hodge Structure on the Fundamental Group, {\em Proc. Symp. Pure Math.} \textbf{48} (1987), 247--282.

\bibitem[O]{Ogg} A. Ogg, On the Weierstrass points of $X_0(N)$, {\em Illinois J. Math.} \textbf{22} (1978), 31--35. 

\bibitem[R]{rad} Rademacher, The Fourier coefficients of the modular ivnariant $J(\tau)$.

\bibitem[S1]{S1} W. Stein, \textit{Generating the Hecke algebra}. \url{http://www.wstein.org/papers/generating_hecke/generating_hecke.ps}.

\bibitem[S2]{S2} W. Stein, \textit{Modular forms: a computational approach}. \url{http://wstein.org/books/modform/modform/}.

\bibitem[St]{Stevens} G. Stevens, \textit{Arithmetic on modular curves}.

% \bibitem[Ya]{Ya} Y. Yang, Defining equations of modular curves, {\em Advances Math.} {\bf 204} (2006), 481-508.

%\bibitem[Zh]{Zh} S. Zhang. {\em Arithmetic of Shimura curves}. Science China Mathematics {\bf 53} (2010), 573--592.



\end{thebibliography}

\end{document}

