This file contains functions useful for computing discrete Fourier
transforms and probability distribution functions for discrete random
variables for sequences of elements of or
, indexed by a
range(N),
, an abelian group, the conjugacy classes
of a permutation group, or the conjugacy classes of a matrix group.
This file implements:
Todo
AUTHORS:
Bases: sage.structure.sage_object.SageObject
An indexed sequence.
INPUT:
This just returns the common parent of the
list
elements. In some applications (say, when computing the
discrete Fourier transform, dft), it is more accurate to think
of the base_ring as the group ring
.
EXAMPLES:
sage: J = range(10)
sage: A = [1/10 for j in J]
sage: s = IndexedSequence(A,J)
sage: s.base_ring()
Rational Field
Convolves two sequences of the same length (automatically expands the shortest one by extending it by 0 if they have different lengths).
If and
are sequences indexed by
,
extended by zero for all
in
, then the convolution is
INPUT:
OUTPUT:
The Dirichlet convolution of self and other.
EXAMPLES:
sage: J = range(5)
sage: A = [ZZ(1) for i in J]
sage: B = [ZZ(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: t = IndexedSequence(B,J)
sage: s.convolution(t)
[1, 2, 3, 4, 5, 4, 3, 2, 1]
AUTHOR: David Joyner (2006-09)
Convolves two collections indexed by a range(...) of the same length (automatically expands the shortest one by extending it by 0 if they have different lengths).
If and
are sequences indexed by
,
extended periodically for all
in
, then the convolution is
INPUT:
OUTPUT:
The Dirichlet convolution of self and other.
EXAMPLES:
sage: I = range(5)
sage: A = [ZZ(1) for i in I]
sage: B = [ZZ(1) for i in I]
sage: s = IndexedSequence(A,I)
sage: t = IndexedSequence(B,I)
sage: s.convolution_periodic(t)
[5, 5, 5, 5, 5, 5, 5, 5, 5]
AUTHOR: David Joyner (2006-09)
A discrete Cosine transform.
EXAMPLES:
sage: J = range(5)
sage: A = [exp(-2*pi*i*I/5) for i in J]
sage: s = IndexedSequence(A,J)
sage: s.dct()
Indexed sequence: [1/4*(sqrt(5) - 1)*e^(-2/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-4/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-6/5*I*pi) + 1/4*(sqrt(5) - 1)*e^(-8/5*I*pi) + 1, 1/4*(sqrt(5) - 1)*e^(-2/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-4/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-6/5*I*pi) + 1/4*(sqrt(5) - 1)*e^(-8/5*I*pi) + 1, 1/4*(sqrt(5) - 1)*e^(-2/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-4/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-6/5*I*pi) + 1/4*(sqrt(5) - 1)*e^(-8/5*I*pi) + 1, 1/4*(sqrt(5) - 1)*e^(-2/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-4/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-6/5*I*pi) + 1/4*(sqrt(5) - 1)*e^(-8/5*I*pi) + 1, 1/4*(sqrt(5) - 1)*e^(-2/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-4/5*I*pi) - 1/4*(sqrt(5) + 1)*e^(-6/5*I*pi) + 1/4*(sqrt(5) - 1)*e^(-8/5*I*pi) + 1]
indexed by [0, 1, 2, 3, 4]
A discrete Fourier transform “over ” using exact
-th roots of unity.
EXAMPLES:
sage: J = range(6)
sage: A = [ZZ(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: s.dft(lambda x:x^2)
Indexed sequence: [6, 0, 0, 6, 0, 0]
indexed by [0, 1, 2, 3, 4, 5]
sage: s.dft()
Indexed sequence: [6, 0, 0, 0, 0, 0]
indexed by [0, 1, 2, 3, 4, 5]
sage: G = SymmetricGroup(3)
sage: J = G.conjugacy_classes_representatives()
sage: s = IndexedSequence([1,2,3],J) # 1,2,3 are the values of a class fcn on G
sage: s.dft() # the "scalar-valued Fourier transform" of this class fcn
Indexed sequence: [8, 2, 2]
indexed by [(), (1,2), (1,2,3)]
sage: J = AbelianGroup(2,[2,3],names='ab')
sage: s = IndexedSequence([1,2,3,4,5,6],J)
sage: s.dft() # the precision of output is somewhat random and architecture dependent.
Indexed sequence: [21.0000000000000, -2.99999999999997 - 1.73205080756885*I, -2.99999999999999 + 1.73205080756888*I, -9.00000000000000 + 0.0000000000000485744257349999*I, -0.00000000000000976996261670137 - 0.0000000000000159872115546022*I, -0.00000000000000621724893790087 - 0.0000000000000106581410364015*I]
indexed by Multiplicative Abelian group isomorphic to C2 x C3
sage: J = CyclicPermutationGroup(6)
sage: s = IndexedSequence([1,2,3,4,5,6],J)
sage: s.dft() # the precision of output is somewhat random and architecture dependent.
Indexed sequence: [21.0000000000000, -2.99999999999997 - 1.73205080756885*I, -2.99999999999999 + 1.73205080756888*I, -9.00000000000000 + 0.0000000000000485744257349999*I, -0.00000000000000976996261670137 - 0.0000000000000159872115546022*I, -0.00000000000000621724893790087 - 0.0000000000000106581410364015*I]
indexed by Cyclic group of order 6 as a permutation group
sage: p = 7; J = range(p); A = [kronecker_symbol(j,p) for j in J]
sage: s = IndexedSequence(A,J)
sage: Fs = s.dft()
sage: c = Fs.list()[1]; [x/c for x in Fs.list()]; s.list()
[0, 1, 1, -1, 1, -1, -1]
[0, 1, 1, -1, 1, -1, -1]
The DFT of the values of the quadratic residue symbol is itself, up to a constant factor (denoted c on the last line above).
Todo
Read the parent of the elements of S; if or
leave as
is; if AbelianGroup, use abelian_group_dual; if some other
implemented Group (permutation, matrix), call .characters()
and test if the index list is the set of conjugacy classes.
Return a python dict of self where the keys are elments in the indexing set.
EXAMPLES:
sage: J = range(10)
sage: A = [1/10 for j in J]
sage: s = IndexedSequence(A,J)
sage: s.dict()
{0: 1/10, 1: 1/10, 2: 1/10, 3: 1/10, 4: 1/10, 5: 1/10, 6: 1/10, 7: 1/10, 8: 1/10, 9: 1/10}
A discrete Sine transform.
EXAMPLES:
sage: J = range(5)
sage: I = CC.0; pi = CC(pi)
sage: A = [exp(-2*pi*i*I/5) for i in J]
sage: s = IndexedSequence(A,J)
sage: s.dst() # discrete sine
Indexed sequence: [1.11022302462516e-16 - 2.50000000000000*I, 1.11022302462516e-16 - 2.50000000000000*I, 1.11022302462516e-16 - 2.50000000000000*I, 1.11022302462516e-16 - 2.50000000000000*I, 1.11022302462516e-16 - 2.50000000000000*I]
indexed by [0, 1, 2, 3, 4]
Wraps the gsl WaveletTransform.forward in dwt (written by Joshua Kantor). Assumes the length of the sample is a power of 2. Uses the GSL function gsl_wavelet_transform_forward().
INPUT:
other – the the name of the type of wavelet; valid choices are:
wavelet_k – For daubechies wavelets, wavelet_k specifies a
daubechie wavelet with vanishing moments.
for
even are the only ones implemented.
For Haar wavelets, wavelet_k must be 2.
For bspline wavelets, wavelet_k equal to will give biorthogonal B-spline wavelets
of order
where wavelet_k equals
.
The wavelet transform uses levels.
EXAMPLES:
sage: J = range(8)
sage: A = [RR(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: t = s.dwt()
sage: t # slightly random output
Indexed sequence: [2.82842712474999, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000]
indexed by [0, 1, 2, 3, 4, 5, 6, 7]
Wraps the gsl FastFourierTransform.forward() in fft.
If the length is a power of 2 then this automatically uses the radix2 method. If the number of sample points in the input is a power of 2 then the wrapper for the GSL function gsl_fft_complex_radix2_forward() is automatically called. Otherwise, gsl_fft_complex_forward() is used.
EXAMPLES:
sage: J = range(5)
sage: A = [RR(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: t = s.fft(); t
Indexed sequence: [5.00000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000]
indexed by [0, 1, 2, 3, 4]
A discrete inverse Fourier transform. Only works over .
EXAMPLES:
sage: J = range(5)
sage: A = [ZZ(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: fs = s.dft(); fs
Indexed sequence: [5, 0, 0, 0, 0]
indexed by [0, 1, 2, 3, 4]
sage: it = fs.idft(); it
Indexed sequence: [1, 1, 1, 1, 1]
indexed by [0, 1, 2, 3, 4]
sage: it == s
True
Implements the gsl WaveletTransform.backward() in dwt.
Assumes the length of the sample is a power of 2. Uses the GSL function gsl_wavelet_transform_backward().
INPUT:
other – Must be one of the following:
wavelet_k – For daubechies wavelets, wavelet_k specifies a
daubechie wavelet with vanishing moments.
for
even are the only ones implemented.
For Haar wavelets, wavelet_k must be 2.
For bspline wavelets, wavelet_k equal to will give biorthogonal B-spline wavelets
of order
where wavelet_k equals
.
EXAMPLES:
sage: J = range(8)
sage: A = [RR(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: t = s.dwt()
sage: t # random arch dependent output
Indexed sequence: [2.82842712474999, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000]
indexed by [0, 1, 2, 3, 4, 5, 6, 7]
sage: t.idwt() # random arch dependent output
Indexed sequence: [1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000]
indexed by [0, 1, 2, 3, 4, 5, 6, 7]
sage: t.idwt() == s
True
sage: J = range(16)
sage: A = [RR(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: t = s.dwt("bspline", 103)
sage: t # random arch dependent output
Indexed sequence: [4.00000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000]
indexed by [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15]
sage: t.idwt("bspline", 103) == s
True
Implements the gsl FastFourierTransform.inverse in fft.
If the number of sample points in the input is a power of 2 then the wrapper for the GSL function gsl_fft_complex_radix2_inverse() is automatically called. Otherwise, gsl_fft_complex_inverse() is used.
EXAMPLES:
sage: J = range(5)
sage: A = [RR(1) for i in J]
sage: s = IndexedSequence(A,J)
sage: t = s.fft(); t
Indexed sequence: [5.00000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000, 0.000000000000000]
indexed by [0, 1, 2, 3, 4]
sage: t.ifft()
Indexed sequence: [1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000, 1.00000000000000]
indexed by [0, 1, 2, 3, 4]
sage: t.ifft() == s
1
Return the indexing object.
EXAMPLES:
sage: J = range(10)
sage: A = [1/10 for j in J]
sage: s = IndexedSequence(A,J)
sage: s.index_object()
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9]
Return the list of self.
EXAMPLES:
sage: J = range(10)
sage: A = [1/10 for j in J]
sage: s = IndexedSequence(A,J)
sage: s.list()
[1/10, 1/10, 1/10, 1/10, 1/10, 1/10, 1/10, 1/10, 1/10, 1/10]
Plot the points of the sequence.
Elements of the sequence are assumed to be real or from a finite field, with a real indexing set I = range(len(self)).
EXAMPLES:
sage: I = range(3)
sage: A = [ZZ(i^2)+1 for i in I]
sage: s = IndexedSequence(A,I)
sage: P = s.plot()
sage: show(P) # Not tested
Plot the histogram plot of the sequence.
The sequence is assumed to be real or from a finite field,
with a real indexing set I coercible into .
Options are clr, which is an RGB value, and eps, which is the spacing between the bars.
EXAMPLES:
sage: J = range(3)
sage: A = [ZZ(i^2)+1 for i in J]
sage: s = IndexedSequence(A,J)
sage: P = s.plot_histogram()
sage: show(P) # Not tested