Python sympy.diff函数代码示例

    def __init__(self, nddistr=None, d=None):
        super(TwoVarsModel, self).__init__(nddistr, [d])
        self.eliminate_other([d])
        self.d = d
        self.vars = []
        self.symvars = []
        for var in nddistr.Vars: #self.free_rvs:
            self.vars.append(var)
            self.symvars.append(var.getSymname())
        #print "=====", self.vars
        #print self.symvars
        #print self.dep_rvs
        #print self.rv_to_equation
        self.symop = self.rv_to_equation[d]

        if len(self.vars) != 2:
            raise Exception("use it with two variables")
        x = self.symvars[0]
        y = self.symvars[1]
        z = sympy.Symbol("z")
        self.fun_alongx = eq_solve(self.symop, z, y)[0]
        self.fun_alongy = eq_solve(self.symop, z, x)[0]

        self.lfun_alongx = my_lambdify([x, z], self.fun_alongx, "numpy")
        self.lfun_alongy = my_lambdify([y, z], self.fun_alongy, "numpy")
        self.Jx = 1 * sympy.diff(self.fun_alongx, z)
        #print "Jx=", self.Jx
        #print "fun_alongx=", self.fun_alongx
        self.Jy = 1 * sympy.diff(self.fun_alongy, z)
        self.lJx = my_lambdify([x, z], self.Jx, "numpy")
        self.lJy = my_lambdify([y, z], self.Jy, "numpy")
        self.z = z

def test_atan2():
    assert atan2(0, 0) == S.NaN
    assert atan2(0, 1) == 0
    assert atan2(1, 1) == pi/4
    assert atan2(1, 0) == pi/2
    assert atan2(1, -1) == 3*pi/4
    assert atan2(0, -1) == pi
    assert atan2(-1, -1) == -3*pi/4
    assert atan2(-1, 0) == -pi/2
    assert atan2(-1, 1) == -pi/4

    u = Symbol("u", positive=True)
    assert atan2(0, u) == 0
    u = Symbol("u", negative=True)
    assert atan2(0, u) == pi

    assert atan2(y, oo) ==  0
    assert atan2(y, -oo)==  2*pi*Heaviside(re(y)) - pi

    assert atan2(y, x).rewrite(log) == -I*log((x + I*y)/sqrt(x**2 + y**2))
    assert atan2(y, x).rewrite(atan) == 2*atan(y/(x + sqrt(x**2 + y**2)))

    assert diff(atan2(y, x), x) == -y/(x**2 + y**2)
    assert diff(atan2(y, x), y) == x/(x**2 + y**2)

    assert simplify(diff(atan2(y, x).rewrite(log), x)) == -y/(x**2 + y**2)
    assert simplify(diff(atan2(y, x).rewrite(log), y)) ==  x/(x**2 + y**2)

    assert isinstance(atan2(2, 3*I).n(), atan2)

def Simple_manifold_with_scalar_function_derivative():
    Print_Function()
    coords = (x,y,z) = symbols('x y z')
    basis = (e1, e2, e3, grad) = MV.setup('e_1 e_2 e_3',metric='[1,1,1]',coords=coords)
    # Define surface
    mfvar = (u,v) = symbols('u v')
    X = u*e1+v*e2+(u**2+v**2)*e3
    print '\\f{X}{u,v} =',X
    MF = Manifold(X,mfvar)
    (eu,ev) = MF.Basis()
    # Define field on the surface.
    g = (v+1)*log(u)

    print '\\f{g}{u,v} =',g

    # Method 1: Using old Manifold routines.
    VectorDerivative = (MF.rbasis[0]/MF.E_sq)*diff(g,u) + (MF.rbasis[1]/MF.E_sq)*diff(g,v)
    print '\\eval{\\nabla g}{u=1,v=0} =', VectorDerivative.subs({u:1,v:0})

    # Method 2: Using new Manifold routines.
    dg = MF.Grad(g)
    print '\\eval{\\f{Grad}{g}}{u=1,v=0} =', dg.subs({u:1,v:0})
    dg = MF.grad*g
    print '\\eval{\\nabla g}{u=1,v=0} =', dg.subs({u:1,v:0})
    return

    def __init__(self):
        self.num_lines, self.num_marks = 5, 2
        self._qp = {}

        function_type = random.choice(['sqrt', 'quadratic', 'product'])

        if function_type == 'sqrt':
            outer_function = all_functions.request_linear(difficulty=2).equation
            inner_function = all_functions.request_linear(difficulty=1).equation
            inner_function = inner_function.replace(lambda expr: expr.is_Symbol, lambda expr: sympy.sqrt(expr))

            self._qp['equation'] = outer_function.replace(x, inner_function)
            self._qp['derivative'] = sympy.diff(self._qp['equation'])

        elif function_type == 'quadratic':
            power_two_coeff = not_named_yet.randint_no_zero(-3, 3)
            power_one_coeff = not_named_yet.randint_no_zero(-5, 5)
            inner_function = power_two_coeff * x ** 2 + power_one_coeff * x
            index = random.randint(3, 5)

            self._qp['equation'] = inner_function ** index
            self._qp['derivative'] = sympy.diff(self._qp['equation'])

        elif function_type == 'product':
            left_function = x ** random.randint(1, 3)
            right_outer_function = random.choice([sympy.sin, sympy.cos, sympy.log, sympy.exp])
            right_inner_function = not_named_yet.randint_no_zero(-3, 3) * x

            self._qp['equation'] = left_function * right_outer_function(right_inner_function)
            self._qp['derivative'] = sympy.diff(self._qp['equation'])

def test_derivative_evaluate():
    assert Derivative(sin(x), x) != diff(sin(x), x)
    assert Derivative(sin(x), x).doit() == diff(sin(x), x)

    assert Derivative(Derivative(f(x), x), x) == diff(f(x), x, x)
    assert Derivative(sin(x), x, 0) == sin(x)
    assert Derivative(sin(x), (x, y), (x, -y)) == sin(x)

 def grad(self, func):
     """
     Calculate the gradient of 'func'.
     """
     return Matrix([diff(func, self.xs[0]) / self.h[0], \
                    diff(func, self.xs[1]) / self.h[1], \
                    diff(func, self.xs[2]) / self.h[2]])

def test_gegenbauer():
    n = Symbol("n")
    a = Symbol("a")

    assert gegenbauer(0, a, x) == 1
    assert gegenbauer(1, a, x) == 2*a*x
    assert gegenbauer(2, a, x) == -a + x**2*(2*a**2 + 2*a)
    assert gegenbauer(3, a, x) == \
        x**3*(4*a**3/3 + 4*a**2 + 8*a/3) + x*(-2*a**2 - 2*a)

    assert gegenbauer(-1, a, x) == 0
    assert gegenbauer(n, S(1)/2, x) == legendre(n, x)
    assert gegenbauer(n, 1, x) == chebyshevu(n, x)
    assert gegenbauer(n, -1, x) == 0

    X = gegenbauer(n, a, x)
    assert isinstance(X, gegenbauer)

    assert gegenbauer(n, a, -x) == (-1)**n*gegenbauer(n, a, x)
    assert gegenbauer(n, a, 0) == 2**n*sqrt(pi) * \
        gamma(a + n/2)/(gamma(a)*gamma(-n/2 + S(1)/2)*gamma(n + 1))
    assert gegenbauer(n, a, 1) == gamma(2*a + n)/(gamma(2*a)*gamma(n + 1))

    assert gegenbauer(n, Rational(3, 4), -1) == zoo

    m = Symbol("m", positive=True)
    assert gegenbauer(m, a, oo) == oo*RisingFactorial(a, m)

    assert conjugate(gegenbauer(n, a, x)) == gegenbauer(n, conjugate(a), conjugate(x))

    assert diff(gegenbauer(n, a, x), n) == Derivative(gegenbauer(n, a, x), n)
    assert diff(gegenbauer(n, a, x), x) == 2*a*gegenbauer(n - 1, a + 1, x)

def Delp2(f):
    """ Laplacian in X-Z
    """
    d2fdx2 = diff(f, x, 2)
    d2fdz2 = diff(f, z, 2)

    return d2fdx2 + d2fdz2

def ItoFormula_2(function,a,b,dt):
	func = compile(function,"",'eval')
	a=float(a)
	b=float(b)
	n=int((b-a)/dt)	
	WP = 0
	t = a
	for i in range(int(a/dt)):
		WP = WP+np.random.normal(0,math.sqrt(dt))
	part_1 = eval(func)
	f = sym.sympify(function)
	WP = sym.Symbol('WP')
	df = str(sym.diff(f,WP))
	ddf = str(sym.diff(df,WP))
	part_2 = 0
	part_3 = 0
	result_2 = 0
	WP = 0
	for i in range(0,n):
		t=dt*((a/dt)+i)
		WP = WP+np.random.normal(0,math.sqrt(dt))
		k2 = eval(df,valueDic)
		k3 = eval(ddf,valueDic)
		part_2 = part_2 + k2*(WT[i+1]-WT[i])
		part_3 = part_3 + k3*dt
	result_2 = part_1 + part_2 + part_3
	print 'ito formula result: '+str(result_2)

def eleq(Lagrangian, Friction = 0, t = Symbol('t')):
    """
    Returns Euler-Lagrange equations of the lagrangian system.

    Examples
    ========

    >>> from sympy import *
    >>> t, k = symbols('t k')
    >>> x = symbols('x', cls=Function)
    >>> eleq(diff(x(t),t)**2/2 - k*x(t)**2/2)
    {a_x: -k*x}

    >>> a = symbols('a')
    >>> eleq(diff(x(t),t)**2/2 - k*x(t)**2/2, a*diff(x(t),t)**2/2)
    {a_x: -*a*v_x - k*x}
    """
    Lagrangian = simplify(Lagrangian)
    var_list = [list(x.atoms(Function))[0] for x in Lagrangian.atoms(Derivative)]
    nvar = len(var_list)
    ecu_list = [ diff(Lagrangian, variable) - diff(Lagrangian, diff(variable,t), t) -  diff(Friction, diff(variable,t)) for variable in var_list ]
    str_list = [ str(variable).replace("("+str(t)+")","") for variable in var_list ]
    a_subs = {diff(var_list[i],t,2): Symbol('a_' + str_list[i]) for i in range(nvar)}
    v_subs = {diff(var_list[i],t): Symbol('v_' + str_list[i]) for i in range(nvar)}
    x_subs = {var_list[i]: Symbol(str_list[i]) for i in range(nvar)}
    for i in range(nvar):
        if hasattr(ecu_list[i], "subs"):
            ecu_list[i] = ecu_list[i].subs(a_subs).subs(v_subs).subs(x_subs)
    a_list = sorted(list(a_subs.values()), key = str)
    return solveswc(ecu_list, a_list)

    def _derive_dw_dpsi():
        """Derive the slope dw_dpsi with sympy
        
        Returns
        -------
        dw_dpsi1 : sympy expression
            Derivative including correction factor.
        dw_dpsi2 : sympy expression
            Derivative with correction factor equal to one.
         """
        psi, a, n, m, psir, ws = sympy.symbols('psi, a, n, m, psir, ws')
        psi, ws, a, b, wr, s1, psir = sympy.symbols('psi, ws, a, b, wr, s1, psir')
        from sympy import log, exp        
        C1 = (1 - log(1 + psi / psir) / 
                          log(1.0 + 1.0e6 / psir))

        C2 = 1      
        
        l10 = sympy.symbols('l10')
        
        w = C1 * ((ws - s1 * log(psi)/l10 - wr) * a / (psi**b + a) + wr)
        dw_dpsi1 = sympy.diff(w, psi)                    
        
        w = C2 * ((ws - s1 * log(psi)/l10 - wr) * a / (psi**b + a) + wr)
        dw_dpsi2 = sympy.diff(w, psi)
        return  dw_dpsi1, dw_dpsi2

def hamiltonian(Lagrangian, t = Symbol('t'), delta = False):
    """
    Returns the Hamiltonian of the Lagrangian.

    Examples
    ========

    >>> from sympy import *
    >>> t, k = symbols('t k')
    >>> x = symbols('x', cls=Function)
    >>> hamiltonian(diff(x(t),t)**2/2 - k*x(t)**2/2)
    k*x**2/2 + v_x**2/2
    """
    Lagrangian = simplify(Lagrangian)
    var_list = [list(x.atoms(Function))[0] for x in Lagrangian.atoms(Derivative)]
    nvar = len(var_list)
    # New variables.
    str_list = [ str(variable).replace("("+str(t)+")","") for variable in var_list ]
    v_subs = {diff(var_list[i],t): Symbol('v_' + str_list[i]) for i in range(nvar)}
    x_subs = {var_list[i]: Symbol(str_list[i]) for i in range(nvar)}
    # Hamiltonian calculus.
    dxdLv = 0
    for variable in var_list:
        dxdLv += diff(variable,t)*diff(Lagrangian, diff(variable,t))
    result = simplify((dxdLv - Lagrangian).subs(v_subs).subs(x_subs))
    if delta:
        v0_subs = {Symbol('v_' + str_list[i]): Symbol('v_' + str_list[i] + "0") for i in range(nvar)}
        x0_subs = {Symbol(str_list[i]): Symbol(str_list[i] + "0") for i in range(nvar)}
        return result - result.subs(v0_subs).subs(x0_subs)
    else:
        return result

def test_legendre():
    raises(ValueError, lambda: legendre(-1, x))
    assert legendre(0, x) == 1
    assert legendre(1, x) == x
    assert legendre(2, x) == ((3*x**2 - 1)/2).expand()
    assert legendre(3, x) == ((5*x**3 - 3*x)/2).expand()
    assert legendre(4, x) == ((35*x**4 - 30*x**2 + 3)/8).expand()
    assert legendre(5, x) == ((63*x**5 - 70*x**3 + 15*x)/8).expand()
    assert legendre(6, x) == ((231*x**6 - 315*x**4 + 105*x**2 - 5)/16).expand()

    assert legendre(10, -1) == 1
    assert legendre(11, -1) == -1
    assert legendre(10, 1) == 1
    assert legendre(11, 1) == 1
    assert legendre(10, 0) != 0
    assert legendre(11, 0) == 0

    assert roots(legendre(4, x), x) == {
        sqrt(Rational(3, 7) - Rational(2, 35)*sqrt(30)): 1,
        -sqrt(Rational(3, 7) - Rational(2, 35)*sqrt(30)): 1,
        sqrt(Rational(3, 7) + Rational(2, 35)*sqrt(30)): 1,
        -sqrt(Rational(3, 7) + Rational(2, 35)*sqrt(30)): 1,
    }

    n = Symbol("n")

    X = legendre(n, x)
    assert isinstance(X, legendre)

    assert legendre(-n, x) == legendre(n - 1, x)
    assert legendre(n, -x) == (-1)**n*legendre(n, x)
    assert diff(legendre(n, x), x) == \
        n*(x*legendre(n, x) - legendre(n - 1, x))/(x**2 - 1)
    assert diff(legendre(n, x), n) == Derivative(legendre(n, x), n)

def make_qdiff_to_q012(x_names):
    qdiff_to_q012 = {}
    qdiff_0 = {}
    qdiff_1 = {}
    qdiff_2 = {}

    for date_str in ['tp1', 't', 'tm1']:
        dated_names = [name + date_str for name in x_names]
        dated_q_fun_sym = [sympy.Function(x)(q) for x in dated_names]

        dated_names_0 = [name + '_0' for name in dated_names]
        dated_names_0_sym = [sympy.Symbol(x) for x in dated_names_0]
        qdiff_0.update(dict(zip(dated_q_fun_sym, dated_names_0_sym)))

        dated_qdiffs_1 = [sympy.diff(x, q, 1) for x in dated_q_fun_sym]
        dated_names_1 = [name + '_1' for name in dated_names]
        dated_names_1_sym = [sympy.Symbol(x) for x in dated_names_1]
        qdiff_1.update(dict(zip(dated_qdiffs_1, dated_names_1_sym)))

        dated_qdiffs_2 = [sympy.diff(x, q, 2) for x in dated_q_fun_sym]
        dated_names_2 = [name + '_2' for name in dated_names]
        dated_names_2_sym = [sympy.Symbol(x) for x in dated_names_2]
        qdiff_2.update(dict(zip(dated_qdiffs_2, dated_names_2_sym)))

    qdiff_to_q012.update(qdiff_0)
    qdiff_to_q012.update(qdiff_1)
    qdiff_to_q012.update(qdiff_2)

    return qdiff_to_q012

    def __call__(self, scalar_field):
        """
        Represents the gradient of the given scalar field.

        Parameters
        ==========

        scalar_field : SymPy expression
            The scalar field to calculate the gradient of.

        Examples
        ========

        >>> from sympy.vector import CoordSysCartesian
        >>> C = CoordSysCartesian('C')
        >>> C.delop(C.x*C.y*C.z)
        C.y*C.z*C.i + C.x*C.z*C.j + C.x*C.y*C.k

        """

        scalar_field = express(scalar_field, self.system,
                               variables = True)
        vx = diff(scalar_field, self._x)
        vy = diff(scalar_field, self._y)
        vz = diff(scalar_field, self._z)

        return vx*self._i + vy*self._j + vz*self._k