moment_functionals¶
Exact moment functionals for circles, ellipses and other varieties with respect to uniform measure.
EXAMPLES:
sage: from momentproblems.moment_functionals import *
sage: moment_ellipse = mk_moment_ellipse(center=(1,1), radius=(5,4))
sage: h0 = hankel(moment_circle, dim=2, d=2)
sage: h1 = hankel(moment_ellipse, dim=2, d=2)
sage: h0.right_kernel_matrix()
[ 1 0 0 -1 0 -1]
sage: h1.right_kernel_matrix()
[ 1 32/359 50/359 -16/359 0 -25/359]
sage: mom_a = mk_moment_torus_line(ratio=(1, 1), translate=e^(i*pi/4))
sage: mom_b = mk_moment_torus_line(ratio=(-1, 1), translate=1)
sage: mom_c = mk_moment_torus_line(ratio=(-2, 1), translate=1)
sage: mom_d = mk_moment_torus_line(ratio=(2, 1), translate=e^(i*1))
sage: mom_e = mk_moment_torus_line(ratio=(2, 1), translate=2)
sage: T_c = toeplitz(mom_c, dim=2, d=3)
sage: T_c.is_hermitian()
True
sage: T_e = toeplitz(mom_e, dim=2, d=3)
sage: T_e.is_hermitian()
False
The Toeplitz moment matrices have kernels which allow to recover the expected ideals:
sage: basis = vector(monomial_basis_laurent(dim=2, d=3))
sage: ker = T_c.right_kernel_matrix() * basis
sage: basis.base_ring().ideal(list(ker)).groebner_basis()
(z0^2*z1 - 1,)
sage: ker = T_e.right_kernel_matrix() * basis
sage: basis.base_ring().ideal(list(ker)).groebner_basis()
(z0^2 - 1/2*z1,)
Affine line segments:
sage: mk_moment_affine_line_segment((0,0,0), (1,2,3), normalize=False)((1,2,1)).radical_expression()
12/5*sqrt(14)
sage: mk_moment_affine_line_segment((0,0,0), (1,2,3))((1,2,1))
12/5
sage: mk_moment_affine_line_segment((1,0,0), (2,0,0))((2,0,0)) == 7/3
True
sage: mk_moment_affine_line_segment((0,1,0), (0,2,0))((0,2,0)) == 7/3
True
The kernel of the Hankel matrix recovers the ideal corresponding to the line defined by start and end points:
sage: start, end = (2,3,1), (3,4,0)
sage: mom = mk_moment_affine_line_segment(start, end)
sage: basis = vector(monomial_basis(dim=3, d=2))
sage: H = hankel(mom, dim=3, d=2)
sage: ker = H.right_kernel_matrix() * basis
sage: gb = basis.base_ring().ideal(list(ker)).groebner_basis(); gb
[y0 + y2 - 3, y1 + y2 - 4]
sage: all(f(*pnt).is_zero() for f in gb for pnt in [start, end])
True
We construct a moment matrix corresponding to points lying on a trigonometric curve. We normalize by the number of points, so that the 0-th moment becomes 1:
sage: R.<x1,x2> = QQ[]
sage: g = -1/4 + x1 + x2 + x1*x2
sage: points = list(generate_points_on_trigonometric_variety(R.ideal(g), gridsize=6, parent=RDF))
sage: basis = monomial_basis_laurent(dim=R.ngens(), d=2)
sage: A = vandermonde(points, basis, trigonometric=True)
sage: A.dimensions()
(18, 9)
sage: T = A.H * A / len(points)
sage: T.is_hermitian()
True
sage: T.zero_at(1e-8)[:4,:4].n(30)
[ 1.0000000 0.069444444 -0.55381944 0.37037037 + 0.11111111*I]
[ 0.069444444 1.0000000 0.069444444 -0.18981481 + 0.027777778*I]
[ -0.55381944 0.069444444 1.0000000 -0.40219907 - 0.097222222*I]
[ 0.37037037 - 0.11111111*I -0.18981481 - 0.027777778*I -0.40219907 + 0.097222222*I 1.0000000]
sage: T.zero_at(1e-8)[-4:,-4:].n(30)
[ 1.0000000 -0.40219907 - 0.097222222*I -0.18981481 + 0.027777778*I 0.37037037 + 0.11111111*I]
[-0.40219907 + 0.097222222*I 1.0000000 0.069444444 -0.55381944]
[-0.18981481 - 0.027777778*I 0.069444444 1.0000000 0.069444444]
[ 0.37037037 - 0.11111111*I -0.55381944 0.069444444 1.0000000]
sage: T.zero_at(1e-12).diagonal() # abs tol 1e-12
[1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]
The matrix has one-dimensional kernel, as the relation defined by g is in the kernel:
sage: from numpy.linalg import matrix_rank
sage: T.ncols() - matrix_rank(T, tol=1e-13)
1
sage: ker = T.SVD()[0][:,-1]
sage: (ker / ker[0,0]).zero_at(1e-12) # abs tol 1e-12
[ 1.0]
[ 2.0]
[ 1.0]
[ 2.0]
[-1.0]
[ 2.0]
[ 1.0]
[ 2.0]
[ 1.0]
- moment_functionals.eval_functional(moment_functional, p)¶
- moment_functionals.mk_moment_ellipse(center, radius)¶
- moment_functionals.monomial_basis(dim, d=None)¶
- moment_functionals.hankel(moment_functional, dim=None, d=None, basis=None)¶
- moment_functionals.mk_moment_torus_line(ratio, translate)¶
Construct a moment functional of a straight line on the two-dimensional torus.
ratio = (q1, q2) – coprime translate = c1 – non-zero scalar (if on unit circle, then the Toeplitz moment matrix is Hermitian)
where:
q2*t2 == q1*t1 + b1 | exp(i*_) exp(i*q2*t2) == exp(i*q1*t1) * exp(i*b1) | c1 = exp(i*b1) z2^q2 == z1^q1 * c1
- moment_functionals.monomial_basis_laurent(dim, d)¶
- moment_functionals.toeplitz(moment_functional, dim, d)¶
- moment_functionals.mk_moment_affine_line_segment(start, end, *, normalize=True)¶
- moment_functionals.mk_moment_trigonometric_curve_example(c1, c2, *, prec=80)¶
Return the moment functional for the uniform measure on a trigonometric example curve.
The curve depends on two parameters, defining a family of curves. The curve is defined in degrees -1..+1.
If c1 > c2, the curve is connected, if c1 < c2 it has two disconnected parts. The parameters must satisfy c1 != c2 and c1 < c2 + 2.
EXAMPLES:
sage: from momentproblems.moment_functionals import * sage: mom = mk_moment_trigonometric_curve_example(5/8, 0) sage: [mom(a) for a in [(0,0), (1,0), (0,1)]] [[12.61257092367025016329 +/- ...], [2.753059781992467887758 +/- ...] + [+/- ...]*I, [2.753059781992467887758 +/- ...] + [+/- ...]*I]
sage: mom = mk_moment_trigonometric_curve_example(1, 5/4) sage: [mom(a) for a in [(0,0), (1,0), (0,1)]] [[13.448678364566850580435 +/- ...], [-0.408926940531123205236 +/- ...] + [+/- ...]*I, [0.060306710833167660509 +/- ...] + [+/- ...]*I]
- moment_functionals.mk_moment_trigonometric_curve_example1(c1, c2, *, prec=80)¶
A parametric family of curves of max-degree 1 on the torus.
Requirements: c1,c2 real such that c1+c2 != 0.
Orientation changes at c1==c2.
- moment_functionals.modulated(functional, shift)¶
Apply a modulation to a moment functional.
EXAMPLES:
sage: from momentproblems.moment_functionals import * sage: mom = mk_moment_trigonometric_curve_example(5/8, 0) sage: K = mom((0,0)).parent() sage: mom2 = modulated(mom, vector(K, (exp(2*pi*I*0.1), exp(2*pi*I*0.3)))) sage: (mom((0,0)) - mom2((0,0))).diameter() < 1e-13 True sage: (mom((3,2)) * K(exp(2*pi*I*(0.1*3 + 0.3*2))) - mom2((3,2))).diameter() < 1e-12 True
- moment_functionals.generate_points_on_trigonometric_variety(I, gridsize, parent=None)¶
Find points on trigonometric variety by cutting with (deterministic) hyperplanes.
We assume that the variety can be defined purely in terms of coordinates x_k, i.e. coordinates of the real part, where
\[z_k = exp(i t_k) = x_k + i y_k = cos(t_k) + i sin(t_k).\]INPUT:
gridsize– number of grid points -1 at which the hyperplanes are evaluatedparent– ring of the coordinates of the result (default:SR)
OUTPUT:
Iterator of some points on the variety, given with respect to the coordinates t_k.
EXAMPLES:
sage: from momentproblems.moment_functionals import * sage: P.<x0,x1,x2> = QQ[] sage: g = -1/4 + x1 + x2 + x1*x2 + x0 sage: points = list(generate_points_on_trigonometric_variety(P.ideal(g), gridsize=4)) sage: len(points) 102 sage: points[:5] [(pi, 1/3*pi, 1/3*pi), (pi, 1/3*pi, -1/3*pi), (pi, -1/3*pi, 1/3*pi), (pi, -1/3*pi, -1/3*pi), (-pi, 1/3*pi, 1/3*pi)]
sage: from itertools import islice sage: list(islice(generate_points_on_trigonometric_variety(P.ideal(g, x2-1/2), gridsize=5, parent=RealBallField(100)), 5)) [([3.14159265358979323846264338328 +/- 2.25e-30], [1.04719755119659774615421446109 +/- 5.86e-30], [1.04719755119659774615421446109 +/- 5.86e-30]), ([3.14159265358979323846264338328 +/- 2.25e-30], [1.04719755119659774615421446109 +/- 5.86e-30], [-1.04719755119659774615421446109 +/- 5.86e-30]), ([3.14159265358979323846264338328 +/- 2.25e-30], [-1.04719755119659774615421446109 +/- 5.86e-30], [1.04719755119659774615421446109 +/- 5.86e-30]), ([3.14159265358979323846264338328 +/- 2.25e-30], [-1.04719755119659774615421446109 +/- 5.86e-30], [-1.04719755119659774615421446109 +/- 5.86e-30]), ([-3.14159265358979323846264338328 +/- 2.25e-30], [1.04719755119659774615421446109 +/- 5.86e-30], [1.04719755119659774615421446109 +/- 5.86e-30])]
- moment_functionals.vandermonde(points, basis, trigonometric=False)¶
Construct the Vandermonde matrix of size
len(points) × len(basis).EXAMPLES:
sage: from momentproblems.moment_functionals import * sage: points = [(1/2, 1/2), (1/3, 2/5)] sage: basis = monomial_basis_laurent(dim=2, d=1) sage: V = vandermonde(points, basis, trigonometric=True); V.T [ 1.0 1.0] [ 0.8775825618903728 + 0.479425538604203*I 0.9210609940028851 + 0.3894183423086505*I] [ 0.8775825618903728 + 0.479425538604203*I 0.9449569463147377 + 0.3271946967961522*I] [0.5403023058681398 + 0.8414709848078965*I 0.742947367824044 + 0.6693498402504662*I]
- moment_functionals.random_positive_definite_matrix(K, n)¶
EXAMPLES:
sage: from momentproblems.moment_functionals import * sage: random_positive_definite_matrix(QQ, 5).is_positive_semidefinite() True sage: random_positive_definite_matrix(RDF, 5).is_positive_definite() True sage: random_positive_definite_matrix(CDF, 5).is_positive_definite() True
- moment_functionals.random_moment_functional(ideal, *, random_fun=None)¶
Return a random functional on a coordinate ring of a variety defined by some ideal.
EXAMPLES:
sage: from momentproblems.moment_functionals import * sage: R.<x,y,z> = QQ[] sage: J = R.ideal(y^2 - x^3 + 1, z^4 - z) sage: moment_fun = random_moment_functional(J, random_fun=lambda: ZZ.random_element(1, 10)) sage: H = hankel(moment_fun, dim=R.ngens(), d=4) sage: J == R.ideal(list(H.right_kernel_matrix() * vector(monomial_basis(dim=R.ngens(), d=4)))) True
Note that the matrices are not positive-semidefinite:
sage: H.is_positive_semidefinite() False
Laurent case:
sage: L.<x,y,z> = LaurentPolynomialRing(QQ) sage: J = L.ideal([y^2 - x^3 + 1, z^4 - z]) sage: moment_fun = random_moment_functional(J, random_fun=lambda: ZZ.random_element(1, 10)) sage: H = toeplitz(moment_fun, dim=L.ngens(), d=3) sage: J == L.ideal(list(H.right_kernel_matrix() * vector(monomial_basis_laurent(dim=L.ngens(), d=3)))) True sage: H.is_positive_semidefinite() False
