Python basic.C类代码示例

    def eval(cls, arg):
        arg = sympify(arg)

        if arg.is_Number:
            if arg is S.NaN:
                return S.NaN
            elif arg is S.Infinity:
                return S.Infinity
            elif arg is S.NegativeInfinity:
                return S.NegativeInfinity
            elif arg is S.Zero:
                return S.Zero
            elif arg is S.One:
                return C.log(2 ** S.Half + 1)
            elif arg is S.NegativeOne:
                return C.log(2 ** S.Half - 1)
            elif arg.is_negative:
                return -cls(-arg)
        else:
            i_coeff = arg.as_coefficient(S.ImaginaryUnit)

            if i_coeff is not None:
                return S.ImaginaryUnit * C.asin(i_coeff)
            else:
                coeff, terms = arg.as_coeff_terms()

                if coeff.is_negative:
                    return -cls(-arg)

def monomial_count(V, N):
    r"""
    Computes the number of monomials.

    The number of monomials is given by the following formula:

    .. math::

        \frac{(\#V + N)!}{\#V! N!}

    where `N` is a total degree and `V` is a set of variables.

    **Examples**

    >>> from sympy import monomials, monomial_count
    >>> from sympy.abc import x, y

    >>> monomial_count(2, 2)
    6

    >>> M = monomials([x, y], 2)

    >>> sorted(M)
    [1, x, y, x**2, y**2, x*y]
    >>> len(M)
    6

    """
    return C.factorial(V + N) / C.factorial(V) / C.factorial(N)

    def _eval_expand_trig(self, **hints):
        from sympy import expand_mul
        arg = self.args[0]
        x = None
        if arg.is_Add:  # TODO, implement more if deep stuff here
            # TODO: Do this more efficiently for more than two terms
            x, y = arg.as_two_terms()
            sx = sin(x, evaluate=False)._eval_expand_trig()
            sy = sin(y, evaluate=False)._eval_expand_trig()
            cx = cos(x, evaluate=False)._eval_expand_trig()
            cy = cos(y, evaluate=False)._eval_expand_trig()
            return sx*cy + sy*cx
        else:
            n, x = arg.as_coeff_Mul(rational=True)
            if n.is_Integer:  # n will be positive because of .eval
                # canonicalization

                # See http://mathworld.wolfram.com/Multiple-AngleFormulas.html
                if n.is_odd:
                    return (-1)**((n - 1)/2)*C.chebyshevt(n, sin(x))
                else:
                    return expand_mul((-1)**(n/2 - 1)*cos(x)*C.chebyshevu(n -
                        1, sin(x)), deep=False)
            pi_coeff = _pi_coeff(arg)
            if pi_coeff is not None:
                if pi_coeff.is_Rational:
                    return self.rewrite(sqrt)
        return sin(arg)

 def vertices(self):
     points = []
     c, r, n = self
     v = 2*S.Pi/n
     for k in xrange(0, n):
         points.append( Point(c[0] + r*C.cos(k*v), c[1] + r*C.sin(k*v)) )
     return points

 def _eval_expand_complex(self, *args):
     if self.args[0].is_real:
         return self
     re, im = self.args[0].as_real_imag()
     denom = sin(re)**2 + C.sinh(im)**2
     return (sin(re)*cos(re) - \
         S.ImaginaryUnit*C.sinh(im)*C.cosh(im))/denom

def solve_ODE_first_order(eq, f):
    """
    solves many kinds of first order odes, different methods are used
    depending on the form of the given equation. Now the linear
    and Bernoulli cases are implemented.
    """
    from sympy.integrals.integrals import integrate
    x = f.args[0]
    f = f.func

    #linear case: a(x)*f'(x)+b(x)*f(x)+c(x) = 0
    a = Wild('a', exclude=[f(x)])
    b = Wild('b', exclude=[f(x)])
    c = Wild('c', exclude=[f(x)])

    r = eq.match(a*diff(f(x),x) + b*f(x) + c)
    if r:
        t = C.exp(integrate(r[b]/r[a], x))
        tt = integrate(t*(-r[c]/r[a]), x)
        return (tt + Symbol("C1"))/t

    #Bernoulli case: a(x)*f'(x)+b(x)*f(x)+c(x)*f(x)^n = 0
    n = Wild('n', exclude=[f(x)])

    r = eq.match(a*diff(f(x),x) + b*f(x) + c*f(x)**n)
    if r:
        t = C.exp((1-r[n])*integrate(r[b]/r[a],x))
        tt = (r[n]-1)*integrate(t*r[c]/r[a],x)
        return ((tt + Symbol("C1"))/t)**(1/(1-r[n]))

    #other cases of first order odes will be implemented here

    raise NotImplementedError("solve_ODE_first_order: Cannot solve " + str(eq))

 def _eval_expand_complex(self, deep=True, **hints):
     re, im = self.args[0].as_real_imag()
     if deep:
         re = re.expand(deep, **hints)
         im = im.expand(deep, **hints)
     cos, sin = C.cos(im), C.sin(im)
     return exp(re) * cos + S.ImaginaryUnit * exp(re) * sin

    def _eval_expand_trig(self, **hints):
        arg = self.args[0]
        x = None
        if arg.is_Add:
            from sympy import symmetric_poly

            n = len(arg.args)
            CX = []
            for x in arg.args:
                cx = cot(x, evaluate=False)._eval_expand_trig()
                CX.append(cx)

            Yg = numbered_symbols("Y")
            Y = [Yg.next() for i in xrange(n)]

            p = [0, 0]
            for i in xrange(n, -1, -1):
                p[(n - i) % 2] += symmetric_poly(i, Y) * (-1) ** (((n - i) % 4) // 2)
            return (p[0] / p[1]).subs(zip(Y, CX))
        else:
            coeff, terms = arg.as_coeff_Mul(rational=True)
            if coeff.is_Integer and coeff > 1:
                I = S.ImaginaryUnit
                z = C.Symbol("dummy", real=True)
                P = ((z + I) ** coeff).expand()
                return (C.re(P) / C.im(P)).subs([(z, cot(terms))])
        return cot(arg)

 def _eval_rewrite_as_polynomial(self, n, m, x):
     k = C.Dummy("k")
     kern = (
         C.factorial(2 * n - 2 * k)
         / (2 ** n * C.factorial(n - k) * C.factorial(k) * C.factorial(n - 2 * k - m))
         * (-1) ** k
         * x ** (n - m - 2 * k)
     )
     return (1 - x ** 2) ** (m / 2) * C.Sum(kern, (k, 0, C.floor((n - m) * S.Half)))

 def _eval_expand_func(self, **hints):
     n, m, theta, phi = self.args
     rv = (
         sqrt((2 * n + 1) / (4 * pi) * C.factorial(n - m) / C.factorial(n + m))
         * C.exp(I * m * phi)
         * assoc_legendre(n, m, C.cos(theta))
     )
     # We can do this because of the range of theta
     return rv.subs(sqrt(-cos(theta) ** 2 + 1), sin(theta))

def solve_ODE_second_order(eq, f):
    """
    solves many kinds of second order odes, different methods are used
    depending on the form of the given equation. So far the constants
    coefficients case and a special case are implemented.
    """
    x = f.args[0]
    f = f.func

    #constant coefficients case: af''(x)+bf'(x)+cf(x)=0
    a = Wild('a', exclude=[x])
    b = Wild('b', exclude=[x])
    c = Wild('c', exclude=[x])

    r = eq.match(a*f(x).diff(x,x) + c*f(x))
    if r:
        return Symbol("C1")*C.sin(sqrt(r[c]/r[a])*x)+Symbol("C2")*C.cos(sqrt(r[c]/r[a])*x)

    r = eq.match(a*f(x).diff(x,x) + b*diff(f(x),x) + c*f(x))
    if r:
        r1 = solve(r[a]*x**2 + r[b]*x + r[c], x)
        if r1[0].is_real:
            if len(r1) == 1:
                return (Symbol("C1") + Symbol("C2")*x)*exp(r1[0]*x)
            else:
                return Symbol("C1")*exp(r1[0]*x) + Symbol("C2")*exp(r1[1]*x)
        else:
            r2 = abs((r1[0] - r1[1])/(2*S.ImaginaryUnit))
            return (Symbol("C2")*C.cos(r2*x) + Symbol("C1")*C.sin(r2*x))*exp((r1[0] + r1[1])*x/2)

    #other cases of the second order odes will be implemented here

    #special equations, that we know how to solve
    a = Wild('a')
    t = x*exp(f(x))
    tt = a*t.diff(x, x)/t
    r = eq.match(tt.expand())
    if r:
        return -solve_ODE_1(f(x), x)

    t = x*exp(-f(x))
    tt = a*t.diff(x, x)/t
    r = eq.match(tt.expand())
    if r:
        #check, that we've rewritten the equation correctly:
        #assert ( r[a]*t.diff(x,2)/t ) == eq.subs(f, t)
        return solve_ODE_1(f(x), x)

    neq = eq*exp(f(x))/exp(-f(x))
    r = neq.match(tt.expand())
    if r:
        #check, that we've rewritten the equation correctly:
        #assert ( t.diff(x,2)*r[a]/t ).expand() == eq
        return solve_ODE_1(f(x), x)

    raise NotImplementedError("solve_ODE_second_order: cannot solve " + str(eq))

    def eval(cls, r, k):
        r, k = map(sympify, (r, k))

        if k.is_Number:
            if k is S.Zero:
                return S.One
            elif k.is_Integer:
                if k.is_negative:
                    return S.Zero
                else:
                    if r.is_Integer and r.is_nonnegative:
                        r, k = int(r), int(k)

                        if k > r:
                            return S.Zero
                        elif k > r // 2:
                            k = r - k

                        M, result = int(sqrt(r)), 1

                        for prime in sieve.primerange(2, r+1):
                            if prime > r - k:
                                result *= prime
                            elif prime > r // 2:
                                continue
                            elif prime > M:
                                if r % prime < k % prime:
                                    result *= prime
                            else:
                                R, K = r, k
                                exp = a = 0

                                while R > 0:
                                    a = int((R % prime) < (K % prime + a))
                                    R, K = R // prime, K // prime
                                    exp = a + exp

                                if exp > 0:
                                    result *= prime**exp

                        return C.Integer(result)
                    else:
                        result = r - k + 1

                        for i in xrange(2, k+1):
                            result *= r-k+i
                            result /= i

                        return result

        if k.is_integer:
            if k.is_negative:
                return S.Zero
        else:
            return C.gamma(r+1)/(C.gamma(r-k+1)*C.gamma(k+1))

 def _eval_expand_complex(self, deep=True, **hints):
     if self.args[0].is_real:
         if deep:
             hints['complex'] = False
             return self.expand(deep, **hints)
         else:
             return self
     if deep:
         re, im = self.args[0].expand(deep, **hints).as_real_imag()
     else:
         re, im = self.args[0].as_real_imag()
     return sin(re)*C.cosh(im) + S.ImaginaryUnit*cos(re)*C.sinh(im)

    def taylor_term(n, x, *previous_terms):
        if n == 0:
            return 1 / sympify(x)
        elif n < 0 or n % 2 == 0:
            return S.Zero
        else:
            x = sympify(x)

            B = C.bernoulli(n+1)
            F = C.factorial(n+1)

            return (-1)**((n+1)//2) * 2**(n+1) * B/F * x**n

 def as_real_imag(self, deep=True, **hints):
     if self.args[0].is_real:
         if deep:
             return (self.expand(deep, **hints), S.Zero)
         else:
             return (self, S.Zero)
     if deep:
         re, im = self.args[0].expand(deep, **hints).as_real_imag()
     else:
         re, im = self.args[0].as_real_imag()
     denom = sinh(re) ** 2 + C.sin(im) ** 2
     return (sinh(re) * cosh(re) / denom, -C.sin(im) * C.cos(im) / denom)