acb_theta_ql_jet_fd, acb_theta_ql_jet

The finite difference method evaluates \theta_{a,b}(z + h_n, \tau) at h_n = (\varepsilon \zeta^{n_0},\ldots, \varepsilon \zeta^{n_{g-1}}) for a primitive m-th root of unity \zeta, m = \mathit{ord} + 1, and all n \in \{0,\ldots, m - 1\}^g; a discrete Fourier transform gives m^g \sum_{k \equiv p \bmod m} a_k \varepsilon^{|k|} for each p, where a_k are the Taylor coefficients. By the Cauchy formula with the bound of acb_theta_ql_local_bound, the error is at most 2cg\varepsilon^{|p|+m}/\rho^{m} before division by \varepsilon^{|p|}. The function chooses \varepsilon, strips (z, \tau) of their radii, raises the precision to compensate for the division by \varepsilon^{\mathit{ord}} (\mathit{ord}+1)^g, and restores certified radii with acb_theta_ql_jet_error and derivatives to order \mathit{ord} + 2 computed at low precision; see Derivatives.