From c266d14a9d100d76369faefc22aa390407e8fd7d Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Thu, 18 Apr 2024 13:57:44 -0400 Subject: [PATCH 01/11] initial commit to acb_theta via rst_to_pxd.py --- src/flint/flintlib/acb_theta.pxd | 121 +++++++++++++++++++++++++++++++ 1 file changed, 121 insertions(+) create mode 100644 src/flint/flintlib/acb_theta.pxd diff --git a/src/flint/flintlib/acb_theta.pxd b/src/flint/flintlib/acb_theta.pxd new file mode 100644 index 00000000..9f47de39 --- /dev/null +++ b/src/flint/flintlib/acb_theta.pxd @@ -0,0 +1,121 @@ +from flint.flintlib.acb cimport acb_t, acb_srcptr, acb_ptr +from flint.flintlib.acb_poly cimport acb_poly_struct, acb_poly_t +from flint.flintlib.arb cimport arb_t, arb_ptr, arb_srcptr +from flint.flintlib.flint cimport ulong, slong, flint_rand_t +from flint.flintlib.fmpz_mat cimport fmpz_mat_t +from flint.flintlib.acb_mat cimport acb_mat_t +from flint.flintlib.arb_mat cimport arb_mat_t +from flint.flintlib.arf cimport arf_t + + +# unimported types {'acb_theta_eld_t', 'acb_theta_naive_worker_t', 'acb_theta_ql_worker_t'} + +cdef extern from "flint/acb_theta.h": + void acb_theta_all(acb_ptr th, acb_srcptr z, const acb_mat_t tau, int sqr, slong prec) + void acb_theta_naive_fixed_ab(acb_ptr th, ulong ab, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) + void acb_theta_naive_all(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) + void acb_theta_jet_all(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec) + void acb_theta_jet_naive_fixed_ab(acb_ptr dth, ulong ab, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec) + void acb_theta_jet_naive_all(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec) + slong sp2gz_dim(const fmpz_mat_t mat) + void sp2gz_set_blocks(fmpz_mat_t mat, const fmpz_mat_t alpha, const fmpz_mat_t beta, const fmpz_mat_t gamma, const fmpz_mat_t delta) + void sp2gz_j(fmpz_mat_t mat) + void sp2gz_block_diag(fmpz_mat_t mat, const fmpz_mat_t U) + void sp2gz_trig(fmpz_mat_t mat, const fmpz_mat_t S) + void sp2gz_embed(fmpz_mat_t res, const fmpz_mat_t mat) + void sp2gz_restrict(fmpz_mat_t res, const fmpz_mat_t mat) + slong sp2gz_nb_fundamental(slong g) + void sp2gz_fundamental(fmpz_mat_t mat, slong j) + int sp2gz_is_correct(const fmpz_mat_t mat) + int sp2gz_is_j(const fmpz_mat_t mat) + int sp2gz_is_block_diag(const fmpz_mat_t mat) + int sp2gz_is_trig(const fmpz_mat_t mat) + int sp2gz_is_embedded(fmpz_mat_t res, const fmpz_mat_t mat) + void sp2gz_inv(fmpz_mat_t inv, const fmpz_mat_t mat) + fmpz_mat_struct * sp2gz_decompose(slong * nb, const fmpz_mat_t mat) + void sp2gz_randtest(fmpz_mat_t mat, flint_rand_t state, slong bits) + void acb_siegel_cocycle(acb_mat_t c, const fmpz_mat_t mat, const acb_mat_t tau, slong prec) + void acb_siegel_transform_cocycle_inv(acb_mat_t w, acb_mat_t c, acb_mat_t cinv, const fmpz_mat_t mat, const acb_mat_t tau, slong prec) + void acb_siegel_transform(acb_mat_t w, const fmpz_mat_t mat, const acb_mat_t tau, slong prec) + void acb_siegel_transform_z(acb_ptr r, acb_mat_t w, const fmpz_mat_t mat, acb_srcptr z, const acb_mat_t tau, slong prec) + void acb_siegel_cho(arb_mat_t C, const acb_mat_t tau, slong prec) + void acb_siegel_yinv(arb_mat_t Yinv, const acb_mat_t tau, slong prec) + void acb_siegel_reduce(fmpz_mat_t mat, const acb_mat_t tau, slong prec) + int acb_siegel_is_reduced(const acb_mat_t tau, slong tol_exp, slong prec) + void acb_siegel_randtest(acb_mat_t tau, flint_rand_t state, slong prec, slong mag_bits) + void acb_siegel_randtest_reduced(acb_mat_t tau, flint_rand_t state, slong prec, slong mag_bits) + void acb_siegel_randtest_vec(acb_ptr z, flint_rand_t state, slong g, slong prec) + void acb_theta_char_get_slong(slong * n, ulong a, slong g) + ulong acb_theta_char_get_a(const slong * n, slong g) + void acb_theta_char_get_arb(arb_ptr v, ulong a, slong g) + void acb_theta_char_get_acb(acb_ptr v, ulong a, slong g) + slong acb_theta_char_dot(ulong a, ulong b, slong g) + slong acb_theta_char_dot_slong(ulong a, const slong * n, slong g) + void acb_theta_char_dot_acb(acb_t x, ulong a, acb_srcptr z, slong g, slong prec) + int acb_theta_char_is_even(ulong ab, slong g) + int acb_theta_char_is_goepel(ulong ch1, ulong ch2, ulong ch3, ulong ch4, slong g) + int acb_theta_char_is_syzygous(ulong ch1, ulong ch2, ulong ch3, slong g) + void acb_theta_eld_init(acb_theta_eld_t E, slong d, slong g) + void acb_theta_eld_clear(acb_theta_eld_t E) + int acb_theta_eld_set(acb_theta_eld_t E, const arb_mat_t C, const arf_t R2, arb_srcptr v) + void acb_theta_eld_points(slong * pts, const acb_theta_eld_t E) + void acb_theta_eld_border(slong * pts, const acb_theta_eld_t E) + int acb_theta_eld_contains(const acb_theta_eld_t E, slong * pt) + void acb_theta_eld_print(const acb_theta_eld_t E) + void acb_theta_naive_radius(arf_t R2, arf_t eps, const arb_mat_t C, slong ord, slong prec) + void acb_theta_naive_reduce(arb_ptr v, acb_ptr new_zs, arb_ptr as, acb_ptr cs, arb_ptr us, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) + void acb_theta_naive_term(acb_t res, acb_srcptr z, const acb_mat_t tau, slong * tup, slong * n, slong prec) + void acb_theta_naive_worker(acb_ptr th, slong len, acb_srcptr zs, slong nb, const acb_mat_t tau, const acb_theta_eld_t E, slong ord, slong prec, acb_theta_naive_worker_t worker) + void acb_theta_naive_00(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) + void acb_theta_naive_0b(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) + void acb_theta_naive_fixed_a(acb_ptr th, ulong a, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) + slong acb_theta_jet_nb(slong ord, slong g) + slong acb_theta_jet_total_order(const slong * tup, slong g) + void acb_theta_jet_tuples(slong * tups, slong ord, slong g) + slong acb_theta_jet_index(const slong * tup, slong g) + void acb_theta_jet_mul(acb_ptr res, acb_srcptr v1, acb_srcptr v2, slong ord, slong g, slong prec) + void acb_theta_jet_compose(acb_ptr res, acb_srcptr v, const acb_mat_t N, slong ord, slong prec) + void acb_theta_jet_exp_pi_i(acb_ptr res, arb_srcptr a, slong ord, slong g, slong prec) + void acb_theta_jet_naive_radius(arf_t R2, arf_t eps, arb_srcptr v, const arb_mat_t C, slong ord, slong prec) + void acb_theta_jet_naive_00(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec) + void acb_theta_jet_error_bounds(arb_ptr err, acb_srcptr z, const acb_mat_t tau, acb_srcptr dth, slong ord, slong prec) + void acb_theta_dist_pt(arb_t d, arb_srcptr v, const arb_mat_t C, slong * n, slong prec) + void acb_theta_dist_lat(arb_t d, arb_srcptr v, const arb_mat_t C, slong prec) + void acb_theta_dist_a0(arb_ptr d, acb_srcptr z, const acb_mat_t tau, slong prec) + slong acb_theta_dist_addprec(const arb_t d) + void acb_theta_agm_hadamard(acb_ptr res, acb_srcptr a, slong g, slong prec) + void acb_theta_agm_sqrt(acb_ptr res, acb_srcptr a, acb_srcptr rts, slong nb, slong prec) + void acb_theta_agm_mul(acb_ptr res, acb_srcptr a1, acb_srcptr a2, slong g, slong prec) + void acb_theta_agm_mul_tight(acb_ptr res, acb_srcptr a0, acb_srcptr a, arb_srcptr d0, arb_srcptr d, slong g, slong prec) + int acb_theta_ql_a0_naive(acb_ptr th, acb_srcptr t, acb_srcptr z, arb_srcptr d0, arb_srcptr d, const acb_mat_t tau, slong guard, slong prec) + int acb_theta_ql_a0_split(acb_ptr th, acb_srcptr t, acb_srcptr z, arb_srcptr d, const acb_mat_t tau, slong s, slong guard, slong prec, acb_theta_ql_worker_t worker) + int acb_theta_ql_a0_steps(acb_ptr th, acb_srcptr t, acb_srcptr z, arb_srcptr d0, arb_srcptr d, const acb_mat_t tau, slong nb_steps, slong s, slong guard, slong prec, acb_theta_ql_worker_t worker) + slong acb_theta_ql_a0_nb_steps(const arb_mat_t C, slong s, slong prec) + int acb_theta_ql_a0(acb_ptr th, acb_srcptr t, acb_srcptr z, arb_srcptr d0, arb_srcptr d, const acb_mat_t tau, slong guard, slong prec) + slong acb_theta_ql_reduce(acb_ptr new_z, acb_t c, arb_t u, slong * n1, acb_srcptr z, const acb_mat_t tau, slong prec) + void acb_theta_ql_all(acb_ptr th, acb_srcptr z, const acb_mat_t tau, int sqr, slong prec) + void acb_theta_jet_ql_bounds(arb_t c, arb_t rho, acb_srcptr z, const acb_mat_t tau, slong ord) + void acb_theta_jet_ql_radius(arf_t eps, arf_t err, const arb_t c, const arb_t rho, slong ord, slong g, slong prec) + void acb_theta_jet_ql_finite_diff(acb_ptr dth, const arf_t eps, const arf_t err, acb_srcptr val, slong ord, slong g, slong prec) + void acb_theta_jet_ql_all(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec) + ulong acb_theta_transform_char(slong * e, const fmpz_mat_t mat, ulong ab) + void acb_theta_transform_sqrtdet(acb_t res, const acb_mat_t tau, slong prec) + slong acb_theta_transform_kappa(acb_t sqrtdet, const fmpz_mat_t mat, const acb_mat_t tau, slong prec) + slong acb_theta_transform_kappa2(const fmpz_mat_t mat) + void acb_theta_transform_proj(acb_ptr res, const fmpz_mat_t mat, acb_srcptr th, int sqr, slong prec) + void acb_theta_g2_jet_naive_1(acb_ptr dth, const acb_mat_t tau, slong prec) + void acb_theta_g2_detk_symj(acb_poly_t res, const acb_mat_t m, const acb_poly_t f, slong k, slong j, slong prec) + void acb_theta_g2_transvectant(acb_poly_t res, const acb_poly_t g, const acb_poly_t h, slong m, slong n, slong k, slong prec) + void acb_theta_g2_transvectant_lead(acb_t res, const acb_poly_t g, const acb_poly_t h, slong m, slong n, slong k, slong prec) + slong acb_theta_g2_character(const fmpz_mat_t mat) + void acb_theta_g2_psi4(acb_t res, acb_srcptr th2, slong prec) + void acb_theta_g2_psi6(acb_t res, acb_srcptr th2, slong prec) + void acb_theta_g2_chi10(acb_t res, acb_srcptr th2, slong prec) + void acb_theta_g2_chi12(acb_t res, acb_srcptr th2, slong prec) + void acb_theta_g2_chi5(acb_t res, acb_srcptr th, slong prec) + void acb_theta_g2_chi35(acb_t res, acb_srcptr th, slong prec) + void acb_theta_g2_chi3_6(acb_poly_t res, acb_srcptr dth, slong prec) + void acb_theta_g2_sextic(acb_poly_t res, const acb_mat_t tau, slong prec) + void acb_theta_g2_sextic_chi5(acb_poly_t res, acb_t chi5, const acb_mat_t tau, slong prec) + void acb_theta_g2_covariants(acb_poly_struct * res, const acb_poly_t f, slong prec) + void acb_theta_g2_covariants_lead(acb_ptr res, const acb_poly_t f, slong prec) From bf2bca747796c19abc0c85b9c4150cffd8690687 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Thu, 13 Jun 2024 17:05:12 -0400 Subject: [PATCH 02/11] Wrapping acb_theta_all --- src/flint/flintlib/acb_theta.pxd | 21 +++++++++++++++++--- src/flint/types/acb_mat.pyx | 33 ++++++++++++++++++++++++++++++++ 2 files changed, 51 insertions(+), 3 deletions(-) diff --git a/src/flint/flintlib/acb_theta.pxd b/src/flint/flintlib/acb_theta.pxd index 9f47de39..2ef2d405 100644 --- a/src/flint/flintlib/acb_theta.pxd +++ b/src/flint/flintlib/acb_theta.pxd @@ -2,15 +2,30 @@ from flint.flintlib.acb cimport acb_t, acb_srcptr, acb_ptr from flint.flintlib.acb_poly cimport acb_poly_struct, acb_poly_t from flint.flintlib.arb cimport arb_t, arb_ptr, arb_srcptr from flint.flintlib.flint cimport ulong, slong, flint_rand_t -from flint.flintlib.fmpz_mat cimport fmpz_mat_t +from flint.flintlib.fmpz_mat cimport fmpz_mat_t, fmpz_mat_struct from flint.flintlib.acb_mat cimport acb_mat_t from flint.flintlib.arb_mat cimport arb_mat_t from flint.flintlib.arf cimport arf_t -# unimported types {'acb_theta_eld_t', 'acb_theta_naive_worker_t', 'acb_theta_ql_worker_t'} - cdef extern from "flint/acb_theta.h": + ctypedef struct acb_theta_eld_struct: + slong dim + slong ambient_dim + slong * last_coords + slong min + slong mid + slong max + slong nr + slong nl + acb_theta_eld_struct * rchildren + acb_theta_eld_struct * lchildren + slong nb_pts, nb_border + slong * box + ctypedef acb_theta_eld_struct acb_theta_eld_t[1] + ctypedef void (*acb_theta_naive_worker_t)(acb_ptr, acb_srcptr, acb_srcptr, const slong *, slong, const acb_t, const slong *, slong, slong, slong, slong) + ctypedef int (*acb_theta_ql_worker_t)(acb_ptr, acb_srcptr, acb_srcptr, arb_srcptr, arb_srcptr, const acb_mat_t, slong, slong) + void acb_theta_all(acb_ptr th, acb_srcptr z, const acb_mat_t tau, int sqr, slong prec) void acb_theta_naive_fixed_ab(acb_ptr th, ulong ab, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) void acb_theta_naive_all(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec) diff --git a/src/flint/types/acb_mat.pyx b/src/flint/types/acb_mat.pyx index 533bdc72..8119f153 100644 --- a/src/flint/types/acb_mat.pyx +++ b/src/flint/types/acb_mat.pyx @@ -19,6 +19,7 @@ from flint.flintlib.arb_mat cimport * from flint.flintlib.arf cimport * from flint.flintlib.acb cimport * from flint.flintlib.acb_mat cimport * +from flint.flintlib.acb_theta cimport * cimport cython @@ -814,3 +815,35 @@ cdef class acb_mat(flint_mat): if left: return Elist, L return Elist, R + + def theta(tau, z, square=False): + r""" + FIXME + """ + g = tau.nrows() + assert g > 1 + assert tau.ncols() == g + assert z.nrows() == g + assert z.ncols() == 1 + + # convert input + cdef acb_ptr zvec + zvec = _acb_vec_init(g) + cdef long i + for i in range(g): + acb_set(zvec + i, acb_mat_entry((z).val, i, 0)) + + # initialize the output + cdef slong nb = 1 << (2 * g) + cdef acb_ptr theta = _acb_vec_init(nb) + + acb_theta_all(theta, zvec, tau.val, square, getprec()) + _acb_vec_clear(zvec, g) + # copy the output + res = tuple() + for i in range(nb): + r = acb.__new__(acb) + acb_set((r).val, theta + i) + res += (r,) + _acb_vec_clear(theta, nb) + return res From 0450b6fb75370ebf323cd159fb494db024a1bbd8 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Fri, 14 Jun 2024 14:34:59 -0400 Subject: [PATCH 03/11] adding documentation --- src/flint/types/acb_mat.pyx | 19 +++++++++++++++++-- 1 file changed, 17 insertions(+), 2 deletions(-) diff --git a/src/flint/types/acb_mat.pyx b/src/flint/types/acb_mat.pyx index 8119f153..04d7eb42 100644 --- a/src/flint/types/acb_mat.pyx +++ b/src/flint/types/acb_mat.pyx @@ -818,10 +818,25 @@ cdef class acb_mat(flint_mat): def theta(tau, z, square=False): r""" - FIXME + Computes the vector valued Riemann theta function `(\theta_{a,b}{z, tau) : a, b \in \{0,1\}^{g}\)` or its squares. + This is a wrapper for the function `acb_theta_all` and it follows the same conventions for the ordering of the theta characteristics. + + >>> from flint import acb, acb_mat, showgood + >>> z = acb(1+1j); tau = acb(1.25+3j) + >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])) + >>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])]) + [+/- 3.82e-14] + >>> for i in range(4):showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]]))[i], dps=25) + ... + 0.9694430387796704100046143 - 0.03055696120816803328582847j + 1.030556961196006476576271 + 0.03055696120816803328582847j + -1.220790267576967690128359 - 1.827055516791154669091679j + -1.820235910124989594900076 + 1.216251950154477951760042j + >>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])) + ([1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j, [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [0.451696791791346 +/- 5.46e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j) + """ g = tau.nrows() - assert g > 1 assert tau.ncols() == g assert z.nrows() == g assert z.ncols() == 1 From 8379c5ebdb49a2f289b658c825cb4c0bc80fca37 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Fri, 14 Jun 2024 16:56:19 -0400 Subject: [PATCH 04/11] making a new module for acb_theta --- meson.build | 3 ++ src/flint/types/acb_mat.pyx | 58 ++++++++--------------------------- src/flint/types/acb_theta.pyx | 56 +++++++++++++++++++++++++++++++++ src/flint/types/meson.build | 4 +++ 4 files changed, 76 insertions(+), 45 deletions(-) create mode 100644 src/flint/types/acb_theta.pyx diff --git a/meson.build b/meson.build index 7e941dbb..d5f4338b 100644 --- a/meson.build +++ b/meson.build @@ -22,6 +22,9 @@ endif # flint.pc was missing -lflint until Flint 3.1.0 if flint_dep.version().version_compare('<3.1') flint_dep = cc.find_library('flint') + have_acb_theta = false +else + have_acb_theta = true endif pyflint_deps = [dep_py, gmp_dep, mpfr_dep, flint_dep] diff --git a/src/flint/types/acb_mat.pyx b/src/flint/types/acb_mat.pyx index 04d7eb42..6066e47d 100644 --- a/src/flint/types/acb_mat.pyx +++ b/src/flint/types/acb_mat.pyx @@ -19,7 +19,6 @@ from flint.flintlib.arb_mat cimport * from flint.flintlib.arf cimport * from flint.flintlib.acb cimport * from flint.flintlib.acb_mat cimport * -from flint.flintlib.acb_theta cimport * cimport cython @@ -818,47 +817,16 @@ cdef class acb_mat(flint_mat): def theta(tau, z, square=False): r""" - Computes the vector valued Riemann theta function `(\theta_{a,b}{z, tau) : a, b \in \{0,1\}^{g}\)` or its squares. - This is a wrapper for the function `acb_theta_all` and it follows the same conventions for the ordering of the theta characteristics. - - >>> from flint import acb, acb_mat, showgood - >>> z = acb(1+1j); tau = acb(1.25+3j) - >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])) - >>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])]) - [+/- 3.82e-14] - >>> for i in range(4):showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]]))[i], dps=25) - ... - 0.9694430387796704100046143 - 0.03055696120816803328582847j - 1.030556961196006476576271 + 0.03055696120816803328582847j - -1.220790267576967690128359 - 1.827055516791154669091679j - -1.820235910124989594900076 + 1.216251950154477951760042j - >>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])) - ([1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j, [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [0.451696791791346 +/- 5.46e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j) - - """ - g = tau.nrows() - assert tau.ncols() == g - assert z.nrows() == g - assert z.ncols() == 1 - - # convert input - cdef acb_ptr zvec - zvec = _acb_vec_init(g) - cdef long i - for i in range(g): - acb_set(zvec + i, acb_mat_entry((z).val, i, 0)) - - # initialize the output - cdef slong nb = 1 << (2 * g) - cdef acb_ptr theta = _acb_vec_init(nb) - - acb_theta_all(theta, zvec, tau.val, square, getprec()) - _acb_vec_clear(zvec, g) - # copy the output - res = tuple() - for i in range(nb): - r = acb.__new__(acb) - acb_set((r).val, theta + i) - res += (r,) - _acb_vec_clear(theta, nb) - return res + Computes the vector valued Riemann theta function + `(\theta_{a,b}{z, tau) : a, b \in \{0,1\}^{g}\)` or its squares. + This is a wrapper for the C-function `acb_theta_all` and it follows the + same conventions for the ordering of the theta characteristics. + + This is a wrapper for :meth:`.acb_theta.acb_mat_theta`; see the + documentation for that method for details for examples. + """ + try: + from .acb_theta import acb_mat_theta + except ImportError: + raise NotImplementedError("acb_mat.theta needs Flint >= 3.1.0") + return acb_mat_theta(z, tau, square=square) diff --git a/src/flint/types/acb_theta.pyx b/src/flint/types/acb_theta.pyx new file mode 100644 index 00000000..d2d8f116 --- /dev/null +++ b/src/flint/types/acb_theta.pyx @@ -0,0 +1,56 @@ +from flint.flint_base.flint_context cimport getprec +from flint.types.acb cimport acb +from flint.types.acb_mat cimport acb_mat +from flint.flintlib.acb cimport * +from flint.flintlib.acb_mat cimport * +from flint.flintlib.acb_theta cimport * + +def acb_mat_theta(acb_mat z, acb_mat tau, ulong square=False): + r""" + Computes the vector valued Riemann theta function `(\theta_{a,b}{z, tau) : a, b \in \{0,1\}^{g}\)` or its squares. + This is a wrapper for the function `acb_theta_all` and it follows the same conventions for the ordering of the theta characteristics. + + This should be used via method `acb_mat.theta` with the order of `z` and `tau` swapped, + + >>> from flint import acb, acb_mat, showgood + >>> z = acb(1+1j); tau = acb(1.25+3j) + >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])) + >>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])]) + [+/- 3.82e-14] + >>> for i in range(4):showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]]))[i], dps=25) + ... + 0.9694430387796704100046143 - 0.03055696120816803328582847j + 1.030556961196006476576271 + 0.03055696120816803328582847j + -1.220790267576967690128359 - 1.827055516791154669091679j + -1.820235910124989594900076 + 1.216251950154477951760042j + >>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])) + ([1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j, [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [0.451696791791346 +/- 5.46e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j) + + """ + g = tau.nrows() + assert tau.ncols() == g + assert z.nrows() == g + assert z.ncols() == 1 + + # convert input + cdef acb_ptr zvec + zvec = _acb_vec_init(g) + cdef long i + for i in range(g): + acb_set(zvec + i, acb_mat_entry(z.val, i, 0)) + + # initialize the output + cdef slong nb = 1 << (2 * g) + cdef acb_ptr theta = _acb_vec_init(nb) + + acb_theta_all(theta, zvec, tau.val, square, getprec()) + _acb_vec_clear(zvec, g) + # copy the output + res = tuple() + cdef acb r + for i in range(nb): + r = acb.__new__(acb) + acb_set(r.val, theta + i) + res += (r,) + _acb_vec_clear(theta, nb) + return res diff --git a/src/flint/types/meson.build b/src/flint/types/meson.build index bf8caec5..7a90bb41 100644 --- a/src/flint/types/meson.build +++ b/src/flint/types/meson.build @@ -41,6 +41,10 @@ exts = [ 'fmpz_mpoly', ] +if have_acb_theta + exts += ['acb_theta'] +endif + py.install_sources( pyfiles, pure: false, From f59c5fe67a3b1eb6063610601e7daf7c5c963b78 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Fri, 14 Jun 2024 17:11:46 -0400 Subject: [PATCH 05/11] add acb_theta to doctests --- src/flint/test/__main__.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/flint/test/__main__.py b/src/flint/test/__main__.py index d16b0aa5..1bd13dec 100644 --- a/src/flint/test/__main__.py +++ b/src/flint/test/__main__.py @@ -80,6 +80,8 @@ def run_doctests(verbose=None): flint.types.acb_series, flint.types.dirichlet, flint.functions.showgood] + if hasattr(flint.types, 'acb_theta'): + modules.append(flint.types.acb_theta) results = [doctest.testmod(x) for x in modules] # ffmpz, tfmpz = doctest.testmod(flint._fmpz, verbose=verbose) # failed, total = doctest.testmod(flint.pyflint, verbose=verbose) From 8c74c77fe197e8242f6a4db35b83121c61e3ead1 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Fri, 14 Jun 2024 17:55:37 -0400 Subject: [PATCH 06/11] return a matrix instead of a tuple --- src/flint/types/acb_theta.pyx | 29 ++++++++++++++++++++++------- 1 file changed, 22 insertions(+), 7 deletions(-) diff --git a/src/flint/types/acb_theta.pyx b/src/flint/types/acb_theta.pyx index d2d8f116..88a2e330 100644 --- a/src/flint/types/acb_theta.pyx +++ b/src/flint/types/acb_theta.pyx @@ -14,17 +14,32 @@ def acb_mat_theta(acb_mat z, acb_mat tau, ulong square=False): >>> from flint import acb, acb_mat, showgood >>> z = acb(1+1j); tau = acb(1.25+3j) - >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])) + >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])).entries() >>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])]) [+/- 3.82e-14] - >>> for i in range(4):showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]]))[i], dps=25) + >>> for i in range(4):showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]])).entries()[i], dps=25) ... 0.9694430387796704100046143 - 0.03055696120816803328582847j 1.030556961196006476576271 + 0.03055696120816803328582847j -1.220790267576967690128359 - 1.827055516791154669091679j -1.820235910124989594900076 + 1.216251950154477951760042j - >>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])) - ([1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j, [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [0.451696791791346 +/- 5.46e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j, [0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j, [+/- 2.43e-16] + [+/- 2.43e-16]j) + >>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])).transpose() + [ [1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j] + [ [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j] + [[0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j] + [[0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j] + [[0.451696791791346 +/- 5.46e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + [[0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + [[0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j] + [[0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + [[0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + [ [+/- 2.43e-16] + [+/- 2.43e-16]j] """ g = tau.nrows() @@ -46,11 +61,11 @@ def acb_mat_theta(acb_mat z, acb_mat tau, ulong square=False): acb_theta_all(theta, zvec, tau.val, square, getprec()) _acb_vec_clear(zvec, g) # copy the output - res = tuple() + res = [] cdef acb r for i in range(nb): r = acb.__new__(acb) acb_set(r.val, theta + i) - res += (r,) + res.append(r) _acb_vec_clear(theta, nb) - return res + return acb_mat([res]) From fbabf32fdf2e6d7d08545c8a815e957cbb547e3c Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Sat, 15 Jun 2024 11:40:45 -0400 Subject: [PATCH 07/11] fix doctests --- src/flint/test/__main__.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/flint/test/__main__.py b/src/flint/test/__main__.py index 1bd13dec..3a874dee 100644 --- a/src/flint/test/__main__.py +++ b/src/flint/test/__main__.py @@ -80,11 +80,14 @@ def run_doctests(verbose=None): flint.types.acb_series, flint.types.dirichlet, flint.functions.showgood] - if hasattr(flint.types, 'acb_theta'): - modules.append(flint.types.acb_theta) + try: + from flint.types import acb_theta + modules.append(acb_theta) + except ImportError: + pass results = [doctest.testmod(x) for x in modules] # ffmpz, tfmpz = doctest.testmod(flint._fmpz, verbose=verbose) - # failed, total = doctest.testmod(flint.pyflint, verbose=verbose) +# failed, total = doctest.testmod(flint.pyflint, verbose=verbose) return tuple(sum(res) for res in zip(*results)) From 0cae1607436cf75dbe8d677fa4ae35efbe19ec8a Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Sat, 15 Jun 2024 17:00:43 -0400 Subject: [PATCH 08/11] documentation --- doc/source/acb_theta.rst | 5 +++++ doc/source/index.rst | 1 + src/flint/types/acb_mat.pyx | 17 ++++++++++------- src/flint/types/acb_theta.pyx | 12 ++++++++---- 4 files changed, 24 insertions(+), 11 deletions(-) create mode 100644 doc/source/acb_theta.rst diff --git a/doc/source/acb_theta.rst b/doc/source/acb_theta.rst new file mode 100644 index 00000000..f9ccfdd2 --- /dev/null +++ b/doc/source/acb_theta.rst @@ -0,0 +1,5 @@ +**acb_theta** -- Riemann theta functions +=============================================================================== + +.. autofunction :: flint.types.acb_theta.acb_theta + diff --git a/doc/source/index.rst b/doc/source/index.rst index 0a1acb33..8e95e722 100644 --- a/doc/source/index.rst +++ b/doc/source/index.rst @@ -61,6 +61,7 @@ Matrix types fmpz_mod_mat.rst arb_mat.rst acb_mat.rst + acb_theta.rst Polynomial types ................ diff --git a/src/flint/types/acb_mat.pyx b/src/flint/types/acb_mat.pyx index 6066e47d..d33a65f1 100644 --- a/src/flint/types/acb_mat.pyx +++ b/src/flint/types/acb_mat.pyx @@ -818,15 +818,18 @@ cdef class acb_mat(flint_mat): def theta(tau, z, square=False): r""" Computes the vector valued Riemann theta function - `(\theta_{a,b}{z, tau) : a, b \in \{0,1\}^{g}\)` or its squares. - This is a wrapper for the C-function `acb_theta_all` and it follows the - same conventions for the ordering of the theta characteristics. + `(\theta_{a,b}(z, \tau) : a, b \in \{0,1\}^{g})` or its squares, + where `\tau` is given by ``self``. + + This is a wrapper for :func:`.acb_theta.acb_theta`; see the + documentation for that method for details and examples. + This follows the same conventions of the C-function + `acb_theta_all `_ + for the ordering of the theta characteristics. - This is a wrapper for :meth:`.acb_theta.acb_mat_theta`; see the - documentation for that method for details for examples. """ try: - from .acb_theta import acb_mat_theta + from .acb_theta import acb_theta except ImportError: raise NotImplementedError("acb_mat.theta needs Flint >= 3.1.0") - return acb_mat_theta(z, tau, square=square) + return acb_theta(z, tau, square=square) diff --git a/src/flint/types/acb_theta.pyx b/src/flint/types/acb_theta.pyx index 88a2e330..ea8e4907 100644 --- a/src/flint/types/acb_theta.pyx +++ b/src/flint/types/acb_theta.pyx @@ -5,12 +5,16 @@ from flint.flintlib.acb cimport * from flint.flintlib.acb_mat cimport * from flint.flintlib.acb_theta cimport * -def acb_mat_theta(acb_mat z, acb_mat tau, ulong square=False): +def acb_theta(acb_mat z, acb_mat tau, ulong square=False): r""" - Computes the vector valued Riemann theta function `(\theta_{a,b}{z, tau) : a, b \in \{0,1\}^{g}\)` or its squares. - This is a wrapper for the function `acb_theta_all` and it follows the same conventions for the ordering of the theta characteristics. + Computes the vector valued Riemann theta function + `(\theta_{a,b}(z, \tau) : a, b \in \{0,1\}^{g})` or its squares. - This should be used via method `acb_mat.theta` with the order of `z` and `tau` swapped, + This is a wrapper for the C-function + `acb_theta_all `_ + and it follows the same conventions for the ordering of the theta characteristics. + + This should be used via the method :meth:`.acb_mat.theta`, explicitly ``tau.theta(z)``. >>> from flint import acb, acb_mat, showgood >>> z = acb(1+1j); tau = acb(1.25+3j) From d9a0273c1c74cf0744452534d63776e35cab5d40 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Sat, 15 Jun 2024 17:51:12 -0400 Subject: [PATCH 09/11] one more example, and now using showgood --- src/flint/types/acb_theta.pyx | 29 +++++++++++++++++++++++------ 1 file changed, 23 insertions(+), 6 deletions(-) diff --git a/src/flint/types/acb_theta.pyx b/src/flint/types/acb_theta.pyx index ea8e4907..b3294bb7 100644 --- a/src/flint/types/acb_theta.pyx +++ b/src/flint/types/acb_theta.pyx @@ -21,12 +21,11 @@ def acb_theta(acb_mat z, acb_mat tau, ulong square=False): >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])).entries() >>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])]) [+/- 3.82e-14] - >>> for i in range(4):showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]])).entries()[i], dps=25) - ... - 0.9694430387796704100046143 - 0.03055696120816803328582847j - 1.030556961196006476576271 + 0.03055696120816803328582847j - -1.220790267576967690128359 - 1.827055516791154669091679j - -1.820235910124989594900076 + 1.216251950154477951760042j + >>> showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]])).transpose(), dps=25) + [0.9694430387796704100046143 - 0.03055696120816803328582847j] + [ 1.030556961196006476576271 + 0.03055696120816803328582847j] + [ -1.220790267576967690128359 - 1.827055516791154669091679j] + [ -1.820235910124989594900076 + 1.216251950154477951760042j] >>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])).transpose() [ [1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j] [ [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j] @@ -44,6 +43,24 @@ def acb_theta(acb_mat z, acb_mat tau, ulong square=False): [ [+/- 2.43e-16] + [+/- 2.43e-16]j] [ [+/- 2.43e-16] + [+/- 2.43e-16]j] [ [+/- 2.43e-16] + [+/- 2.43e-16]j] + >>> ctx.prec = 10000 + >>> print(acb_mat([[1j, 0],[0,1j]]).theta(acb_mat([[0],[0]])).transpose().str(25)) + [ [1.180340599016096226045338 +/- 5.95e-26] + [+/- 1.23e-3010]j] + [[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j] + [[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j] + [[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j] + [[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + [[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + [[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j] + [[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + [[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] """ g = tau.nrows() From 920d48d436f5bb089c30d02cec5ee4d9ffe21c8f Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Sat, 15 Jun 2024 17:56:44 -0400 Subject: [PATCH 10/11] fix doctest --- src/flint/types/acb_theta.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/flint/types/acb_theta.pyx b/src/flint/types/acb_theta.pyx index b3294bb7..3e99cc3a 100644 --- a/src/flint/types/acb_theta.pyx +++ b/src/flint/types/acb_theta.pyx @@ -16,7 +16,7 @@ def acb_theta(acb_mat z, acb_mat tau, ulong square=False): This should be used via the method :meth:`.acb_mat.theta`, explicitly ``tau.theta(z)``. - >>> from flint import acb, acb_mat, showgood + >>> from flint import acb, acb_mat, showgood, ctx >>> z = acb(1+1j); tau = acb(1.25+3j) >>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])).entries() >>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])]) From 778de4b42104c829800ae8b143243ec3933ee840 Mon Sep 17 00:00:00 2001 From: Edgar Costa Date: Sat, 15 Jun 2024 18:25:10 -0400 Subject: [PATCH 11/11] use ctx.default() --- src/flint/types/acb_theta.pyx | 1 + 1 file changed, 1 insertion(+) diff --git a/src/flint/types/acb_theta.pyx b/src/flint/types/acb_theta.pyx index 3e99cc3a..52cca748 100644 --- a/src/flint/types/acb_theta.pyx +++ b/src/flint/types/acb_theta.pyx @@ -61,6 +61,7 @@ def acb_theta(acb_mat z, acb_mat tau, ulong square=False): [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] [ [+/- 1.23e-3010] + [+/- 1.23e-3010]j] + >>> ctx.default() """ g = tau.nrows()