Skip to content

Commit 33a4daa

Browse files
oscarbenjaminoscarbenjamin
authored
Merge pull request #142 from edgarcosta/acb_theta
Exposing acb_theta_all
2 parents ee239b1 + 778de4b commit 33a4daa

File tree

8 files changed

+267
-1
lines changed

8 files changed

+267
-1
lines changed

‎doc/source/acb_theta.rst‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,5 @@
1+
**acb_theta** -- Riemann theta functions
2+
===============================================================================
3+
4+
.. autofunction :: flint.types.acb_theta.acb_theta
5+

‎doc/source/index.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -61,6 +61,7 @@ Matrix types
6161
fmpz_mod_mat.rst
6262
arb_mat.rst
6363
acb_mat.rst
64+
acb_theta.rst
6465

6566
Polynomial types
6667
................

‎meson.build‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -22,6 +22,9 @@ endif
2222
# flint.pc was missing -lflint until Flint 3.1.0
2323
if flint_dep.version().version_compare('<3.1')
2424
flint_dep = cc.find_library('flint')
25+
have_acb_theta = false
26+
else
27+
have_acb_theta = true
2528
endif
2629

2730
pyflint_deps = [dep_py, gmp_dep, mpfr_dep, flint_dep]

‎src/flint/flintlib/acb_theta.pxd‎

Lines changed: 136 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,136 @@
1+
from flint.flintlib.acb cimport acb_t, acb_srcptr, acb_ptr
2+
from flint.flintlib.acb_poly cimport acb_poly_struct, acb_poly_t
3+
from flint.flintlib.arb cimport arb_t, arb_ptr, arb_srcptr
4+
from flint.flintlib.flint cimport ulong, slong, flint_rand_t
5+
from flint.flintlib.fmpz_mat cimport fmpz_mat_t, fmpz_mat_struct
6+
from flint.flintlib.acb_mat cimport acb_mat_t
7+
from flint.flintlib.arb_mat cimport arb_mat_t
8+
from flint.flintlib.arf cimport arf_t
9+
10+
11+
cdef extern from "flint/acb_theta.h":
12+
ctypedef struct acb_theta_eld_struct:
13+
slong dim
14+
slong ambient_dim
15+
slong * last_coords
16+
slong min
17+
slong mid
18+
slong max
19+
slong nr
20+
slong nl
21+
acb_theta_eld_struct * rchildren
22+
acb_theta_eld_struct * lchildren
23+
slong nb_pts, nb_border
24+
slong * box
25+
ctypedef acb_theta_eld_struct acb_theta_eld_t[1]
26+
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)
27+
ctypedef int (*acb_theta_ql_worker_t)(acb_ptr, acb_srcptr, acb_srcptr, arb_srcptr, arb_srcptr, const acb_mat_t, slong, slong)
28+
29+
void acb_theta_all(acb_ptr th, acb_srcptr z, const acb_mat_t tau, int sqr, slong prec)
30+
void acb_theta_naive_fixed_ab(acb_ptr th, ulong ab, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec)
31+
void acb_theta_naive_all(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec)
32+
void acb_theta_jet_all(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec)
33+
void acb_theta_jet_naive_fixed_ab(acb_ptr dth, ulong ab, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec)
34+
void acb_theta_jet_naive_all(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec)
35+
slong sp2gz_dim(const fmpz_mat_t mat)
36+
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)
37+
void sp2gz_j(fmpz_mat_t mat)
38+
void sp2gz_block_diag(fmpz_mat_t mat, const fmpz_mat_t U)
39+
void sp2gz_trig(fmpz_mat_t mat, const fmpz_mat_t S)
40+
void sp2gz_embed(fmpz_mat_t res, const fmpz_mat_t mat)
41+
void sp2gz_restrict(fmpz_mat_t res, const fmpz_mat_t mat)
42+
slong sp2gz_nb_fundamental(slong g)
43+
void sp2gz_fundamental(fmpz_mat_t mat, slong j)
44+
int sp2gz_is_correct(const fmpz_mat_t mat)
45+
int sp2gz_is_j(const fmpz_mat_t mat)
46+
int sp2gz_is_block_diag(const fmpz_mat_t mat)
47+
int sp2gz_is_trig(const fmpz_mat_t mat)
48+
int sp2gz_is_embedded(fmpz_mat_t res, const fmpz_mat_t mat)
49+
void sp2gz_inv(fmpz_mat_t inv, const fmpz_mat_t mat)
50+
fmpz_mat_struct * sp2gz_decompose(slong * nb, const fmpz_mat_t mat)
51+
void sp2gz_randtest(fmpz_mat_t mat, flint_rand_t state, slong bits)
52+
void acb_siegel_cocycle(acb_mat_t c, const fmpz_mat_t mat, const acb_mat_t tau, slong prec)
53+
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)
54+
void acb_siegel_transform(acb_mat_t w, const fmpz_mat_t mat, const acb_mat_t tau, slong prec)
55+
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)
56+
void acb_siegel_cho(arb_mat_t C, const acb_mat_t tau, slong prec)
57+
void acb_siegel_yinv(arb_mat_t Yinv, const acb_mat_t tau, slong prec)
58+
void acb_siegel_reduce(fmpz_mat_t mat, const acb_mat_t tau, slong prec)
59+
int acb_siegel_is_reduced(const acb_mat_t tau, slong tol_exp, slong prec)
60+
void acb_siegel_randtest(acb_mat_t tau, flint_rand_t state, slong prec, slong mag_bits)
61+
void acb_siegel_randtest_reduced(acb_mat_t tau, flint_rand_t state, slong prec, slong mag_bits)
62+
void acb_siegel_randtest_vec(acb_ptr z, flint_rand_t state, slong g, slong prec)
63+
void acb_theta_char_get_slong(slong * n, ulong a, slong g)
64+
ulong acb_theta_char_get_a(const slong * n, slong g)
65+
void acb_theta_char_get_arb(arb_ptr v, ulong a, slong g)
66+
void acb_theta_char_get_acb(acb_ptr v, ulong a, slong g)
67+
slong acb_theta_char_dot(ulong a, ulong b, slong g)
68+
slong acb_theta_char_dot_slong(ulong a, const slong * n, slong g)
69+
void acb_theta_char_dot_acb(acb_t x, ulong a, acb_srcptr z, slong g, slong prec)
70+
int acb_theta_char_is_even(ulong ab, slong g)
71+
int acb_theta_char_is_goepel(ulong ch1, ulong ch2, ulong ch3, ulong ch4, slong g)
72+
int acb_theta_char_is_syzygous(ulong ch1, ulong ch2, ulong ch3, slong g)
73+
void acb_theta_eld_init(acb_theta_eld_t E, slong d, slong g)
74+
void acb_theta_eld_clear(acb_theta_eld_t E)
75+
int acb_theta_eld_set(acb_theta_eld_t E, const arb_mat_t C, const arf_t R2, arb_srcptr v)
76+
void acb_theta_eld_points(slong * pts, const acb_theta_eld_t E)
77+
void acb_theta_eld_border(slong * pts, const acb_theta_eld_t E)
78+
int acb_theta_eld_contains(const acb_theta_eld_t E, slong * pt)
79+
void acb_theta_eld_print(const acb_theta_eld_t E)
80+
void acb_theta_naive_radius(arf_t R2, arf_t eps, const arb_mat_t C, slong ord, slong prec)
81+
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)
82+
void acb_theta_naive_term(acb_t res, acb_srcptr z, const acb_mat_t tau, slong * tup, slong * n, slong prec)
83+
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)
84+
void acb_theta_naive_00(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec)
85+
void acb_theta_naive_0b(acb_ptr th, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec)
86+
void acb_theta_naive_fixed_a(acb_ptr th, ulong a, acb_srcptr zs, slong nb, const acb_mat_t tau, slong prec)
87+
slong acb_theta_jet_nb(slong ord, slong g)
88+
slong acb_theta_jet_total_order(const slong * tup, slong g)
89+
void acb_theta_jet_tuples(slong * tups, slong ord, slong g)
90+
slong acb_theta_jet_index(const slong * tup, slong g)
91+
void acb_theta_jet_mul(acb_ptr res, acb_srcptr v1, acb_srcptr v2, slong ord, slong g, slong prec)
92+
void acb_theta_jet_compose(acb_ptr res, acb_srcptr v, const acb_mat_t N, slong ord, slong prec)
93+
void acb_theta_jet_exp_pi_i(acb_ptr res, arb_srcptr a, slong ord, slong g, slong prec)
94+
void acb_theta_jet_naive_radius(arf_t R2, arf_t eps, arb_srcptr v, const arb_mat_t C, slong ord, slong prec)
95+
void acb_theta_jet_naive_00(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec)
96+
void acb_theta_jet_error_bounds(arb_ptr err, acb_srcptr z, const acb_mat_t tau, acb_srcptr dth, slong ord, slong prec)
97+
void acb_theta_dist_pt(arb_t d, arb_srcptr v, const arb_mat_t C, slong * n, slong prec)
98+
void acb_theta_dist_lat(arb_t d, arb_srcptr v, const arb_mat_t C, slong prec)
99+
void acb_theta_dist_a0(arb_ptr d, acb_srcptr z, const acb_mat_t tau, slong prec)
100+
slong acb_theta_dist_addprec(const arb_t d)
101+
void acb_theta_agm_hadamard(acb_ptr res, acb_srcptr a, slong g, slong prec)
102+
void acb_theta_agm_sqrt(acb_ptr res, acb_srcptr a, acb_srcptr rts, slong nb, slong prec)
103+
void acb_theta_agm_mul(acb_ptr res, acb_srcptr a1, acb_srcptr a2, slong g, slong prec)
104+
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)
105+
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)
106+
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)
107+
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)
108+
slong acb_theta_ql_a0_nb_steps(const arb_mat_t C, slong s, slong prec)
109+
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)
110+
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)
111+
void acb_theta_ql_all(acb_ptr th, acb_srcptr z, const acb_mat_t tau, int sqr, slong prec)
112+
void acb_theta_jet_ql_bounds(arb_t c, arb_t rho, acb_srcptr z, const acb_mat_t tau, slong ord)
113+
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)
114+
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)
115+
void acb_theta_jet_ql_all(acb_ptr dth, acb_srcptr z, const acb_mat_t tau, slong ord, slong prec)
116+
ulong acb_theta_transform_char(slong * e, const fmpz_mat_t mat, ulong ab)
117+
void acb_theta_transform_sqrtdet(acb_t res, const acb_mat_t tau, slong prec)
118+
slong acb_theta_transform_kappa(acb_t sqrtdet, const fmpz_mat_t mat, const acb_mat_t tau, slong prec)
119+
slong acb_theta_transform_kappa2(const fmpz_mat_t mat)
120+
void acb_theta_transform_proj(acb_ptr res, const fmpz_mat_t mat, acb_srcptr th, int sqr, slong prec)
121+
void acb_theta_g2_jet_naive_1(acb_ptr dth, const acb_mat_t tau, slong prec)
122+
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)
123+
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)
124+
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)
125+
slong acb_theta_g2_character(const fmpz_mat_t mat)
126+
void acb_theta_g2_psi4(acb_t res, acb_srcptr th2, slong prec)
127+
void acb_theta_g2_psi6(acb_t res, acb_srcptr th2, slong prec)
128+
void acb_theta_g2_chi10(acb_t res, acb_srcptr th2, slong prec)
129+
void acb_theta_g2_chi12(acb_t res, acb_srcptr th2, slong prec)
130+
void acb_theta_g2_chi5(acb_t res, acb_srcptr th, slong prec)
131+
void acb_theta_g2_chi35(acb_t res, acb_srcptr th, slong prec)
132+
void acb_theta_g2_chi3_6(acb_poly_t res, acb_srcptr dth, slong prec)
133+
void acb_theta_g2_sextic(acb_poly_t res, const acb_mat_t tau, slong prec)
134+
void acb_theta_g2_sextic_chi5(acb_poly_t res, acb_t chi5, const acb_mat_t tau, slong prec)
135+
void acb_theta_g2_covariants(acb_poly_struct * res, const acb_poly_t f, slong prec)
136+
void acb_theta_g2_covariants_lead(acb_ptr res, const acb_poly_t f, slong prec)

‎src/flint/test/__main__.py‎

Lines changed: 6 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -80,9 +80,14 @@ def run_doctests(verbose=None):
8080
flint.types.acb_series,
8181
flint.types.dirichlet,
8282
flint.functions.showgood]
83+
try:
84+
from flint.types import acb_theta
85+
modules.append(acb_theta)
86+
except ImportError:
87+
pass
8388
results = [doctest.testmod(x) for x in modules]
8489
# ffmpz, tfmpz = doctest.testmod(flint._fmpz, verbose=verbose)
85-
# failed, total = doctest.testmod(flint.pyflint, verbose=verbose)
90+
# failed, total = doctest.testmod(flint.pyflint, verbose=verbose)
8691
return tuple(sum(res) for res in zip(*results))
8792

8893

‎src/flint/types/acb_mat.pyx‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -814,3 +814,22 @@ cdef class acb_mat(flint_mat):
814814
if left:
815815
return Elist, L
816816
return Elist, R
817+
818+
def theta(tau, z, square=False):
819+
r"""
820+
Computes the vector valued Riemann theta function
821+
`(\theta_{a,b}(z, \tau) : a, b \in \{0,1\}^{g})` or its squares,
822+
where `\tau` is given by ``self``.
823+
824+
This is a wrapper for :func:`.acb_theta.acb_theta`; see the
825+
documentation for that method for details and examples.
826+
This follows the same conventions of the C-function
827+
`acb_theta_all <https://flintlib.org/doc/acb_theta.html#c.acb_theta_all>`_
828+
for the ordering of the theta characteristics.
829+
830+
"""
831+
try:
832+
from .acb_theta import acb_theta
833+
except ImportError:
834+
raise NotImplementedError("acb_mat.theta needs Flint >= 3.1.0")
835+
return acb_theta(z, tau, square=square)

‎src/flint/types/acb_theta.pyx‎

Lines changed: 93 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,93 @@
1+
from flint.flint_base.flint_context cimport getprec
2+
from flint.types.acb cimport acb
3+
from flint.types.acb_mat cimport acb_mat
4+
from flint.flintlib.acb cimport *
5+
from flint.flintlib.acb_mat cimport *
6+
from flint.flintlib.acb_theta cimport *
7+
8+
def acb_theta(acb_mat z, acb_mat tau, ulong square=False):
9+
r"""
10+
Computes the vector valued Riemann theta function
11+
`(\theta_{a,b}(z, \tau) : a, b \in \{0,1\}^{g})` or its squares.
12+
13+
This is a wrapper for the C-function
14+
`acb_theta_all <https://flintlib.org/doc/acb_theta.html#c.acb_theta_all>`_
15+
and it follows the same conventions for the ordering of the theta characteristics.
16+
17+
This should be used via the method :meth:`.acb_mat.theta`, explicitly ``tau.theta(z)``.
18+
19+
>>> from flint import acb, acb_mat, showgood, ctx
20+
>>> z = acb(1+1j); tau = acb(1.25+3j)
21+
>>> t0, t1, t2, t3 = acb_mat([[tau]]).theta(acb_mat([[z]])).entries()
22+
>>> sum([abs(x) for x in acb_mat([z.modular_theta(tau)]) - acb_mat([[-t3,t2,t0,t1]])])
23+
[+/- 3.82e-14]
24+
>>> showgood(lambda: acb_mat([[tau]]).theta(acb_mat([[z]])).transpose(), dps=25)
25+
[0.9694430387796704100046143 - 0.03055696120816803328582847j]
26+
[ 1.030556961196006476576271 + 0.03055696120816803328582847j]
27+
[ -1.220790267576967690128359 - 1.827055516791154669091679j]
28+
[ -1.820235910124989594900076 + 1.216251950154477951760042j]
29+
>>> acb_mat([[1j,0],[0,2*1j]]).theta(acb_mat([[0],[0]])).transpose()
30+
[ [1.09049252082308 +/- 3.59e-15] + [+/- 2.43e-16]j]
31+
[ [1.08237710165638 +/- 4.15e-15] + [+/- 2.43e-16]j]
32+
[[0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j]
33+
[[0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j]
34+
[[0.451696791791346 +/- 5.46e-16] + [+/- 2.43e-16]j]
35+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
36+
[[0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j]
37+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
38+
[[0.916991251621117 +/- 6.30e-16] + [+/- 2.43e-16]j]
39+
[[0.910167024735558 +/- 7.93e-16] + [+/- 2.43e-16]j]
40+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
41+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
42+
[[0.379830212998946 +/- 4.47e-16] + [+/- 2.43e-16]j]
43+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
44+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
45+
[ [+/- 2.43e-16] + [+/- 2.43e-16]j]
46+
>>> ctx.prec = 10000
47+
>>> print(acb_mat([[1j, 0],[0,1j]]).theta(acb_mat([[0],[0]])).transpose().str(25))
48+
[ [1.180340599016096226045338 +/- 5.95e-26] + [+/- 1.23e-3010]j]
49+
[[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j]
50+
[[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j]
51+
[[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j]
52+
[[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j]
53+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
54+
[[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j]
55+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
56+
[[0.9925441784910574194770081 +/- 3.15e-26] + [+/- 1.23e-3010]j]
57+
[[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j]
58+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
59+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
60+
[[0.8346268416740731862814297 +/- 3.29e-26] + [+/- 1.23e-3010]j]
61+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
62+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
63+
[ [+/- 1.23e-3010] + [+/- 1.23e-3010]j]
64+
>>> ctx.default()
65+
66+
"""
67+
g = tau.nrows()
68+
assert tau.ncols() == g
69+
assert z.nrows() == g
70+
assert z.ncols() == 1
71+
72+
# convert input
73+
cdef acb_ptr zvec
74+
zvec = _acb_vec_init(g)
75+
cdef long i
76+
for i in range(g):
77+
acb_set(zvec + i, acb_mat_entry(z.val, i, 0))
78+
79+
# initialize the output
80+
cdef slong nb = 1 << (2 * g)
81+
cdef acb_ptr theta = _acb_vec_init(nb)
82+
83+
acb_theta_all(theta, zvec, tau.val, square, getprec())
84+
_acb_vec_clear(zvec, g)
85+
# copy the output
86+
res = []
87+
cdef acb r
88+
for i in range(nb):
89+
r = acb.__new__(acb)
90+
acb_set(r.val, theta + i)
91+
res.append(r)
92+
_acb_vec_clear(theta, nb)
93+
return acb_mat([res])

‎src/flint/types/meson.build‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -41,6 +41,10 @@ exts = [
4141
'fmpz_mpoly',
4242
]
4343

44+
if have_acb_theta
45+
exts += ['acb_theta']
46+
endif
47+
4448
py.install_sources(
4549
pyfiles,
4650
pure: false,

0 commit comments

Comments
 (0)