Python arith.divisors函数代码示例

    def submodule_generated_by_images(self, M):
        """
        Return the submodule of this ambient modular symbols space
        generated by the images under all degeneracy maps of M. The space M
        must have the same weight, sign, and group or character as this
        ambient space.

        EXAMPLES::

            sage: ModularSymbols(6, 12).submodule_generated_by_images(ModularSymbols(1,12))
            Modular Symbols subspace of dimension 12 of Modular Symbols space of dimension 22 for Gamma_0(6) of weight 12 with sign 0 over Rational Field
        """
        S = self.zero_submodule()
        if self.level() % M.level() == 0:
            D = arith.divisors(self.level() // M.level())
        elif M.level() % self.level() == 0:
            D = arith.divisors(M.level() // self.level())
        else:
            D = []
        for t in D:
            d = M.degeneracy_map(self, t)
            if d.codomain() != self:
                raise ArithmeticError, "incompatible spaces of modular symbols"
            S += d.image()

        if self.is_full_hecke_module(compute=False):
            S._is_full_hecke_module = True

        return S

def find_product_decomposition(g, k, lmbda=1):
    r"""
    Try to find a product decomposition construction for difference matrices.

    INPUT:

    - ``g,k,lmbda`` -- integers, parameters of the difference matrix

    OUTPUT:

    A pair of pairs ``(g1,lmbda),(g2,lmbda2)`` if Sage knows how to build
    `(g1,k,lmbda1)` and `(g2,k,lmbda2)` difference matrices and ``False``
    otherwise.

    EXAMPLES::

        sage: from sage.combinat.designs.difference_matrices import find_product_decomposition
        sage: find_product_decomposition(77,6)
        ((7, 1), (11, 1))
        sage: find_product_decomposition(616,7)
        ((7, 1), (88, 1))
        sage: find_product_decomposition(24,10)
        False
    """
    for lmbda1 in divisors(lmbda):
        lmbda2 = lmbda//lmbda1

        # To avoid infinite loop:
        # if lmbda1 == lmbda, then g1 should not be g
        # if lmbda2 == lmbda, then g2 should not be g
        if lmbda1 == lmbda:
            if lmbda2 == lmbda:
                div = divisors(g)[1:-1]
            else:
                div = divisors(g)[:-1]
        else:
            if lmbda2 == lmbda:
                div = divisors(g)[1:]
            else:
                div = divisors(g)

        for g1 in div:
            g2 = g//g1
            if g1 > g2:
                break
            if (difference_matrix(g1,k,lmbda1,existence=True) and
                difference_matrix(g2,k,lmbda2,existence=True)):
                return (g1,lmbda1),(g2,lmbda2)

    return False

def test_Hecke_relations(a,b,C):
    r"""Testing Hecke relations for the Fourier coefficients in C



    INPUT:
    -''C'' -- dictionary of complex (Fourier coefficients)
    -''a'' -- integer
    -''b'' -- integer

    OUTPUT:
    -''diff'' -- real : |C(a)C(b)-C(ab)| if (a,b)=1

    EXAMPLE::

    
    sage: S=MaassWaveForms(Gamma0(1))
    sage: R=mpmath.mpf(9.53369526135355755434423523592877032382125639510725198237579046413534)
    sage: Y=mpmath.mpf(0.85)
    sage: C=coefficients_for_Maass_waveforms(S,R,Y,10,20,12)
    sage: d=test_Hecke_relations(C,2,3); mppr(d)
    '9.29e-8'
    sage: C=coefficients_for_Maass_waveforms(S,R,Y,30,50,20)
    sage: d=test_Hecke_relations(C,2,3); mppr(d)
    '3.83e-43'
    
    
    """
    c=gcd(Integer(a),Integer(b))
    lhs=C[0][a]*C[0][b]
    rhs=mpmath.mpf(0)
    for d in divisors(c):
        rhs=rhs+C[0][Integer(a*b/d/d)]
    return abs(rhs-lhs)

    def _coset_reduction_data_second_coord(G):
        """
        Compute data used for determining the canonical coset
        representative of an element of SL_2(Z) modulo G. This
        function specifically returns data needed for the second part
        of the reduction step (the second coordinate).

        INPUT:
            self

        OUTPUT:
            a dictionary v with keys the divisors of N such that v[d]
            is the subgroup {h in H : h = 1 (mod N/d)}.

        EXAMPLES::

            sage: G = GammaH(240,[7,239])
            sage: G._coset_reduction_data_second_coord()
            {1: [1], 2: [1], 3: [1], 4: [1], 5: [1, 49], 6: [1], 48: [1, 191], 8: [1], 80: [1, 7, 49, 103], 10: [1, 49], 12: [1], 15: [1, 49], 240: [1, 7, 49, 103, 137, 191, 233, 239], 40: [1, 7, 49, 103], 20: [1, 49], 24: [1, 191], 120: [1, 7, 49, 103, 137, 191, 233, 239], 60: [1, 49, 137, 233], 30: [1, 49, 137, 233], 16: [1]}
            sage: G = GammaH(1200,[-1,7]); G
            Congruence Subgroup Gamma_H(1200) with H generated by [7, 1199]
            sage: K = G._coset_reduction_data_second_coord().keys() ; K.sort()
            sage: K == divisors(1200)
            True
        """
        H = G._list_of_elements_in_H()
        N = G.level()
        v = { 1: [1] , N: H }
        for d in [x for x in divisors(N) if x > 1 and x < N ]:
            N_over_d = N // d
            v[d] = [x for x in H if x % N_over_d == 1]
        return v

def number_of_Gamma0_NFCusps(N):
    """
    Returns the total number of orbits of cusps under the action of the
    congruence subgroup `\\Gamma_0(N)`.

    INPUT:

    - ``N`` -- a number field ideal.

    OUTPUT:

    ingeter -- the number of orbits of cusps under Gamma0(N)-action.

    EXAMPLES::

        sage: k.<a> = NumberField(x^3 + 11)
        sage: N = k.ideal(2, a+1)
        sage: from sage.modular.cusps_nf import number_of_Gamma0_NFCusps
        sage: number_of_Gamma0_NFCusps(N)
        4
        sage: L = Gamma0_NFCusps(N)
        sage: len(L) == number_of_Gamma0_NFCusps(N)
        True
    """
    k = N.number_field()
    # The number of Gamma0(N)-sub-orbits for each Gamma-orbit:
    from sage.rings.arith import divisors
    s = sum([len(list((d+N/d).invertible_residues_mod(k.unit_group().gens()))) \
                                                for d in divisors(N)])
    # There are h Gamma-orbits, with h class number of underlying number field.
    return s*k.class_number()

def AllCusps(N):
    r"""
    Return a list of CuspFamily objects corresponding to the cusps of
    `X_0(N)`.

    INPUT:

    -  ``N`` - (integer): the level


    EXAMPLES::

        sage: AllCusps(18)
        [(Inf), (c_{2}), (c_{3,1}), (c_{3,2}), (c_{6,1}), (c_{6,2}), (c_{9}), (0)]
    """
    try:
        N = ZZ(N)
        assert N>0
    except TypeError:
        raise TypeError, "N must be an integer"
    except AssertionError:
        raise AssertionError, "N must be positive"
    c = []
    for d in divisors(N):
        n = num_cusps_of_width(N, d)
        if n == 1:
            c.append(CuspFamily(N, d))
        elif n > 1:
            for i in xrange(n):
                c.append(CuspFamily(N, d, label=str(i+1)))
    return c

    def _b_power_k(self, k):
        r"""
        An expression involving moebius inversion in the powersum generators.

        For a positive value of ``k``, this expression is

        .. MATH::

            \frac{1}{k} \sum_{d|k} \mu(d/k) p_d.

        INPUT:

        - ``k`` -- a positive integer

        OUTPUT:

        - an expression in the powersum basis of the symmetric functions

        EXAMPLES::

            sage: st = SymmetricFunctions(QQ).st()
            sage: st._b_power_k(1)
            p[1]
            sage: st._b_power_k(2)
            -1/2*p[1] + 1/2*p[2]
            sage: st._b_power_k(6)
            1/6*p[1] - 1/6*p[2] - 1/6*p[3] + 1/6*p[6]

        """
        if k == 1:
            return self._p([1])
        if k > 0:
            return ~k * self._p.sum(moebius(k/d)*self._p([d])
                                    for d in divisors(k))

def TD_find_product_decomposition(k,n):
    r"""
    Attempts to find a factorization of `n` in order to build a `TD(k,n)`.

    If Sage can build a `TD(k,n_1)` and a `TD(k,n_2)` such that `n=n_1\times
    n_2` then a `TD(k,n)` can be built (from the function
    :func:`transversal_design`). This method returns such a pair of integers if
    it exists, and ``None`` otherwise.

    INPUT:

    - ``k,n`` (integers) -- see above.

    .. SEEALSO::

        :func:`TD_product` that actually build a product

    EXAMPLES::

        sage: from sage.combinat.designs.orthogonal_arrays import TD_find_product_decomposition
        sage: TD_find_product_decomposition(6, 84)
        (7, 12)

        sage: TD1 = designs.transversal_design(6, 7)
        sage: TD2 = designs.transversal_design(6, 12)
        sage: from sage.combinat.designs.orthogonal_arrays import TD_product
        sage: TD = TD_product(6, TD1, 7, TD2, 12)
    """
    from sage.rings.arith import divisors
    for n1 in divisors(n)[1:-1]: # we ignore 1 and n
        n2 = n//n1
        if transversal_design(k, n1, existence = True) and transversal_design(k, n2, existence = True):
            return n1,n2
    return None

def AllCusps(N):
    r"""
    Return a list of CuspFamily objects corresponding to the cusps of
    `X_0(N)`.

    INPUT:

    -  ``N`` - (integer): the level


    EXAMPLES::

        sage: AllCusps(18)
        [(Inf), (c_{2}), (c_{3,1}), (c_{3,2}), (c_{6,1}), (c_{6,2}), (c_{9}), (0)]
        sage: AllCusps(0)
        Traceback (most recent call last):
        ...
        ValueError: N must be positive
    """
    N = ZZ(N)
    if N <= 0:
        raise ValueError("N must be positive")

    c = []
    for d in divisors(N):
        n = num_cusps_of_width(N, d)
        if n == 1:
            c.append(CuspFamily(N, d))
        elif n > 1:
            for i in xrange(n):
                c.append(CuspFamily(N, d, label=str(i+1)))
    return c

    def _Weyl_law_consts(self):
        r"""
        Compute constants for the Weyl law on self._G

        OUTPUT:

        - tuple of real numbers

        EXAMPLES::


            sage: M=MaassWaveForms(MySubgroup(Gamma0(1)))
            sage: M._Weyl_law_consts()
            (0, 2/pi, (log(pi) - log(2) + 2)/pi, 0, -2)
        """
        import mpmath
        pi=mpmath.fp.pi
        ix=Integer(self._G.index())
        nc=Integer(len(self._G.cusps()))
        if(self._G.is_congruence()):
            lvl=Integer(self._G.level())
        else:
            lvl=0
        n2=Integer(self._G.nu2())
        n3=Integer(self._G.nu3())
        c1=ix/Integer(12)
        c2=Integer(2)*nc/pi
        c3=nc*(Integer(2)-ln(Integer(2))+ln(pi))/pi
        if(lvl<>0):
            A=1
            for q in divisors(lvl):
                num_prim_dc=0
                DG=DirichletGroup(q)
                for chi in DG.list():
                    if(chi.is_primitive()):
                        num_prim_dc=num_prim_dc+1
                for m in divisors(lvl):
                    if(lvl % (m*q) == 0   and m % q ==0 ): 
                        fak=(q*lvl)/gcd(m,lvl/m)
                        A=A*Integer(fak)**num_prim_dc        
            c4=-ln(A)/pi
        else:
            c4=Integer(0)
        # constant term
        c5=-ix/144+n2/8+n3*2/9-nc/4-1
        return (c1,c2,c3,c4,c5)

    def reduce_basis(self, long_etas):
        r"""
        Produce a more manageable basis via LLL-reduction.

        INPUT:


        - ``long_etas`` -  a list of EtaGroupElement objects (which
          should all be of the same level)


        OUTPUT:


        - a new list of EtaGroupElement objects having
          hopefully smaller norm


        ALGORITHM: We define the norm of an eta-product to be the
        `L^2` norm of its divisor (as an element of the free
        `\ZZ`-module with the cusps as basis and the
        standard inner product). Applying LLL-reduction to this gives a
        basis of hopefully more tractable elements. Of course we'd like to
        use the `L^1` norm as this is just twice the degree, which
        is a much more natural invariant, but `L^2` norm is easier
        to work with!

        EXAMPLES::

            sage: EtaGroup(4).reduce_basis([ EtaProduct(4, {1:8,2:24,4:-32}), EtaProduct(4, {1:8, 4:-8})])
            [Eta product of level 4 : (eta_1)^8 (eta_4)^-8,
            Eta product of level 4 : (eta_1)^-8 (eta_2)^24 (eta_4)^-16]
        """
        N = self.level()
        cusps = AllCusps(N)
        r = matrix(ZZ, [[et.order_at_cusp(c) for c in cusps] for et in long_etas])
        V = FreeModule(ZZ, r.ncols())
        A = V.submodule_with_basis([V(rw) for rw in r.rows()])
        rred = r.LLL()
        short_etas = []
        for shortvect in rred.rows():
            bv = A.coordinates(shortvect)
            dict = {}
            for d in divisors(N):
                dict[d] = sum( [bv[i]*long_etas[i].r(d) for i in xrange(r.nrows())])
            short_etas.append(self(dict))
        return short_etas

def p1NFlist(N):
    """
    Returns a list of the normalized elements of `\\mathbb{P}^1(R/N)`, where
    `N` is an integral ideal.

    INPUT:

    -  ``N`` - integral ideal (the level or modulus).

    EXAMPLES::

        sage: k.<a> = NumberField(x^2 + 23)
        sage: N = k.ideal(3)
        sage: from sage.modular.modsym.p1list_nf import p1NFlist, psi
        sage: len(p1NFlist(N))==psi(N)
        True
    """
    k = N.number_field()

    L = [MSymbol(N, k(0),k(1), check=False)]
    #N.residues() = iterator through the residues mod N
    L = L+[MSymbol(N, k(1), r, check=False) for r in N.residues()]

    from sage.rings.arith import divisors
    for D in divisors(N):
        if not D.is_trivial() and D!=N:
            #we find Dp ideal coprime to N, in inverse class to D
            if D.is_principal():
                Dp = k.ideal(1)
                c = D.gens_reduced()[0]
            else:
                it = k.primes_of_degree_one_iter()
                Dp = it.next()
                while not Dp.is_coprime(N) or not (Dp*D).is_principal():
                    Dp = it.next()
                c = (D*Dp).gens_reduced()[0]
            #now we find all the (c,d)'s which have associated divisor D
            I = D + N/D
            for d in (N/D).residues():
                if I.is_coprime(d):
                    M = D.prime_to_idealM_part(N/D)
                    u = (Dp*M).element_1_mod(N/D)
                    d1 = u*d + (1-u)
                    L.append(MSymbol(N, c, d1, check=False).normalize())
    return L

    def cardinality(self):
        r"""
        Return the number of integer necklaces with the evaluation ``content``.

        The formula for the number of necklaces of content `\alpha`
        a composition of `n` is:

        .. MATH::

            \sum_{d|gcd(\alpha)} \phi(d)
            \binom{n/d}{\alpha_1/d, \ldots, \alpha_\ell/d},

        where `\phi(d)` is the Euler `\phi` function.

        EXAMPLES::

            sage: Necklaces([]).cardinality()
            0
            sage: Necklaces([2,2]).cardinality()
            2
            sage: Necklaces([2,3,2]).cardinality()
            30
            sage: Necklaces([0,3,2]).cardinality()
            2

        Check to make sure that the count matches up with the number of
        necklace words generated.

        ::

            sage: comps = [[],[2,2],[3,2,7],[4,2],[0,4,2],[2,0,4]]+Compositions(4).list()
            sage: ns = [ Necklaces(comp) for comp in comps]
            sage: all( [ n.cardinality() == len(n.list()) for n in ns] )
            True
        """
        evaluation = self._content
        le = list(evaluation)
        if not le:
            return 0

        n = sum(le)

        return sum(euler_phi(j)*factorial(n/j) / prod(factorial(ni/j)
                    for ni in evaluation) for j in divisors(gcd(le))) / n

def summand(part, n):
    """
    Create the summand used in the Harrison count for a given partition.

    Args:
        part (tuple): A partition of `n` represented as a tuple.
        n (int): The integer for which `part` is a partition.

    Returns:
        int: The summand corresponding to the partition `part` of `n`.
    """

    t = 1
    count = list(cycle_count(part, n)) + (factorial(n)-n)*[0]
    for i in range(1,n+1):
        for j in range(1,n+1):
            s = sum([d*(count[d-1]) for d in divisors(lcm(i,j))])
            t = t*(s**(count[i-1]*count[j-1]*gcd(i,j)))
    t = t*factorial(n)/(prod(factorial(count[d-1])*(d**(count[d-1])) for d in range(1,n+1)))
    return t

    def _find_cusps(self):
        r"""
        Return an ordered list of inequivalent cusps for self, i.e. a
        set of representatives for the orbits of self on
        `\mathbb{P}^1(\QQ)`.  These are returned in a reduced
        form; see self.reduce_cusp for the definition of reduced.

        ALGORITHM:
            Uses explicit formulae specific to `\Gamma_0(N)`: a reduced cusp on
            `\Gamma_0(N)` is always of the form `a/d` where `d | N`, and `a_1/d
            \sim a_2/d` if and only if `a_1 \cong a_2 \bmod {\rm gcd}(d,
            N/d)`.

        EXAMPLES::

            sage: Gamma0(90)._find_cusps()
            [0, 1/45, 1/30, 1/18, 1/15, 1/10, 1/9, 2/15, 1/6, 1/5, 1/3, 11/30, 1/2, 2/3, 5/6, Infinity]
            sage: Gamma0(1).cusps()
            [Infinity]
            sage: Gamma0(180).cusps() == Gamma0(180).cusps(algorithm='modsym')
            True
        """
        N = self.level()
        s = []

        for d in arith.divisors(N):
            w = arith.gcd(d, N//d)
            if w == 1:
                if d == 1:
                    s.append(Cusp(1,0))
                elif d == N:
                    s.append(Cusp(0,1))
                else:
                    s.append(Cusp(1,d))
            else:
                for a in xrange(1, w):
                    if arith.gcd(a, w) == 1:
                        while arith.gcd(a, d//w) != 1:
                            a += w
                        s.append(Cusp(a,d))
        return sorted(s)