SymPy-1-13-中文文档-十四-
SymPy 1.13 中文文档(十四)
矢量积分的应用
原文:
docs.sympy.org/latest/modules/vector/vector_integration.html
要在区域上积分一个标量或矢量场,我们必须首先定义一个区域。SymPy 提供了三种定义区域的方法:
-
使用带有
ParametricRegion
的参数方程。 -
使用带有
ImplicitRegion
的隐式方程。 -
使用几何模块的对象。
vector_integrate()
函数用于在任何类型的区域上积分标量或矢量场。它根据对象的性质自动确定积分的类型(线、面或体)。
我们定义一个坐标系并为示例做必要的导入。
>>> from sympy import sin, cos, exp, pi, symbols
>>> from sympy.vector import CoordSys3D, ParametricRegion, ImplicitRegion, vector_integrate
>>> from sympy.abc import r, x, y, z, theta, phi
>>> C = CoordSys3D('C')
周长、表面积和体积的计算
要计算圆的周长,我们需要定义它。让我们使用其参数方程定义它。
>>> param_circle = ParametricRegion((4*cos(theta), 4*sin(theta)), (theta, 0, 2*pi))
我们也可以使用其隐式方程定义一个圆。
>>> implicit_circle = ImplicitRegion((x, y), x**2 + y**2 - 4)
图形的周长等于它在单位标量场上的积分的绝对值。
>>> vector_integrate(1, param_circle)
8*pi
>>> vector_integrate(1, implicit_circle)
4*pi
假设用户想要计算三角形的周长。确定三角形的参数表示可能很困难。相反,用户可以使用几何模块中 Polygon
类的对象。
>>> from sympy.geometry import Point, Polygon
>>> triangle = Polygon(Point(1, 2), (3, 5), (1,6))
>>> vector_integrate(1, triangle)
sqrt(5) + sqrt(13) + 4
要定义一个实心球,我们需要使用三个参数(r,theta 和 phi)。对于ParametricRegion
对象来说,限制的顺序决定积分的符号。
>>> solidsphere = ParametricRegion((r*sin(phi)*cos(theta),r*sin(phi)*sin(theta), r*cos(phi)),
... (phi, 0, pi), (theta, 0, 2*pi), (r, 0, 3))
>>> vector_integrate(1, solidsphere)
36*pi
物体质量的计算
考虑一个三角形片段 𝑅,其顶点为 (0,0),(0, 5),(5,0),密度为 (\rho(x, y) = xy:kg/m²)。找到总质量。
>>> triangle = ParametricRegion((x, y), (x, 0, 5), (y, 0, 5 - x))
>>> vector_integrate(C.x*C.y, triangle)
625/24
找到以 z 轴为中心的圆柱体的质量,其高度为 h,半径为 a,密度为 (\rho = x² + y²:kg/m²)。
>>> a, h = symbols('a h', positive=True)
>>> cylinder = ParametricRegion((r*cos(theta), r*sin(theta), z),
... (theta, 0, 2*pi), (z, 0, h), (r, 0, a))
>>> vector_integrate(C.x**2 + C.y**2, cylinder)
pi*a**4*h/2
通量的计算
1. 考虑空间中存在一个恒定向量场 (E(x, y, z) = a\mathbf{\hat{k}})。半径为 r 的半球位于 x-y 平面上。该场通过球体的通量是多少?
>>> semisphere = ParametricRegion((r*sin(phi)*cos(theta), r*sin(phi)*sin(theta), r*cos(phi)),\
... (phi, 0, pi/2), (theta, 0, 2*pi))
>>> flux = vector_integrate(a*C.k, semisphere)
>>> flux
pi*a*r**2
2. 考虑空间中存在向量场 (E(x, y, z) = x² \mathbf{\hat{k}}) 在 x-y 平面上方,并且在 x-y 平面下方存在一个场 (E(x, y, z) = y² \mathbf{\hat{k}})。该向量场穿过边长为 L 且中心在原点的立方体的通量是多少?
该场是沿着 z 轴平行的,因此只有箱子的顶部和底部会对通量有贡献。
>>> L = symbols('L', positive=True)
>>> top_face = ParametricRegion((x, y, L/2), (x, -L/2, L/2), (y, -L/2, L/2))
>>> bottom_face = ParametricRegion((x, y, -L/2), (x, -L/2, L/2), (y, -L/2, L/2))
>>> flux = vector_integrate(C.x**2*C.k, top_face) + vector_integrate(C.y**2*C.k, bottom_face)
>>> flux
L**4/6
验证斯托克斯定理
参见斯托克斯定理
示例 1
>>> from sympy.vector import curl
>>> curve = ParametricRegion((cos(theta), sin(theta)), (theta, 0, pi/2))
>>> surface = ParametricRegion((r*cos(theta), r*sin(theta)), (r, 0, 1), (theta, 0, pi/2))
>>> F = C.y*C.i + C.z*C.k + C.x*C.k
>>>
>>> vector_integrate(F, curve)
-pi/4
>>> vector_integrate(curl(F), surface)
-pi/4
示例 2
>>> circle = ParametricRegion((cos(theta), sin(theta), 1), (theta, 0, 2*pi))
>>> cone = ParametricRegion((r*cos(theta), r*sin(theta), r), (r, 0, 1), (theta, 0, 2*pi))
>>> cone = ParametricRegion((r*cos(theta), r*sin(theta), r), (r, 0, 1), (theta, 0, 2*pi))
>>> f = (-C.y**3/3 + sin(C.x))*C.i + (C.x**3/3 + cos(C.y))*C.j + C.x*C.y*C.z*C.k
>>> vector_integrate(f, circle)
pi/2
>>> vector_integrate(curl(f), cone)
pi/2
验证散度定理
参见散度定理
示例 1
>>> from sympy.vector import divergence
>>> sphere = ParametricRegion((4*sin(phi)*cos(theta),4*sin(phi)*sin(theta), 4*cos(phi)),
... (phi, 0, pi), (theta, 0, 2*pi))
>>> solidsphere = ParametricRegion((r*sin(phi)*cos(theta),r*sin(phi)*sin(theta), r*cos(phi)),
... (r, 0, 4),(phi, 0, pi), (theta, 0, 2*pi))
>>> field = C.x**3*C.i + C.y**3*C.j + C.z**3*C.k
>>> vector_integrate(field, sphere)
12288*pi/5
>>> vector_integrate(divergence(field), solidsphere)
12288*pi/5
示例 2
>>> cube = ParametricRegion((x, y, z), (x, 0, 1), (y, 0, 1), (z, 0, 1))
>>> field = 2*C.x*C.y*C.i + 3*C.x*C.y*C.j + C.z*exp(C.x + C.y)*C.k
>>> vector_integrate(divergence(field), cube)
-E + 7/2 + E*(-1 + E)
向量 API
-
sympy.vector 中的关键类(文档字符串)
-
定向器类(文档字符串)
-
sympy.vector 中的关键函数(文档字符串)
sympy.vector 中的基本类(文档字符串)
class sympy.vector.coordsysrect.CoordSys3D(name, transformation=None, parent=None, location=None, rotation_matrix=None, vector_names=None, variable_names=None)
代表三维空间中的坐标系。
__init__(name, location=None, rotation_matrix=None, parent=None, vector_names=None, variable_names=None, latex_vects=None, pretty_vects=None, latex_scalars=None, pretty_scalars=None, transformation=None)
如果此系统在某个方向或位置相对于另一个定义,则方向/位置参数是必需的。
参数:
name : str
新 CoordSys3D 实例的名称。
transformation : Lambda, Tuple, str
根据变换方程定义的转换或从预定义的转换中选择的转换。
location : Vector
新系统原点相对于父实例的位置向量。
rotation_matrix : SymPy ImmutableMatrix
新坐标系的旋转矩阵,相对于父坐标系。换句话说,这是 new_system.rotation_matrix(parent) 的输出。
parent : CoordSys3D
相对于其方向/位置(或两者)正在定义的坐标系。
vector_names, variable_names : iterable(optional)
每个都是包含 3 个字符串的迭代器,分别用于新系统的基本向量和基本标量的自定义名称。用于简单的字符串打印。
create_new(name, transformation, variable_names=None, vector_names=None)
返回一个通过变换与自身连接的 CoordSys3D。
参数:
name : str
新 CoordSys3D 实例的名称。
transformation : Lambda, Tuple, str
根据变换方程定义的转换或从预定义的转换中选择的转换。
vector_names, variable_names : iterable(optional)
每个都是包含 3 个字符串的迭代器,分别用于新系统的基本向量和基本标量的自定义名称。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> a = CoordSys3D('a')
>>> b = a.create_new('b', transformation='spherical')
>>> b.transformation_to_parent()
(b.r*sin(b.theta)*cos(b.phi), b.r*sin(b.phi)*sin(b.theta), b.r*cos(b.theta))
>>> b.transformation_from_parent()
(sqrt(a.x**2 + a.y**2 + a.z**2), acos(a.z/sqrt(a.x**2 + a.y**2 + a.z**2)), atan2(a.y, a.x))
locate_new(name, position, vector_names=None, variable_names=None)
返回一个 CoordSys3D,其原点位于给定位置相对于此坐标系原点的位置。
参数:
name : str
新 CoordSys3D 实例的名称。
position : Vector
相对于此系统原点的新系统原点的位置向量。
vector_names, variable_names : iterable(optional)
每个都是包含 3 个字符串的迭代器,分别用于新系统的基本向量和基本标量的自定义名称。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> A = CoordSys3D('A')
>>> B = A.locate_new('B', 10 * A.i)
>>> B.origin.position_wrt(A.origin)
10*A.i
orient_new(name, orienters, location=None, vector_names=None, variable_names=None)
使用用户指定的方式创建一个与此系统相关的新 CoordSys3D。
请参阅有关定向程序的定向器类文档以获取更多信息。
参数:
name : str
新 CoordSys3D 实例的名称。
orienters : iterable/Orienter
一个 Orienter 或 Orienter 的迭代器,用于定向新坐标系。如果提供了一个 Orienter,则应用它以获得新系统。如果提供了一个可迭代对象,则按照它们出现的顺序应用定向器。
location : Vector(optional)
新坐标系原点相对于此系统原点的位置。如果未指定,则认为原点重合。
vector_names, variable_names : iterable(optional)
每个都是包含 3 个字符串的迭代器,分别用于新系统的基本向量和基本标量的自定义名称。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q0, q1, q2, q3 = symbols('q0 q1 q2 q3')
>>> N = CoordSys3D('N')
使用 AxisOrienter
>>> from sympy.vector import AxisOrienter
>>> axis_orienter = AxisOrienter(q1, N.i + 2 * N.j)
>>> A = N.orient_new('A', (axis_orienter, ))
使用 BodyOrienter
>>> from sympy.vector import BodyOrienter
>>> body_orienter = BodyOrienter(q1, q2, q3, '123')
>>> B = N.orient_new('B', (body_orienter, ))
使用 SpaceOrienter
>>> from sympy.vector import SpaceOrienter
>>> space_orienter = SpaceOrienter(q1, q2, q3, '312')
>>> C = N.orient_new('C', (space_orienter, ))
使用 QuaternionOrienter
>>> from sympy.vector import QuaternionOrienter
>>> q_orienter = QuaternionOrienter(q0, q1, q2, q3)
>>> D = N.orient_new('D', (q_orienter, ))
orient_new_axis(name, angle, axis, location=None, vector_names=None, variable_names=None)
轴旋转是围绕任意轴的旋转,旋转角度由 SymPy 表达式标量提供,轴由矢量提供。
参数:
name : 字符串
新坐标系的名称
angle : 表达式
新系统旋转的角度
axis : 矢量
执行旋转的轴
location : 矢量(可选)
新坐标系的原点位置相对于该系统的原点。如果未指定,则认为原点重合。
vector_names, variable_names : 可迭代对象(可选)
每个具有 3 个字符串的可迭代对象,分别为新系统的基向量和基标量的自定义名称。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q1 = symbols('q1')
>>> N = CoordSys3D('N')
>>> B = N.orient_new_axis('B', q1, N.i + 2 * N.j)
orient_new_body(name, angle1, angle2, angle3, rotation_order, location=None, vector_names=None, variable_names=None)
Body orientation 通过连续三个简单旋转带此坐标系。
Body fixed rotations 包括欧拉角和 Tait-Bryan 角,请参见 zh.wikipedia.org/wiki/欧拉角
。
参数:
name : 字符串
新坐标系的名称
angle1, angle2, angle3 : 表达式
连续三个角度旋转坐标系
rotation_order : 字符串
定义旋转轴顺序的字符串
location : 矢量(可选)
新坐标系的原点位置相对于该系统的原点。如果未指定,则认为原点重合。
vector_names, variable_names : 可迭代对象(可选)
每个具有 3 个字符串的可迭代对象,分别为新系统的基向量和基标量的自定义名称。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q1, q2, q3 = symbols('q1 q2 q3')
>>> N = CoordSys3D('N')
'Body' 固定旋转由三个角度和三个固定于身体的旋转轴描述。为了将坐标系 D 相对于 N 定向,每个连续的旋转总是围绕固定在 D 上的正交单位向量进行的。例如,'123' 旋转将指定关于 N.i、然后 D.j、然后 D.k 的旋转。(最初,D.i 与 N.i 相同)因此,
>>> D = N.orient_new_body('D', q1, q2, q3, '123')
与之相同
>>> D = N.orient_new_axis('D', q1, N.i)
>>> D = D.orient_new_axis('D', q2, D.j)
>>> D = D.orient_new_axis('D', q3, D.k)
可接受的旋转顺序长度为 3,表达为 XYZ 或 123,并且不能连续两次围绕同一轴旋转。
>>> B = N.orient_new_body('B', q1, q2, q3, '123')
>>> B = N.orient_new_body('B', q1, q2, 0, 'ZXZ')
>>> B = N.orient_new_body('B', 0, 0, 0, 'XYX')
orient_new_quaternion(name, q0, q1, q2, q3, location=None, vector_names=None, variable_names=None)
四元数方向用四元数使新的 CoordSys3D 定向,由 lambda,一个单位向量,以某个量 theta 进行有限旋转定义。
此方向由四个参数描述:
q0 = cos(theta/2)
q1 = lambda_x sin(theta/2)
q2 = lambda_y sin(theta/2)
q3 = lambda_z sin(theta/2)
四元数不接受旋转顺序。
参数:
name : 字符串
新坐标系的名称
q0, q1, q2, q3 : 表达式
用于旋转坐标系的四元数
location : 矢量(可选)
新坐标系的原点位置相对于该系统的原点。如果未指定,则认为原点重合。
vector_names, variable_names : 可迭代对象(可选)
每个具有新系统的基向量和基标量的自定义名称的 3 个字符串的可迭代。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q0, q1, q2, q3 = symbols('q0 q1 q2 q3')
>>> N = CoordSys3D('N')
>>> B = N.orient_new_quaternion('B', q0, q1, q2, q3)
orient_new_space(name, angle1, angle2, angle3, rotation_order, location=None, vector_names=None, variable_names=None)
空间旋转类似于体旋转,但是旋转顺序相反。
参数:
name:字符串
新坐标系的名称
angle1, angle2, angle3:Expr
用于旋转坐标系的三个连续角度
rotation_order:字符串
定义旋转轴顺序的字符串
location:向量(可选)
新坐标系原点的位置相对于此系统原点的位置。如果未指定,则假定原点重合。
vector_names, variable_names:可迭代(可选)
每个具有新系统的基向量和基标量的自定义名称的 3 个字符串的可迭代。用于简单的字符串打印。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q1, q2, q3 = symbols('q1 q2 q3')
>>> N = CoordSys3D('N')
要将坐标系 D 定向到 N,每个顺序旋转始终围绕 N 的正交单位向量进行。例如,'123'旋转将指定围绕 N.i,然后 N.j,然后 N.k 的旋转。因此,
>>> D = N.orient_new_space('D', q1, q2, q3, '312')
与...相同
>>> B = N.orient_new_axis('B', q1, N.i)
>>> C = B.orient_new_axis('C', q2, N.j)
>>> D = C.orient_new_axis('D', q3, N.k)
另请参阅
CoordSys3D.orient_new_body
通过欧拉角定向的方法
position_wrt(other)
返回此坐标系的原点位置向量与另一个点/CoordSys3D 的原点的位置向量之间的位置向量。
参数:
other:点/CoordSys3D
如果 other 是一个点,则返回此系统原点相对于其的位置。如果其是 CoordSyRect 的实例,则返回相对于其原点的位置。
示例
>>> from sympy.vector import CoordSys3D
>>> N = CoordSys3D('N')
>>> N1 = N.locate_new('N1', 10 * N.i)
>>> N.position_wrt(N1)
(-10)*N.i
rotation_matrix(other)
返回此坐标系与另一个系统之间的方向余弦矩阵(DCM),也称为‘旋转矩阵’。
如果 v_a 是在系统‘A’中定义的向量(以矩阵格式),v_b 是在系统‘B’中定义的相同向量,则 v_a = A.rotation_matrix(B) * v_b。
返回一个 SymPy 矩阵。
参数:
other:CoordSys3D
生成 DCM 的系统。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q1 = symbols('q1')
>>> N = CoordSys3D('N')
>>> A = N.orient_new_axis('A', q1, N.i)
>>> N.rotation_matrix(A)
Matrix([
[1, 0, 0],
[0, cos(q1), -sin(q1)],
[0, sin(q1), cos(q1)]])
scalar_map(other)
返回一个表达此框架的坐标变量(基标量)与 otherframe 的变量相关的字典。
参数:
otherframe:CoordSys3D
映射变量到其他系统。
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import Symbol
>>> A = CoordSys3D('A')
>>> q = Symbol('q')
>>> B = A.orient_new_axis('B', q, A.k)
>>> A.scalar_map(B)
{A.x: B.x*cos(q) - B.y*sin(q), A.y: B.x*sin(q) + B.y*cos(q), A.z: B.z}
class sympy.vector.vector.Vector(*args)
所有向量类的超类。理想情况下,用户不应该实例化此类或其任何子类。
property components
返回此向量的分量,以 Python 字典形式将 BaseVector 实例映射到相应的测量数。
示例
>>> from sympy.vector import CoordSys3D
>>> C = CoordSys3D('C')
>>> v = 3*C.i + 4*C.j + 5*C.k
>>> v.components
{C.i: 3, C.j: 4, C.k: 5}
cross(other)
返回此向量与另一个向量或二重实例的叉积。如果‘other’是向量,则叉积是一个向量。如果‘other’是二重,这将返回一个二重实例。
参数:
other:向量/二重
我们正在交叉的向量或二重的。
示例
>>> from sympy.vector import CoordSys3D
>>> C = CoordSys3D('C')
>>> C.i.cross(C.j)
C.k
>>> C.i ^ C.i
0
>>> v = 3*C.i + 4*C.j + 5*C.k
>>> v ^ C.i
5*C.j + (-4)*C.k
>>> d = C.i.outer(C.i)
>>> C.j.cross(d)
(-1)*(C.k|C.i)
dot(other)
返回此向量与另一个向量、二阶张量或梯度算子的点积。如果‘other’是一个向量,则返回点积标量(SymPy 表达式)。如果‘other’是一个二阶张量,则返回点积作为一个向量。如果‘other’是 Del 的实例,则返回 Python 函数形式的方向导数算子。如果将此函数应用于标量表达式,则返回标量场相对于此向量的方向导数。
参数:
other: 向量/二阶张量/梯度算子
我们正在与之点乘的向量或二阶张量,或者是一个梯度算子。
示例
>>> from sympy.vector import CoordSys3D, Del
>>> C = CoordSys3D('C')
>>> delop = Del()
>>> C.i.dot(C.j)
0
>>> C.i & C.i
1
>>> v = 3*C.i + 4*C.j + 5*C.k
>>> v.dot(C.k)
5
>>> (C.i & delop)(C.x*C.y*C.z)
C.y*C.z
>>> d = C.i.outer(C.i)
>>> C.i.dot(d)
C.i
magnitude()
返回此向量的大小。
normalize()
返回此向量的归一化版本。
outer(other)
返回此向量与另一个向量的外积,以一个二阶张量实例的形式。
参数:
other : 向量
与之进行外积计算的向量。
示例
>>> from sympy.vector import CoordSys3D
>>> N = CoordSys3D('N')
>>> N.i.outer(N.j)
(N.i|N.j)
projection(other, scalar=False)
返回‘self’上‘other’的向量或标量投影。
示例
>>> from sympy.vector.coordsysrect import CoordSys3D
>>> C = CoordSys3D('C')
>>> i, j, k = C.base_vectors()
>>> v1 = i + j + k
>>> v2 = 3*i + 4*j
>>> v1.projection(v2)
7/3*C.i + 7/3*C.j + 7/3*C.k
>>> v1.projection(v2, scalar=True)
7/3
separate()
这个向量在不同坐标系中的成分,根据其定义。
返回一个字典,将每个 CoordSys3D 映射到相应的成分向量。
示例
>>> from sympy.vector import CoordSys3D
>>> R1 = CoordSys3D('R1')
>>> R2 = CoordSys3D('R2')
>>> v = R1.i + R2.i
>>> v.separate() == {R1: R1.i, R2: R2.i}
True
to_matrix(system)
返回此向量相对于指定坐标系的矩阵形式。
参数:
system : 三维坐标系
计算矩阵形式的系统
示例
>>> from sympy.vector import CoordSys3D
>>> C = CoordSys3D('C')
>>> from sympy.abc import a, b, c
>>> v = a*C.i + b*C.j + c*C.k
>>> v.to_matrix(C)
Matrix([
[a],
[b],
[c]])
class sympy.vector.dyadic.Dyadic(*args)
所有二阶张量类的超类。
参考文献
[R1074]
[R1075]
Kane, T., Levinson, D. 动力学理论与应用. 1985 McGraw-Hill
property components
返回此二阶张量的分量,以 Python 字典形式映射 BaseDyadic 实例到相应的测量数。
cross(other)
返回此二阶张量与一个向量的叉乘,作为一个向量实例。
参数:
other : 向量
我们正在与此二阶张量进行叉乘的向量。
示例
>>> from sympy.vector import CoordSys3D
>>> N = CoordSys3D('N')
>>> d = N.i.outer(N.i)
>>> d.cross(N.j)
(N.i|N.k)
dot(other)
返回此二阶张量与另一个二阶张量或向量的点积(也称为内积)。如果‘other’是一个二阶张量,则返回一个二阶张量。否则,返回一个向量(除非出现错误)。
参数:
other : 二阶张量/向量
与之进行内积运算的其他二阶张量或向量
示例
>>> from sympy.vector import CoordSys3D
>>> N = CoordSys3D('N')
>>> D1 = N.i.outer(N.j)
>>> D2 = N.j.outer(N.j)
>>> D1.dot(D2)
(N.i|N.j)
>>> D1.dot(N.j)
N.i
to_matrix(system, second_system=None)
返回与一个或两个坐标系相关的二阶张量的矩阵形式。
参数:
system : 三维坐标系
矩阵的行和列对应的坐标系。如果提供第二个系统,则仅对应矩阵的行。
second_system : 三维坐标系,可选,默认为 None
矩阵列对应的坐标系。
示例
>>> from sympy.vector import CoordSys3D
>>> N = CoordSys3D('N')
>>> v = N.i + 2*N.j
>>> d = v.outer(N.i)
>>> d.to_matrix(N)
Matrix([
[1, 0, 0],
[2, 0, 0],
[0, 0, 0]])
>>> from sympy import Symbol
>>> q = Symbol('q')
>>> P = N.orient_new_axis('P', q, N.k)
>>> d.to_matrix(N, P)
Matrix([
[ cos(q), -sin(q), 0],
[2*cos(q), -2*sin(q), 0],
[ 0, 0, 0]])
class sympy.vector.deloperator.Del
表示向量微分算子,通常在数学表达式中表示为‘nabla’符号。
cross(vect, doit=False)
表示此算子与给定向量的叉乘 - 等同于向量场的旋度。
参数:
vect : 向量
要计算其旋度的向量。
doit : 布尔值
如果为 True,则在调用每个分量的.doit()后返回结果。否则,返回的表达式包含 Derivative 实例。
示例
>>> from sympy.vector import CoordSys3D, Del
>>> C = CoordSys3D('C')
>>> delop = Del()
>>> v = C.x*C.y*C.z * (C.i + C.j + C.k)
>>> delop.cross(v, doit = True)
(-C.x*C.y + C.x*C.z)*C.i + (C.x*C.y - C.y*C.z)*C.j +
(-C.x*C.z + C.y*C.z)*C.k
>>> (delop ^ C.i).doit()
0
dot(vect, doit=False)
表示该运算符与给定向量的点积,等于向量场的散度。
参数:
vect : 向量
要计算其散度的向量。
doit : bool
如果为 True,则在调用每个分量的.doit()后返回结果。否则,返回的表达式包含 Derivative 实例。
示例
>>> from sympy.vector import CoordSys3D, Del
>>> delop = Del()
>>> C = CoordSys3D('C')
>>> delop.dot(C.x*C.i)
Derivative(C.x, C.x)
>>> v = C.x*C.y*C.z * (C.i + C.j + C.k)
>>> (delop & v).doit()
C.x*C.y + C.x*C.z + C.y*C.z
gradient(scalar_field, doit=False)
返回给定标量场的梯度,作为 Vector 实例。
参数:
scalar_field : SymPy 表达式
要计算其梯度的标量场。
doit : bool
如果为 True,则在调用每个分量的.doit()后返回结果。否则,返回的表达式包含 Derivative 实例。
示例
>>> from sympy.vector import CoordSys3D, Del
>>> C = CoordSys3D('C')
>>> delop = Del()
>>> delop.gradient(9)
0
>>> delop(C.x*C.y*C.z).doit()
C.y*C.z*C.i + C.x*C.z*C.j + C.x*C.y*C.k
class sympy.vector.parametricregion.ParametricRegion(definition, *bounds)
表示空间中的参数区域。
参数:
definition : 用于根据参数定义基础标量的元组。
bounds : 用于定义参数及其相应下限和上限的参数或长度为 3 的元组。
示例
>>> from sympy import cos, sin, pi
>>> from sympy.abc import r, theta, t, a, b, x, y
>>> from sympy.vector import ParametricRegion
>>> ParametricRegion((t, t**2), (t, -1, 2))
ParametricRegion((t, t**2), (t, -1, 2))
>>> ParametricRegion((x, y), (x, 3, 4), (y, 5, 6))
ParametricRegion((x, y), (x, 3, 4), (y, 5, 6))
>>> ParametricRegion((r*cos(theta), r*sin(theta)), (r, -2, 2), (theta, 0, pi))
ParametricRegion((r*cos(theta), r*sin(theta)), (r, -2, 2), (theta, 0, pi))
>>> ParametricRegion((a*cos(t), b*sin(t)), t)
ParametricRegion((a*cos(t), b*sin(t)), t)
>>> circle = ParametricRegion((r*cos(theta), r*sin(theta)), r, (theta, 0, pi))
>>> circle.parameters
(r, theta)
>>> circle.definition
(r*cos(theta), r*sin(theta))
>>> circle.limits
{theta: (0, pi)}
参数化区域的维数确定区域是曲线、曲面还是体积区域。它不表示空间中的维数。
>>> circle.dimensions
1
class sympy.vector.implicitregion.ImplicitRegion(variables, equation)
表示空间中的隐式区域。
参数:
variables : 用于将隐式方程中的变量映射到基础标量的元组。
equation : 表示区域隐式方程的表达式或等式。
示例
>>> from sympy import Eq
>>> from sympy.abc import x, y, z, t
>>> from sympy.vector import ImplicitRegion
>>> ImplicitRegion((x, y), x**2 + y**2 - 4)
ImplicitRegion((x, y), x**2 + y**2 - 4)
>>> ImplicitRegion((x, y), Eq(y*x, 1))
ImplicitRegion((x, y), x*y - 1)
>>> parabola = ImplicitRegion((x, y), y**2 - 4*x)
>>> parabola.degree
2
>>> parabola.equation
-4*x + y**2
>>> parabola.rational_parametrization(t)
(4/t**2, 4/t)
>>> r = ImplicitRegion((x, y, z), Eq(z, x**2 + y**2))
>>> r.variables
(x, y, z)
>>> r.singular_points()
EmptySet
>>> r.regular_point()
(-10, -10, 200)
multiplicity(point)
返回区域上奇点的多重性。
区域的奇点(x,y)如果所有 m-1 阶偏导数在此处为零,则称其为多重性 m。
示例
>>> from sympy.abc import x, y, z
>>> from sympy.vector import ImplicitRegion
>>> I = ImplicitRegion((x, y, z), x**2 + y**3 - z**4)
>>> I.singular_points()
{(0, 0, 0)}
>>> I.multiplicity((0, 0, 0))
2
rational_parametrization(parameters=('t', 's'), reg_point=None)
返回隐式区域的有理参数化。
示例
>>> from sympy import Eq
>>> from sympy.abc import x, y, z, s, t
>>> from sympy.vector import ImplicitRegion
>>> parabola = ImplicitRegion((x, y), y**2 - 4*x)
>>> parabola.rational_parametrization()
(4/t**2, 4/t)
>>> circle = ImplicitRegion((x, y), Eq(x**2 + y**2, 4))
>>> circle.rational_parametrization()
(4*t/(t**2 + 1), 4*t**2/(t**2 + 1) - 2)
>>> I = ImplicitRegion((x, y), x**3 + x**2 - y**2)
>>> I.rational_parametrization()
(t**2 - 1, t*(t**2 - 1))
>>> cubic_curve = ImplicitRegion((x, y), x**3 + x**2 - y**2)
>>> cubic_curve.rational_parametrization(parameters=(t))
(t**2 - 1, t*(t**2 - 1))
>>> sphere = ImplicitRegion((x, y, z), x**2 + y**2 + z**2 - 4)
>>> sphere.rational_parametrization(parameters=(t, s))
(-2 + 4/(s**2 + t**2 + 1), 4*s/(s**2 + t**2 + 1), 4*t/(s**2 + t**2 + 1))
对于某些圆锥曲线,regular_points()无法找到曲线上的点。在这种情况下,用户需要确定区域上的一个点,并使用 reg_point 传递它以计算参数化表示。
>>> c = ImplicitRegion((x, y), (x - 1/2)**2 + (y)**2 - (1/4)**2)
>>> c.rational_parametrization(reg_point=(3/4, 0))
(0.75 - 0.5/(t**2 + 1), -0.5*t/(t**2 + 1))
参考文献
- Christoph M. Hoffmann,“参数曲线和曲面之间的转换方法”,普渡大学 e-Pubs,1990 年。可查看:
docs.lib.purdue.edu/cgi/viewcontent.cgi?article=1827&context=cstech
regular_point()
返回隐式区域上的一个点。
示例
>>> from sympy.abc import x, y, z
>>> from sympy.vector import ImplicitRegion
>>> circle = ImplicitRegion((x, y), (x + 2)**2 + (y - 3)**2 - 16)
>>> circle.regular_point()
(-2, -1)
>>> parabola = ImplicitRegion((x, y), x**2 - 4*y)
>>> parabola.regular_point()
(0, 0)
>>> r = ImplicitRegion((x, y, z), (x + y + z)**4)
>>> r.regular_point()
(-10, -10, 20)
参考文献
- Erik Hillgarter,“圆锥曲线上的有理点”,学位论文,RISC-Linz,约翰·开普勒林茨大学,1996 年。可查看:
www3.risc.jku.at/publications/download/risc_1355/Rational%20Points%20on%20Conics.pdf
singular_points()
返回区域的奇点集合。
区域上的奇点是区域上所有偏导数均为零的点。
示例
>>> from sympy.abc import x, y
>>> from sympy.vector import ImplicitRegion
>>> I = ImplicitRegion((x, y), (y-1)**2 -x**3 + 2*x**2 -x)
>>> I.singular_points()
{(1, 1)}
class sympy.vector.integrals.ParametricIntegral(field, parametricregion)
表示标量或矢量场在参数区域上的积分。
示例
>>> from sympy import cos, sin, pi
>>> from sympy.vector import CoordSys3D, ParametricRegion, ParametricIntegral
>>> from sympy.abc import r, t, theta, phi
>>> C = CoordSys3D('C')
>>> curve = ParametricRegion((3*t - 2, t + 1), (t, 1, 2))
>>> ParametricIntegral(C.x, curve)
5*sqrt(10)/2
>>> length = ParametricIntegral(1, curve)
>>> length
sqrt(10)
>>> semisphere = ParametricRegion((2*sin(phi)*cos(theta), 2*sin(phi)*sin(theta), 2*cos(phi)), (theta, 0, 2*pi), (phi, 0, pi/2))
>>> ParametricIntegral(C.z, semisphere)
8*pi
>>> ParametricIntegral(C.j + C.k, ParametricRegion((r*cos(theta), r*sin(theta)), r, theta))
0
定向器类(文档字符串)
原文:
docs.sympy.org/latest/modules/vector/api/orienterclasses.html
class sympy.vector.orienters.Orienter(*args)
所有定向器类的超类。
rotation_matrix()
与此定向器实例对应的旋转矩阵。
class sympy.vector.orienters.AxisOrienter(angle, axis)
表示轴定向器的类。
__init__(angle, axis)
轴旋转是围绕任意轴的旋转,角度由 SymPy 表达式标量提供,轴由向量提供。
参数:
角度:Expr
用于旋转顺序的角度
轴:向量
需要执行旋转的轴
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q1 = symbols('q1')
>>> N = CoordSys3D('N')
>>> from sympy.vector import AxisOrienter
>>> orienter = AxisOrienter(q1, N.i + 2 * N.j)
>>> B = N.orient_new('B', (orienter, ))
rotation_matrix(system)
与此定向器实例对应的旋转矩阵。
参数:
系统:CoordSys3D
计算旋转矩阵的坐标系
class sympy.vector.orienters.BodyOrienter(angle1, angle2, angle3, rot_order)
表示体定向器的类。
__init__(angle1, angle2, angle3, rot_order)
体定向将此坐标系通过三个连续的简单旋转。
'Body'固定旋转包括欧拉角和泰特-布莱恩角,参见en.wikipedia.org/wiki/Euler_angles
。
参数:
角度 1,角度 2,角度 3:Expr
三个连续角度来旋转坐标系
旋转顺序:字符串
定义旋转轴顺序的字符串
示例
>>> from sympy.vector import CoordSys3D, BodyOrienter
>>> from sympy import symbols
>>> q1, q2, q3 = symbols('q1 q2 q3')
>>> N = CoordSys3D('N')
'Body'固定旋转由三个角度和三个固定于 D 的体旋转轴描述。为了将坐标系 D 定向到 N,每次连续旋转都是关于固定于 D 的正交单位向量。例如,'123'旋转将指定关于 N.i、然后 D.j、然后 D.k 的旋转。(最初,D.i 与 N.i 相同)因此,
>>> body_orienter = BodyOrienter(q1, q2, q3, '123')
>>> D = N.orient_new('D', (body_orienter, ))
同上
>>> from sympy.vector import AxisOrienter
>>> axis_orienter1 = AxisOrienter(q1, N.i)
>>> D = N.orient_new('D', (axis_orienter1, ))
>>> axis_orienter2 = AxisOrienter(q2, D.j)
>>> D = D.orient_new('D', (axis_orienter2, ))
>>> axis_orienter3 = AxisOrienter(q3, D.k)
>>> D = D.orient_new('D', (axis_orienter3, ))
可接受的旋转顺序长度为 3,表示为 XYZ 或 123,并且不能连续两次围绕同一轴旋转。
>>> body_orienter1 = BodyOrienter(q1, q2, q3, '123')
>>> body_orienter2 = BodyOrienter(q1, q2, 0, 'ZXZ')
>>> body_orienter3 = BodyOrienter(0, 0, 0, 'XYX')
class sympy.vector.orienters.SpaceOrienter(angle1, angle2, angle3, rot_order)
表示空间定向器的类。
__init__(angle1, angle2, angle3, rot_order)
空间旋转类似于体旋转,但是旋转的顺序相反。
参数:
角度 1,角度 2,角度 3:Expr
三个连续角度来旋转坐标系
旋转顺序:字符串
定义旋转轴顺序的字符串
示例
>>> from sympy.vector import CoordSys3D, SpaceOrienter
>>> from sympy import symbols
>>> q1, q2, q3 = symbols('q1 q2 q3')
>>> N = CoordSys3D('N')
为了将坐标系 D 定向到 N,每次连续旋转都是关于 N 的正交单位向量。例如,'123'旋转将指定关于 N.i、然后 N.j、然后 N.k 的旋转。因此,
>>> space_orienter = SpaceOrienter(q1, q2, q3, '312')
>>> D = N.orient_new('D', (space_orienter, ))
同上
>>> from sympy.vector import AxisOrienter
>>> axis_orienter1 = AxisOrienter(q1, N.i)
>>> B = N.orient_new('B', (axis_orienter1, ))
>>> axis_orienter2 = AxisOrienter(q2, N.j)
>>> C = B.orient_new('C', (axis_orienter2, ))
>>> axis_orienter3 = AxisOrienter(q3, N.k)
>>> D = C.orient_new('C', (axis_orienter3, ))
另请参阅
BodyOrienter
相对于欧拉角定向系统的定向器。
class sympy.vector.orienters.QuaternionOrienter(q0, q1, q2, q3)
表示四元数定向器的类。
__init__(angle1, angle2, angle3, rot_order)
四元数定向使用四元数将新的 CoordSys3D 定向,定义为围绕单位向量 lambda 的有限旋转,旋转量为 theta。
这个方向由四个参数描述:
q0 = cos(theta/2)
q1 = lambda_x sin(theta/2)
q2 = lambda_y sin(theta/2)
q3 = lambda_z sin(theta/2)
四元数不接受旋转顺序。
参数:
q0, q1, q2, q3:Expr
用于旋转坐标系的四元数
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import symbols
>>> q0, q1, q2, q3 = symbols('q0 q1 q2 q3')
>>> N = CoordSys3D('N')
>>> from sympy.vector import QuaternionOrienter
>>> q_orienter = QuaternionOrienter(q0, q1, q2, q3)
>>> B = N.orient_new('B', (q_orienter, ))
sympy.vector 中的基本功能(文档字符串)
原文:
docs.sympy.org/latest/modules/vector/api/vectorfunctions.html
sympy.vector.matrix_to_vector(matrix, system)
将矩阵形式的向量转换为 Vector 实例。
假设矩阵的元素表示‘system’的基向量上的向量分量的测量数。
参数:
矩阵:SymPy Matrix,尺寸:(3, 1)
要转换为向量的矩阵
system:CoordSys3D
定义向量的坐标系统
示例
>>> from sympy import ImmutableMatrix as Matrix
>>> m = Matrix([1, 2, 3])
>>> from sympy.vector import CoordSys3D, matrix_to_vector
>>> C = CoordSys3D('C')
>>> v = matrix_to_vector(m, C)
>>> v
C.i + 2*C.j + 3*C.k
>>> v.to_matrix(C) == m
True
sympy.vector.express(expr, system, system2=None, variables=False)
用于‘表达’功能的全局函数。
在给定的坐标系中重新表达向量、二重或标量(可用 sympyfiable)。
如果‘variables’为 True,则将向量/标量场或二重标量中其他坐标系统的坐标变量(基标量)也用给定系统的基标量代换。
参数:
表达式:向量/二重/标量(可用 sympyfiable)
要在 CoordSys3D ‘system’ 中重新表达的表达式
system: CoordSys3D
要表达的表达式的坐标系统
system2: CoordSys3D
重新表达所需的其他坐标系统(仅适用于 Dyadic Expr)
variables:布尔值
指定是否要用参数系统中的表达式中存在的坐标变量替换它们
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy import Symbol, cos, sin
>>> N = CoordSys3D('N')
>>> q = Symbol('q')
>>> B = N.orient_new_axis('B', q, N.k)
>>> from sympy.vector import express
>>> express(B.i, N)
(cos(q))*N.i + (sin(q))*N.j
>>> express(N.x, B, variables=True)
B.x*cos(q) - B.y*sin(q)
>>> d = N.i.outer(N.i)
>>> express(d, B, N) == (cos(q))*(B.i|N.i) + (-sin(q))*(B.j|N.i)
True
sympy.vector.curl(vect, doit=True)
返回相对于给定坐标系的基标量计算的向量场的旋度。
参数:
向量:向量
向量操作数
doit:bool
如果为 True,则在每个分量上调用 .doit() 后返回结果。否则,返回表达式包含 Derivative 实例
示例
>>> from sympy.vector import CoordSys3D, curl
>>> R = CoordSys3D('R')
>>> v1 = R.y*R.z*R.i + R.x*R.z*R.j + R.x*R.y*R.k
>>> curl(v1)
0
>>> v2 = R.x*R.y*R.z*R.i
>>> curl(v2)
R.x*R.y*R.j + (-R.x*R.z)*R.k
sympy.vector.divergence(vect, doit=True)
返回相对于给定坐标系的基标量计算的向量场的散度。
参数:
向量:向量
向量操作数
doit:bool
如果为 True,则在每个分量上调用 .doit() 后返回结果。否则,返回表达式包含 Derivative 实例
示例
>>> from sympy.vector import CoordSys3D, divergence
>>> R = CoordSys3D('R')
>>> v1 = R.x*R.y*R.z * (R.i+R.j+R.k)
>>> divergence(v1)
R.x*R.y + R.x*R.z + R.y*R.z
>>> v2 = 2*R.y*R.z*R.j
>>> divergence(v2)
2*R.z
sympy.vector.gradient(scalar_field, doit=True)
返回相对于给定坐标系的基标量计算的标量场的向量梯度。
参数:
标量场:SymPy Expr
要计算其梯度的标量场
doit:bool
如果为 True,则在每个分量上调用 .doit() 后返回结果。否则,返回表达式包含 Derivative 实例
示例
>>> from sympy.vector import CoordSys3D, gradient
>>> R = CoordSys3D('R')
>>> s1 = R.x*R.y*R.z
>>> gradient(s1)
R.y*R.z*R.i + R.x*R.z*R.j + R.x*R.y*R.k
>>> s2 = 5*R.x**2*R.z
>>> gradient(s2)
10*R.x*R.z*R.i + 5*R.x**2*R.k
sympy.vector.is_conservative(field)
检查场是否保守。
参数:
场:向量
要检查其保守属性的场
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy.vector import is_conservative
>>> R = CoordSys3D('R')
>>> is_conservative(R.y*R.z*R.i + R.x*R.z*R.j + R.x*R.y*R.k)
True
>>> is_conservative(R.z*R.j)
False
sympy.vector.is_solenoidal(field)
检查场是否为旋量场。
参数:
场:向量
要检查旋量性质的场
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy.vector import is_solenoidal
>>> R = CoordSys3D('R')
>>> is_solenoidal(R.y*R.z*R.i + R.x*R.z*R.j + R.x*R.y*R.k)
True
>>> is_solenoidal(R.y * R.j)
False
sympy.vector.scalar_potential(field, coord_sys)
返回给定坐标系中场的标量势函数(不包括添加的积分常数)。
参数:
场:向量
要计算其标量势函数的向量场
coord_sys:CoordSys3D
进行计算的坐标系统
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy.vector import scalar_potential, gradient
>>> R = CoordSys3D('R')
>>> scalar_potential(R.k, R) == R.z
True
>>> scalar_field = 2*R.x**2*R.y*R.z
>>> grad_field = gradient(scalar_field)
>>> scalar_potential(grad_field, R)
2*R.x**2*R.y*R.z
sympy.vector.scalar_potential_difference(field, coord_sys, point1, point2)
返回在给定坐标系中,关于给定场的两点之间的标量势差。
如果提供了标量场,则考虑它在两个点处的值。如果提供了保守向量场,则使用其在两点处的标量势函数的值。
返回(点 2 处的电势)-(点 1 处的电势)
两个点的位置向量是相对于提供的坐标系原点计算的。
参数:
field : 向量/表达式
要计算的场
coord_sys : CoordSys3D
在进行计算时使用的坐标系
point1 : 点
在给定的坐标系中的初始点
position2 : 点
给定坐标系中的第二个点
示例
>>> from sympy.vector import CoordSys3D
>>> from sympy.vector import scalar_potential_difference
>>> R = CoordSys3D('R')
>>> P = R.origin.locate_new('P', R.x*R.i + R.y*R.j + R.z*R.k)
>>> vectfield = 4*R.x*R.y*R.i + 2*R.x**2*R.j
>>> scalar_potential_difference(vectfield, R, R.origin, P)
2*R.x**2*R.y
>>> Q = R.origin.locate_new('O', 3*R.i + R.j + 2*R.k)
>>> scalar_potential_difference(vectfield, R, P, Q)
-2*R.x**2*R.y + 18
sympy.vector.integrals.vector_integrate(field, *region)
计算在区域或一组参数上的向量/标量场的积分。
示例
>>> from sympy.vector import CoordSys3D, ParametricRegion, vector_integrate
>>> from sympy.abc import x, y, t
>>> C = CoordSys3D('C')
>>> region = ParametricRegion((t, t**2), (t, 1, 5))
>>> vector_integrate(C.x*C.i, region)
12
还可以计算几何模块中某些对象上的积分。
>>> from sympy.geometry import Point, Circle, Triangle
>>> c = Circle(Point(0, 2), 5)
>>> vector_integrate(C.x**2 + C.y**2, c)
290*pi
>>> triangle = Triangle(Point(-2, 3), Point(2, 3), Point(0, 5))
>>> vector_integrate(3*C.x**2*C.y*C.i + C.j, triangle)
-8
可以计算一些简单隐式区域上的积分。但在大多数情况下,计算它们需要太长时间。这是因为参数表示的表达式变得很大。
>>> from sympy.vector import ImplicitRegion
>>> c2 = ImplicitRegion((x, y), (x - 2)**2 + (y - 1)**2 - 9)
>>> vector_integrate(1, c2)
6*pi
与基础标量相关的场的积分:
>>> vector_integrate(12*C.y**3, (C.y, 1, 3))
240
>>> vector_integrate(C.x**2*C.z, C.x)
C.x**3*C.z/3
>>> vector_integrate(C.x*C.i - C.y*C.k, C.x)
(Integral(C.x, C.x))*C.i + (Integral(-C.y, C.x))*C.k
>>> _.doit()
C.x**2/2*C.i + (-C.x*C.y)*C.k
数论
原文:
docs.sympy.org/latest/reference/public/numbertheory/index.html
目录
- 数论
数论
Ntheory 类参考
class sympy.ntheory.generate.Sieve(sieve_interval=1000000)
质数列表,实现为动态增长的埃拉托斯特尼筛。当请求涉及尚未筛选的奇数时,筛网自动扩展至该数。实现细节限制了质数数量为 2³²-1
。
示例
>>> from sympy import sieve
>>> sieve._reset() # this line for doctest only
>>> 25 in sieve
False
>>> sieve._list
array('L', [2, 3, 5, 7, 11, 13, 17, 19, 23])
extend(n)
扩展筛以涵盖所有小于等于 n 的质数。
示例
>>> from sympy import sieve
>>> sieve._reset() # this line for doctest only
>>> sieve.extend(30)
>>> sieve[10] == 29
True
extend_to_no(i)
扩展以包括第 i 个质数。
参数:
i:整数
示例
>>> from sympy import sieve
>>> sieve._reset() # this line for doctest only
>>> sieve.extend_to_no(9)
>>> sieve._list
array('L', [2, 3, 5, 7, 11, 13, 17, 19, 23])
注意事项
如果列表过短,会增加 50%,因此它很可能比请求的要长。
mobiusrange(a, b)
生成范围 [a, b) 内的所有莫比乌斯函数数。
参数:
a:整数
范围内的第一个数字
b:整数
超出范围的第一个数字
示例
>>> from sympy import sieve
>>> print([i for i in sieve.mobiusrange(7, 18)])
[-1, 0, 0, 1, -1, 0, -1, 1, 1, 0, -1]
primerange(a, b=None)
生成范围 [2, a) 或 [a, b) 内的所有质数。
示例
>>> from sympy import sieve, prime
小于 19 的所有质数:
>>> print([i for i in sieve.primerange(19)])
[2, 3, 5, 7, 11, 13, 17]
大于等于 7 且小于 19 的所有质数:
>>> print([i for i in sieve.primerange(7, 19)])
[7, 11, 13, 17]
通过第 10 个质数的所有质数
>>> list(sieve.primerange(prime(10) + 1))
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
search(n)
返回界定 n 的质数的索引 i, j。
如果 n 是质数,则 i == j。
虽然 n 可以是一个表达式,但如果 ceiling 不能将其转换为整数,则会引发一个 n 错误。
示例
>>> from sympy import sieve
>>> sieve.search(25)
(9, 10)
>>> sieve.search(23)
(9, 9)
totientrange(a, b)
生成范围 [a, b) 内的所有欧拉函数数。
示例
>>> from sympy import sieve
>>> print([i for i in sieve.totientrange(7, 18)])
[6, 4, 6, 4, 10, 4, 12, 6, 8, 8, 16]
Ntheory 函数参考
sympy.ntheory.generate.prime(nth)
返回第 n 个质数,质数索引为 prime(1) = 2, prime(2) = 3,等等… 第 n 个质数大约是 (n\log(n))。
对于 x 的对数积分是小于等于 x 的质数数的一个相当不错的近似值,即 li(x) ~ pi(x) 实际上,对于我们关心的数字( x<1e11 ),li(x) - pi(x) < 50000
此外,可以安全地假设对于此函数可以评估的数字,li(x) > pi(x)。
在这里,我们使用二分查找找到最小的整数 m,使得 li(m) > n。现在 pi(m-1) < li(m-1) <= n,
我们使用 primepi
函数找到 pi(m - 1)。
从 m 开始,我们必须找到 n - pi(m-1) 更多的质数。
对于此实现可以处理的输入,我们最多需要测试约 10**5 个数的素性以获得答案。
示例
>>> from sympy import prime
>>> prime(10)
29
>>> prime(1)
2
>>> prime(100000)
1299709
另请参阅
sympy.ntheory.primetest.isprime
测试 n 是否为质数
primerange
在给定范围内生成所有质数
primepi
返回小于或等于 n 的质数的数量
参考资料
[R648]
en.wikipedia.org/wiki/Prime_number_theorem#Table_of_.CF.80.28x.29.2C_x_.2F_log_x.2C_and_li.28x.29
[R649]
en.wikipedia.org/wiki/Prime_number_theorem#Approximations_for_the_nth_prime_number
[R650]
en.wikipedia.org/wiki/Skewes%27_number
sympy.ntheory.generate.primepi(n)
表示质数计数函数 pi(n) = 小于或等于 n 的质数的数量。
自版本 1.13 开始不推荐使用:primepi
函数已弃用。请改用sympy.functions.combinatorial.numbers.primepi
。有关更多信息,请参阅其文档。有关详细信息,请参阅从 ntheory 移动到 functions 的符号函数。
算法描述:
在筛法中,我们移除所有质数 p 的倍数,除了 p 本身。
让 phi(i,j)为从小于或等于 j 的质数中筛除后剩余的 2 <= k <= i 的整数数量。显然,pi(n) = phi(n, sqrt(n))
如果 j 不是一个质数,phi(i,j) = phi(i, j - 1)
如果 j 是一个质数,我们删除所有最小质因数为 j 的数(除了 j 本身)。
让(x= j \times a)为这样一个数,其中(2 \le a \le i / j)。现在,在从质数(\le j - 1)中筛选后,a 必须保留(因为 x,因此 a 没有质因数(\le j - 1))。显然,有 phi(i / j, j - 1)个这样的 a 在从质数(\le j - 1)中筛选后保留。
如果 a 是小于等于 j-1 的质数,(x= j \times a) 的最小质因数为 a,并且已经被移除(通过筛法从 a 中)。因此,我们不需要再次移除它。(注意:这样的 x 有 pi(j-1)个)
因此,将被移除的 x 的数量为:phi(i / j, j - 1) - phi(j - 1, j - 1)(注意 pi(j - 1) = phi(j - 1, j - 1))
(\Rightarrow) phi(i,j) = phi(i, j - 1) - phi(i / j, j - 1) + phi(j - 1, j - 1)
因此,使用以下递归并实现 dp:
phi(a, b) = phi(a, b - 1),如果 b 不是质数 phi(a, b) = phi(a, b-1)-phi(a / b, b-1) + phi(b-1, b-1),如果 b 是质数
显然 a 总是形式为 floor(n / k),最多可以取(2\sqrt{n})个值。维护两个数组 arr1,arr2,arr1[i] = phi(i, j),arr2[i] = phi(n // i, j)
最后的答案是 arr2[1]
示例
>>> from sympy import primepi, prime, prevprime, isprime
>>> primepi(25)
9
因此,小于或等于 25 的有 9 个质数。25 是质数吗?
>>> isprime(25)
False
不是。因此,小于 25 的第一个质数必须是第 9 个质数:
>>> prevprime(25) == prime(9)
True
参见
sympy.ntheory.primetest.isprime
测试 n 是否为质数
primerange
在给定范围内生成所有的质数
prime
返回第 n 个质数
sympy.ntheory.generate.nextprime(n, ith=1)
返回大于 n 的第 ith 个质数。
参数:
n : 整数
ith : 正整数
返回:
int : 返回大于 n 的第 ith 个质数
异常:
ValueError
如果
ith <= 0
。如果n
或ith
不是整数。
注意
潜在的质数位于 6*j +/- 1。这个性质在搜索中被使用。
>>> from sympy import nextprime
>>> [(i, nextprime(i)) for i in range(10, 15)]
[(10, 11), (11, 13), (12, 13), (13, 17), (14, 17)]
>>> nextprime(2, ith=2) # the 2nd prime after 2
5
参见
prevprime
返回小于 n 的最大素数
primerange
生成给定范围内的所有素数
sympy.ntheory.generate.prevprime(n)
返回小于 n 的最大素数。
注释
潜在的素数位于 6*j +/- 1。这个性质在搜索时使用。
>>> from sympy import prevprime
>>> [(i, prevprime(i)) for i in range(10, 15)]
[(10, 7), (11, 7), (12, 11), (13, 11), (14, 13)]
另请参阅
nextprime
返回大于 n 的第 i 个素数
primerange
生成给定范围内的所有素数
sympy.ntheory.generate.primerange(a, b=None)
生成范围为 [2, a) 或 [a, b) 的所有素数列表。
如果默认筛选器中存在范围,则将从那里返回值;否则将返回值但不会修改筛选器。
示例
>>> from sympy import primerange, prime
小于 19 的所有素数:
>>> list(primerange(19))
[2, 3, 5, 7, 11, 13, 17]
大于或等于 7 且小于 19 的所有素数:
>>> list(primerange(7, 19))
[7, 11, 13, 17]
所有小于第 10 个素数的素数
>>> list(primerange(prime(10) + 1))
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
筛法方法 primerange 通常更快,但由于筛选器存储值,内存占用更多。可以使用名为 sieve 的默认筛选器实例:
>>> from sympy import sieve
>>> list(sieve.primerange(1, 30))
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
注释
一些关于素数出现频率的著名猜想是 [1]:
-
孪生素数:尽管经常不是,以下将给出两个素数
无限次数:
primerange(6n - 1, 6n + 2)
-
勒让德定理:以下始终至少产生一个素数
primerange(n2, (n+1)2+1)
-
贝特朗定理(已证明):在范围内总是有一个素数
primerange(n, 2*n)
-
布洛卡德定理:在给定范围内至少有四个素数
primerange(prime(n)2, prime(n+1)2)
素数之间的平均间隔是 log(n) [2];素数之间的间隔可以任意大,因为合数序列可以任意大,例如序列 n! + 2, n! + 3 … n! + n 中的数字都是合数。
另请参阅
prime
返回第 n 个素数
nextprime
返回大于 n 的第 i 个素数
prevprime
返回小于 n 的最大素数
randprime
返回给定范围内的随机素数
primorial
基于条件返回素数的乘积
Sieve.primerange
返回从已计算的素数范围或扩展筛选器以包含所请求范围的范围。
参考文献
[R651]
en.wikipedia.org/wiki/Prime_number
[R652]
primes.utm.edu/notes/gaps.html
sympy.ntheory.generate.randprime(a, b)
返回在区间 a, b) 中的随机素数。
贝特朗定理确保 randprime(a, 2*a) 对于 a > 1 总是成功。
请注意,由于实现困难,选择的素数不是均匀随机的。例如,在范围 [112, 128) 中有两个素数 113
和 127
,但 randprime(112, 128)
以 15/17 的概率返回 127
。
示例
>>> from sympy import randprime, isprime
>>> randprime(1, 30)
13
>>> isprime(randprime(1, 30))
True
另请参阅
[primerange
在给定范围内生成所有素数
参考资料
[R653]
sympy.ntheory.generate.primorial(n, nth=True)
返回前 n 个素数的乘积(默认)或小于等于 n 的素数(当 nth=False
时)。
示例
>>> from sympy.ntheory.generate import primorial, primerange
>>> from sympy import factorint, Mul, primefactors, sqrt
>>> primorial(4) # the first 4 primes are 2, 3, 5, 7
210
>>> primorial(4, nth=False) # primes <= 4 are 2 and 3
6
>>> primorial(1)
2
>>> primorial(1, nth=False)
1
>>> primorial(sqrt(101), nth=False)
210
人们可以认为素数是无限的,因为如果你取一组素数并将它们相乘(例如,素数阶乘),然后加或减 1,结果不能被任何原始因子之一整除,因此这个素数的乘积必须被除以 1 或多个新的素数。
在这种情况下,数本身是一个新的素数:
>>> factorint(primorial(4) + 1)
{211: 1}
在这种情况下,两个新的素数是因子:
>>> factorint(primorial(4) - 1)
{11: 1, 19: 1}
这里,乘在一起的素数比这些素数更大或更小:
>>> p = list(primerange(10, 20))
>>> sorted(set(primefactors(Mul(*p) + 1)).difference(set(p)))
[2, 5, 31, 149]
另请参阅
primerange
在给定范围内生成所有素数
sympy.ntheory.generate.cycle_length(f, x0, nmax=None, values=False)
对于给定的迭代序列,返回一个生成器,提供迭代周期的长度(lambda)和循环开始之前的项的长度(mu);如果 values
为 True,则将返回序列的项。序列以值 x0
开始。
注意:可能会返回多于第一个 lambda + mu 项,这是使用布伦特方法进行循环检测的成本;然而,通常比使用弗洛伊德方法确定的正确结束点计算的项要少。
>>> from sympy.ntheory.generate import cycle_length
这将产生 i <– func(i) 的连续值:
>>> def gen(func, i):
... while 1:
... yield i
... i = func(i)
...
函数定义如下:
>>> func = lambda i: (i**2 + 1) % 51
并给定种子 4 和计算的 mu 和 lambda 项:
>>> next(cycle_length(func, 4))
(6, 3)
我们可以看到输出的含义:
>>> iter = cycle_length(func, 4, values=True)
>>> list(iter)
[4, 17, 35, 2, 5, 26, 14, 44, 50, 2, 5, 26, 14]
第一个 3 之后有 6 个重复值。
如果怀疑一个序列的长度比你希望的长,可以使用 nmax
提前退出(mu 将返回为 None):
>>> next(cycle_length(func, 4, nmax = 4))
(4, None)
>>> list(cycle_length(func, 4, nmax = 4, values=True))
[4, 17, 35, 2]
代码修改自:
sympy.ntheory.generate.composite(nth)
返回第 n 个复合数,其中复合数索引为 composite(1) = 4,composite(2) = 6,等等……
示例
>>> from sympy import composite
>>> composite(36)
52
>>> composite(1)
4
>>> composite(17737)
20000
另请参阅
sympy.ntheory.primetest.isprime
测试 n 是否为素数
primerange
在给定范围内生成所有素数
primepi
返回小于或等于 n 的素数的数量
prime
返回第 n 个素数
compositepi
返回小于或等于 n
的正复合数的数量
sympy.ntheory.generate.compositepi(n)
返回小于或等于 n
的正复合数的数量。第一个正复合数是 4
,即 compositepi(4) = 1
。
示例
>>> from sympy import compositepi
>>> compositepi(25)
15
>>> compositepi(1000)
831
另见
sympy.ntheory.primetest.isprime
测试 n
是否为素数
primerange
在给定范围内生成所有素数
prime
返回第 n
个质数
primepi
返回小于或等于 n
的素数的数量
composite
返回第 n
个复合数
sympy.ntheory.factor_.smoothness(n)
返回 n
的 B-平滑和 B-幂平滑值。
n
的平滑性是 n
的最大素因子;幂平滑性是最大除数升至其重数。
示例
>>> from sympy.ntheory.factor_ import smoothness
>>> smoothness(2**7*3**2)
(3, 128)
>>> smoothness(2**4*13)
(13, 16)
>>> smoothness(2)
(2, 2)
另见
factorint
, smoothness_p
sympy.ntheory.factor_.smoothness_p(n, m=-1, power=0, visual=None)
返回列表 [m, (p, (M, sm(p + m), psm(p + m)))…]
其中:
-
p**M
是n
的基数p
除数 -
sm(p + m)
是p + m
的平滑性(默认情况下m = -1
) -
psm(p + m)
是p + m
的幂平滑性
列表按平滑性(默认)或按幂平滑性(如果 power=1
)排序。
左侧(m = -1
)或右侧(m = 1
)的平滑性决定了从 p +/- 1
类型的因式方法获得的结果。
>>> from sympy.ntheory.factor_ import smoothness_p, factorint
>>> smoothness_p(10431, m=1)
(1, [(3, (2, 2, 4)), (19, (1, 5, 5)), (61, (1, 31, 31))])
>>> smoothness_p(10431)
(-1, [(3, (2, 2, 2)), (19, (1, 3, 9)), (61, (1, 5, 5))])
>>> smoothness_p(10431, power=1)
(-1, [(3, (2, 2, 2)), (61, (1, 5, 5)), (19, (1, 3, 9))])
如果 visual=True
,则返回带注释的字符串:
>>> print(smoothness_p(21477639576571, visual=1))
p**i=4410317**1 has p-1 B=1787, B-pow=1787
p**i=4869863**1 has p-1 B=2434931, B-pow=2434931
该字符串也可以直接从因数分解字典生成,反之亦然:
>>> factorint(17*9)
{3: 2, 17: 1}
>>> smoothness_p(_)
'p**i=3**2 has p-1 B=2, B-pow=2\np**i=17**1 has p-1 B=2, B-pow=16'
>>> smoothness_p(_)
{3: 2, 17: 1}
输出逻辑的表格如下:
可视化 输入 真 --- --- dict str str str tuple str n str mul str
另见
factorint
, smoothness
sympy.ntheory.factor_.multiplicity(p, n)
找到最大整数 m
,使得 p**m
可整除 n
。
示例
>>> from sympy import multiplicity, Rational
>>> [multiplicity(5, n) for n in [8, 5, 25, 125, 250]]
[0, 1, 2, 3, 3]
>>> multiplicity(3, Rational(1, 9))
-2
注意:在检查大阶乘中数字的重复性时,最有效的方法是将其作为未评估的阶乘发送或直接调用 multiplicity_in_factorial
:
>>> from sympy.ntheory import multiplicity_in_factorial
>>> from sympy import factorial
>>> p = factorial(25)
>>> n = 2**100
>>> nfac = factorial(n, evaluate=False)
>>> multiplicity(p, nfac)
52818775009509558395695966887
>>> _ == multiplicity_in_factorial(p, n)
True
另见
trailing
sympy.ntheory.factor_.perfect_power(n, candidates=None, big=True, factor=True)
返回 (b, e)
,使得若 n
是具有 e > 1
的唯一完全幂,则 n
== b**e
,否则返回 False
(例如 1
不是完全幂)。如果 n
不是有理数,则引发 ValueError
。
默认情况下,基数将递归分解并收集指数,因此寻找最大可能的 e
。如果 big=False
则选择最小可能的 e
(因此是素数)。
如果factor=True
,则尝试同时因式分解n
,因为找到一个因子表明n
的唯一可能根。这是默认值,因为在搜索完美幂的过程中只会测试几个小因子。
使用candidates
主要是为了内部使用;如果提供,则如果n
不能写成其中一个候选数的幂,则返回 False,并且不会尝试因式分解(超出测试因子 2)。
示例
>>> from sympy import perfect_power, Rational
>>> perfect_power(16)
(2, 4)
>>> perfect_power(16, big=False)
(4, 2)
负数只能有奇数完全幂:
>>> perfect_power(-4)
False
>>> perfect_power(-8)
(-2, 3)
有理数也被识别:
>>> perfect_power(Rational(1, 2)**3)
(1/2, 3)
>>> perfect_power(Rational(-3, 2)**3)
(-3/2, 3)
注释
要知道一个整数是否是 2 的完全幂,使用
>>> is2pow = lambda n: bool(n and not n & (n - 1))
>>> [(i, is2pow(i)) for i in range(5)]
[(0, False), (1, True), (2, True), (3, False), (4, True)]
不需要提供candidates
。提供时,假定它们是整数。第一个大于计算的最大可能指数的候选数将为该程序的常规信号。
>>> perfect_power(3**8, [9])
False
>>> perfect_power(3**8, [2, 4, 8])
(3, 8)
>>> perfect_power(3**8, [4, 8], big=False)
(9, 4)
参见
sympy.core.intfunc.integer_nthroot
,sympy.ntheory.primetest.is_square
sympy.ntheory.factor_.pollard_rho(n, s=2, a=1, retries=5, seed=1234, max_steps=None, F=None)
使用波拉德ρ方法尝试提取n
的非平凡因子。返回的因子可能是复合数。如果找不到因子,则返回None
。
算法使用生成器函数生成 x 的伪随机值,并用 F(x)替换 x。如果未提供 F,则使用函数 x**2 + a
。提供给 F(x)的第一个值是s
。如果失败(如果retries
> 0),将提供新的a
和s
;如果提供了 F,则将忽略a
。
这些函数生成的数字序列通常会引导到某个数字,然后循环回到该数字并开始重复序列,例如 1, 2, 3, 4, 5, 3, 4, 5 – 这种引导和循环看起来有点像希腊字母ρ,因此得名‘rho’。
对于给定的函数,可能会得到非常不同的主循环值,因此允许重试是个好主意:
>>> from sympy.ntheory.generate import cycle_length
>>> n = 16843009
>>> F = lambda x:(2048*pow(x, 2, n) + 32767) % n
>>> for s in range(5):
... print('loop length = %4i; leader length = %3i' % next(cycle_length(F, s)))
...
loop length = 2489; leader length = 43
loop length = 78; leader length = 121
loop length = 1482; leader length = 100
loop length = 1482; leader length = 286
loop length = 1482; leader length = 101
这是一个明确的例子,其中有一个三元素引导到一个序列的 3 个数字(11, 14, 4),然后重复:
>>> x=2
>>> for i in range(9):
... print(x)
... x=(x**2+12)%17
...
2
16
13
11
14
4
11
14
4
>>> next(cycle_length(lambda x: (x**2+12)%17, 2))
(3, 3)
>>> list(cycle_length(lambda x: (x**2+12)%17, 2, values=True))
[2, 16, 13, 11, 14, 4]
不是检查所有生成的值与 n 的差异的替代品,仅检查第 k 和 2*k 个数字,例如第一个和第二个,第二个和第四个,第三个和第六个,直到检测到循环已遍历。在 rho 找到因子或报告失败之前,循环可能会有数千个步骤长。如果指定了max_steps
,则在指定的步数后取消迭代并失败。
示例
>>> from sympy import pollard_rho
>>> n=16843009
>>> F=lambda x:(2048*pow(x,2,n) + 32767) % n
>>> pollard_rho(n, F=F)
257
使用默认设置和糟糕的a
值和无重试:
>>> pollard_rho(n, a=n-2, retries=0)
如果重试> 0,则当为 a 生成新值时,问题可能会得到纠正:
>>> pollard_rho(n, a=n-2, retries=1)
257
参考
[R654]
Richard Crandall & Carl Pomerance(2005 年),《素数:计算视角》,Springer,第二版,229-231
sympy.ntheory.factor_.pollard_pm1(n, B=10, a=2, retries=0, seed=1234)
使用波拉德 p-1 方法尝试提取n
的非平凡因子。返回一个除数(可能是复合数)或None
。
a
的值是在测试 gcd(a**M - 1, n) 中使用的基数。默认值为 2. 如果 retries
> 0,则在第一次尝试未找到因子后,将随机生成一个新的 a
(使用 seed
),并重复该过程。
注意:M 的值是 lcm(1..B) = reduce(ilcm, range(2, B + 1))。
寻找具有比 B
小的幂平滑度的偶数旁边的因子。选择更大的 B 增加找到较大因子的可能性,但需要更长时间。无论是否找到 n 的因子取决于 a
和略小于因子 p 的偶数的幂平滑度(因此称为 p - 1)。
尽管有关什么构成一个好的 a
的讨论,有些描述很难解释。在下面引用的模块化数学网站中指出,如果 gcd(aM - 1, n) = N,则对于 N 的每个素数幂除数,aM % q**r 等于 1。但请考虑以下情况:
>>> from sympy.ntheory.factor_ import smoothness_p, pollard_pm1
>>> n=257*1009
>>> smoothness_p(n)
(-1, [(257, (1, 2, 256)), (1009, (1, 7, 16))])
因此,我们应该(并且可以)找到一个具有 B=16 的根:
>>> pollard_pm1(n, B=16, a=3)
1009
如果我们尝试将 B 增加到 256,我们发现它不起作用:
>>> pollard_pm1(n, B=256)
>>>
但是,如果 a
的值改变,我们发现只有 257 的倍数有效,例如:
>>> pollard_pm1(n, B=256, a=257)
1009
检查不同的 a
值表明所有未能成功的值都具有不等于 n
而等于其中一个因子的 gcd 值:
>>> from sympy import ilcm, igcd, factorint, Pow
>>> M = 1
>>> for i in range(2, 256):
... M = ilcm(M, i)
...
>>> set([igcd(pow(a, M, n) - 1, n) for a in range(2, 256) if
... igcd(pow(a, M, n) - 1, n) != n])
{1009}
但是,对于 n 的每个除数,aM % d 是否都等于 1?
>>> aM = pow(255, M, n)
>>> [(d, aM%Pow(*d.args)) for d in factorint(n, visual=True).args]
[(257**1, 1), (1009**1, 1)]
不,只有其中一个。因此,也许原则是,在给定 B 值的情况下,可以找到一个根,只要:
-
p - 1 值旁边的幂平滑度不超过 B
-
对于每个 n 的除数,a**M % p != 1。
通过尝试多个 a
可能会有一个能够产生一个因子。
例子
使用默认平滑界限,这个数字无法破解:
>>> from sympy.ntheory import pollard_pm1
>>> pollard_pm1(21477639576571)
增加平滑界限有助于:
>>> pollard_pm1(21477639576571, B=2000)
4410317
查看该数字的因子的平滑度时,我们发现:
>>> from sympy.ntheory.factor_ import smoothness_p, factorint
>>> print(smoothness_p(21477639576571, visual=1))
p**i=4410317**1 has p-1 B=1787, B-pow=1787
p**i=4869863**1 has p-1 B=2434931, B-pow=2434931
p - 1 的因子化的除数的 B 和 B-pow 是相同的,因为这些因子化有一个非常大的素因子:
>>> factorint(4410317 - 1)
{2: 2, 617: 1, 1787: 1}
>>> factorint(4869863-1)
{2: 1, 2434931: 1}
注意,直到 B 达到 B-pow 值为 1787,该数字才能被破解;
>>> pollard_pm1(21477639576571, B=1786)
>>> pollard_pm1(21477639576571, B=1787)
4410317
B 值与除数旁边的数字的因子有关,而不是除数本身。最坏的情况是,旁边的数字 p 有一个大的素数除数或者是一个完美的幂。如果这些条件适用,则幂平滑度将约为 p/2 或 p。更现实的情况是,旁边将有一个大的素数因子需要一个大约 p/2 的 B 值。尽管可能已经在这个级别搜索了素数,但 p/2 是 p - 1 的一个因子,我们不知道。下面的模块化数学引用表示,在 1015 到 1515 + 104 范围内的数的 15% 是 106 平滑的,因此在该范围内 B=106 的失败率为 85%。从 108 到 108 + 103,百分比几乎颠倒了…但在该范围内,简单的试除法非常快。
参考文献
[R655]
Richard Crandall & Carl Pomerance (2005), “Prime Numbers: A Computational Perspective”, Springer, 第二版, 236-238
[R656]
[R657]
www.cs.toronto.edu/~yuvalf/Factorization.pdf
sympy.ntheory.factor_.factorint(n, limit=None, use_trial=True, use_rho=True, use_pm1=True, use_ecm=True, verbose=False, visual=None, multiple=False)
给定正整数n
,factorint(n)
返回一个字典,其中包含n
的素因子作为键和它们的重数作为值。例如:
>>> from sympy.ntheory import factorint
>>> factorint(2000) # 2000 = (2**4) * (5**3)
{2: 4, 5: 3}
>>> factorint(65537) # This number is prime
{65537: 1}
对于小于 2 的输入,factorint
的行为如下:
factorint(1)
返回空因式分解{}
factorint(0)
返回{0:1}
factorint(-n)
将-1:1
添加到因子中,然后因式分解n
部分因式分解:
如果指定了limit
(>3),则在执行试除法达到(包括)限制(或进行相应数量的 rho/p-1 步骤)后停止搜索。如果有一个大数,并且只想找到小因子(如果有的话),这是很有用的。注意,设置限制并不阻止较大的因子提前被发现;它仅仅意味着最大因子可能是复合的。由于检查完全幂的成本相对较低,因此无论限制设置如何,都会执行此检查。
例如,此数字有两个小因子和一个难以轻易分解的巨大半素数因子:
>>> from sympy.ntheory import isprime
>>> a = 1407633717262338957430697921446883
>>> f = factorint(a, limit=10000)
>>> f == {991: 1, int(202916782076162456022877024859): 1, 7: 1}
True
>>> isprime(max(f))
False
此数有一个小因子和一个基数大于限制的剩余完全幂:
>>> factorint(3*101**7, limit=5)
{3: 1, 101: 7}
因子列表:
如果设置了multiple
为True
,则返回包含素因子及其重数的列表。
>>> factorint(24, multiple=True)
[2, 2, 2, 3]
视觉因式分解:
如果visual
设置为True
,则会返回整数的视觉因式分解。例如:
>>> from sympy import pprint
>>> pprint(factorint(4200, visual=True))
3 1 2 1
2 *3 *5 *7
注意,这是通过在 Mul 和 Pow 中使用 evaluate=False 标志实现的。如果您对 evaluate=False 的表达式进行其他操作,它可能会计算。因此,如果您希望执行因式分解后的操作,请仅在可视化时使用 visual 选项,并且在 visual=False 时使用正常的字典返回。
您可以通过将它们发送回factorint
轻松地在这两种形式之间切换。
>>> from sympy import Mul
>>> regular = factorint(1764); regular
{2: 2, 3: 2, 7: 2}
>>> pprint(factorint(regular))
2 2 2
2 *3 *7
>>> visual = factorint(1764, visual=True); pprint(visual)
2 2 2
2 *3 *7
>>> print(factorint(visual))
{2: 2, 3: 2, 7: 2}
如果您希望发送一个部分因式分解形式的数字,可以使用字典或未评估的表达式进行发送:
>>> factorint(factorint({4: 2, 12: 3})) # twice to toggle to dict form
{2: 10, 3: 3}
>>> factorint(Mul(4, 12, evaluate=False))
{2: 4, 3: 1}
输出逻辑表格如下:
Input True --- --- dict mul n mul mul mul
注意事项
算法:
该函数在多个算法之间切换。试除法快速找到小因子(约 1-5 位数字的顺序),如果给定足够时间,它会找到所有大因子。波拉德·罗和 p-1 算法用于提前找到大因子;它们通常在几秒钟内找到大约 10 位数字的因子:
>>> factors = factorint(12345678910111213141516)
>>> for base, exp in sorted(factors.items()):
... print('%s %s' % (base, exp))
...
2 2
2507191691 1
1231026625769 1
以上任何方法均可选择性地使用以下布尔参数禁用:
use_trial
:切换使用试除法use_rho
:切换使用波拉德·罗方法use_pm1
:切换使用 Pollard 的 p-1 方法
factorint
还定期检查剩余部分是否为质数或完全幂,并在这些情况下停止。
对于未评估的阶乘,它使用 Legendre 的公式(定理)。
如果 verbose
设置为 True
,则会打印详细进度。
另请参见
smoothness
, smoothness_p
, divisors
sympy.ntheory.factor_.factorrat(rat, limit=None, use_trial=True, use_rho=True, use_pm1=True, verbose=False, visual=None, multiple=False)
给定有理数 r
,factorrat(r)
返回一个包含 r
的质因子为键和它们相应重数为值的字典。例如:
>>> from sympy import factorrat, S
>>> factorrat(S(8)/9) # 8/9 = (2**3) * (3**-2)
{2: 3, 3: -2}
>>> factorrat(S(-1)/987) # -1/789 = -1 * (3**-1) * (7**-1) * (47**-1)
{-1: 1, 3: -1, 7: -1, 47: -1}
请参阅 factorint
的文档字符串,了解以下关键字的详细说明和示例:
limit
:执行试除的整数限制use_trial
:切换使用试除法use_rho
:切换使用 Pollard 的 rho 方法use_pm1
:切换使用 Pollard 的 p-1 方法verbose
:切换详细打印进度multiple
:切换返回因子列表或字典visual
:切换输出的乘积形式
sympy.ntheory.factor_.primefactors(n, limit=None, verbose=False, **kwargs)
返回 n 的质因子的排序列表,忽略重复性和任何如果限制太低而残留的复合因子的完全因子分解。与 factorint()
不同,primefactors()
不返回-1 或 0。
参数:
n:整数
**limit, verbose, kwargs:
要传递给
factorint
的额外关键字参数。由于kwargs
是 1.13 版本中的新内容,因此保留limit
和verbose
以确保兼容性。
返回:
list(int):整数 n 的质数列表
示例
>>> from sympy.ntheory import primefactors, factorint, isprime
>>> primefactors(6)
[2, 3]
>>> primefactors(-5)
[5]
>>> sorted(factorint(123456).items())
[(2, 6), (3, 1), (643, 1)]
>>> primefactors(123456)
[2, 3, 643]
>>> sorted(factorint(10000000001, limit=200).items())
[(101, 1), (99009901, 1)]
>>> isprime(99009901)
False
>>> primefactors(10000000001, limit=300)
[101]
另请参见
factorint
, divisors
sympy.ntheory.factor_.divisors(n, generator=False, proper=False)
默认情况下,返回从 1 到 n 排序的 n 的所有除数。如果 generator
设置为 True
,则返回一个无序生成器。
如果有许多质因子(计算重复因子),则 n 的除数数量可能会很大。如果只需要因子的数量,请使用 divisor_count(n)
。
示例
>>> from sympy import divisors, divisor_count
>>> divisors(24)
[1, 2, 3, 4, 6, 8, 12, 24]
>>> divisor_count(24)
8
>>> list(divisors(120, generator=True))
[1, 2, 4, 8, 3, 6, 12, 24, 5, 10, 20, 40, 15, 30, 60, 120]
注释
这是稍作修改的 Tim Peters 版本,参考自:stackoverflow.com/questions/1010381/python-factorization
另请参见
primefactors
, factorint
, divisor_count
sympy.ntheory.factor_.proper_divisors(n, generator=False)
默认情况下,返回除 n 以外的 n 的所有除数。如果 generator
设置为 True
,则返回一个无序生成器。
示例
>>> from sympy import proper_divisors, proper_divisor_count
>>> proper_divisors(24)
[1, 2, 3, 4, 6, 8, 12]
>>> proper_divisor_count(24)
7
>>> list(proper_divisors(120, generator=True))
[1, 2, 4, 8, 3, 6, 12, 24, 5, 10, 20, 40, 15, 30, 60]
另请参见
factorint
, divisors
, proper_divisor_count
sympy.ntheory.factor_.divisor_count(n, modulus=1, proper=False)
返回 n
的除数的数量。如果 modulus
不为 1,则仅计算可被 modulus
整除的除数。如果 proper
为 True
,则 n
的除数将不被计算。
示例
>>> from sympy import divisor_count
>>> divisor_count(6)
4
>>> divisor_count(6, 2)
2
>>> divisor_count(6, proper=True)
3
另见
factorint
, divisors
, totient
, proper_divisor_count
sympy.ntheory.factor_.proper_divisor_count(n, modulus=1)
返回 n
的真除数的数量。
示例
>>> from sympy import proper_divisor_count
>>> proper_divisor_count(6)
3
>>> proper_divisor_count(6, modulus=2)
1
另见
divisors
, proper_divisors
, divisor_count
sympy.ntheory.factor_.udivisors(n, generator=False)
默认情况下,返回排序后的 n 的所有单位除数。如果 generator
为 True
,则返回一个无序生成器。
如果存在许多质因数,则 n
的单位除数数量可能非常大。如果只需单位除数的数量,请使用 udivisor_count(n)
。
示例
>>> from sympy.ntheory.factor_ import udivisors, udivisor_count
>>> udivisors(15)
[1, 3, 5, 15]
>>> udivisor_count(15)
4
>>> sorted(udivisors(120, generator=True))
[1, 3, 5, 8, 15, 24, 40, 120]
另见
primefactors
, factorint
, divisors
, divisor_count
, udivisor_count
参考文献
[R658]
en.wikipedia.org/wiki/Unitary_divisor
[R659]
mathworld.wolfram.com/UnitaryDivisor.html
sympy.ntheory.factor_.udivisor_count(n)
返回 n
的单位除数的数量。
参数:
n:整数
示例
>>> from sympy.ntheory.factor_ import udivisor_count
>>> udivisor_count(120)
8
另见
factorint
, divisors
, udivisors
, divisor_count
, totient
参考文献
[R660]
mathworld.wolfram.com/UnitaryDivisorFunction.html
sympy.ntheory.factor_.antidivisors(n, generator=False)
默认情况下,返回排序后的 n 的所有反除数。
n 的反除数[R661]是不以最大可能边界除 n 的数。如果 generator 为 True,则返回无序生成器。
示例
>>> from sympy.ntheory.factor_ import antidivisors
>>> antidivisors(24)
[7, 16]
>>> sorted(antidivisors(128, generator=True))
[3, 5, 15, 17, 51, 85]
另见
primefactors
, factorint
, divisors
, divisor_count
, antidivisor_count
参考文献
[R661] (1,2)
定义见oeis.org/A066272/a066272a.html
sympy.ntheory.factor_.antidivisor_count(n)
返回 n
的反除数数目[R662]。
参数:
n:整数
示例
>>> from sympy.ntheory.factor_ import antidivisor_count
>>> antidivisor_count(13)
4
>>> antidivisor_count(27)
5
另见
factorint
, divisors
, antidivisors
, divisor_count
, totient
参考文献
[R662] (1,2)
公式来自oeis.org/A066272
sympy.ntheory.factor_.totient(n)
计算欧拉 totient 函数 phi(n)
自 1.13 版本起不推荐使用:totient
函数已被弃用。请使用sympy.functions.combinatorial.numbers.totient
。更多信息请参见其文档。详细信息请参见将 ntheory 中的符号函数重定位到 functions。
totient(n)
或 (\phi(n)) 是小于等于 n 的与 n 互质的正整数的数量。
参数:
n:整数
示例
>>> from sympy.functions.combinatorial.numbers import totient
>>> totient(1)
1
>>> totient(25)
20
>>> totient(45) == totient(5)*totient(9)
True
另见
divisor_count
参考文献
[R663]
en.wikipedia.org/wiki/Euler%27s_totient_function
[R664]
mathworld.wolfram.com/TotientFunction.html
sympy.ntheory.factor_.reduced_totient(n)
计算 Carmichael 减少的欧拉 totient 函数 lambda(n)
自 1.13 版本起不推荐使用:reduced_totient
函数已被弃用。请使用sympy.functions.combinatorial.numbers.reduced_totient
。更多信息请参见其文档。详细信息请参见将 ntheory 中的符号函数重定位到 functions。
reduced_totient(n)
或 (\lambda(n)) 是最小的 m > 0,使得对于所有与 n 互质的 k,都有 (k^m \equiv 1 \mod n)。
示例
>>> from sympy.functions.combinatorial.numbers import reduced_totient
>>> reduced_totient(1)
1
>>> reduced_totient(8)
2
>>> reduced_totient(30)
4
另见
totient
参考
[R665]
en.wikipedia.org/wiki/Carmichael_function
[R666]
mathworld.wolfram.com/CarmichaelFunction.html
sympy.ntheory.factor_.divisor_sigma(n, k=1)
计算正整数 n 的除数函数 (\sigma_k(n))
自版本 1.13 起不推荐使用:函数divisor_sigma
已不推荐使用。请使用sympy.functions.combinatorial.numbers.divisor_sigma
。有关详细信息,请参阅其文档。有关详细信息,请参见将符号函数从 ntheory 移至函数。
divisor_sigma(n, k)
等于 sum([x**k for x in divisors(n)])
如果 n 的素因数分解为:
[n = \prod_{i=1}^\omega p_i^{m_i},]
则
[\sigma_k(n) = \prod_{i=1}^\omega (1+p_ik+p_i+\cdots + p_i^{m_ik}).]
参数:
n:整数
k:整数,可选
和的除数的幂
对于 k = 0, 1:
divisor_sigma(n, 0)
等于divisor_count(n)
,divisor_sigma(n, 1)
等于sum(divisors(n))
默认情况下,k 为 1。
示例
>>> from sympy.functions.combinatorial.numbers import divisor_sigma
>>> divisor_sigma(18, 0)
6
>>> divisor_sigma(39, 1)
56
>>> divisor_sigma(12, 2)
210
>>> divisor_sigma(37)
38
另见
divisor_count
,totient
,divisors
,factorint
参考
[R667]
en.wikipedia.org/wiki/Divisor_function
sympy.ntheory.factor_.udivisor_sigma(n, k=1)
计算正整数 n 的单位除数函数 (\sigma_k^*(n))
自版本 1.13 起不推荐使用:函数udivisor_sigma
已不推荐使用。请使用sympy.functions.combinatorial.numbers.udivisor_sigma
。有关详细信息,请参阅其文档。有关详细信息,请参见将符号函数从 ntheory 移至函数。
udivisor_sigma(n, k)
等于 sum([x**k for x in udivisors(n)])
如果 n 的素因数分解为:
[n = \prod_{i=1}^\omega p_i^{m_i},]
则
[\sigma_k^*(n) = \prod_{i=1}^\omega (1+ p_i^{m_ik}).]
参数:
k:和的除数的幂
对于 k = 0, 1:
udivisor_sigma(n, 0)
等于udivisor_count(n)
,udivisor_sigma(n, 1)
等于sum(udivisors(n))
默认情况下,k 为 1。
示例
>>> from sympy.functions.combinatorial.numbers import udivisor_sigma
>>> udivisor_sigma(18, 0)
4
>>> udivisor_sigma(74, 1)
114
>>> udivisor_sigma(36, 3)
47450
>>> udivisor_sigma(111)
152
另见
divisor_count
, totient
, divisors
, udivisors
, udivisor_count
, divisor_sigma
, factorint
参考
[R668]
mathworld.wolfram.com/UnitaryDivisorFunction.html
sympy.ntheory.factor_.core(n, t=2)
计算正整数 n 的core_t(n)
= core(n, t)
core_2(n)
等于 n 的无平方数部分
如果 n 的质因数分解为:
[n = \prod_{i=1}^\omega p_i^{m_i},]
则
[core_t(n) = \prod_{i=1}^\omega p_i^{m_i \mod t}.]
参数:
n:整数
t:整数
core(n, t)
计算 n 的第 t 个无幂次方部分
core(n, 2)
是 n 的无平方数部分,core(n, 3)
是 n 的无立方数部分t 的默认值为 2。
例子
>>> from sympy.ntheory.factor_ import core
>>> core(24, 2)
6
>>> core(9424, 3)
1178
>>> core(379238)
379238
>>> core(15**11, 10)
15
参见
factorint
, sympy.solvers.diophantine.diophantine.square_factor
参考
[R669]
en.wikipedia.org/wiki/Square-free_integer#Squarefree_core
sympy.ntheory.factor_.digits(n, b=10, digits=None)
返回 n 在 b 进制下的数字列表。列表中的第一个元素是 b(如果 n 为负,则为-b)。
参数:
n:整数
要返回其数字的数字。
b:整数
计算数字的基数。
digits:整数(或所有数字为 None)
要返回的数字位数(如有必要,用零填充)。
例子
>>> from sympy.ntheory.digits import digits
>>> digits(35)
[10, 3, 5]
如果数字为负,则负号将放置在基数上(即返回列表中的第一个元素):
>>> digits(-35)
[-10, 3, 5]
可以选择除 10 以外(且大于 1)的基数b
:
>>> digits(27, b=2)
[2, 1, 1, 0, 1, 1]
如果需要确定特定位数的数字,请使用digits
关键字:
>>> digits(35, digits=4)
[10, 0, 0, 3, 5]
参见
sympy.core.intfunc.num_digits
, count_digits
sympy.ntheory.factor_.primenu(n)
计算正整数 n 的不同质因子的数量。
自版本 1.13 起已弃用:primenu
函数已弃用。请使用sympy.functions.combinatorial.numbers.primenu
。有关详细信息,请参见其文档。有关详细信息,请参阅从 ntheory 移动符号函数到函数。
如果 n 的素因数分解是:
[n = \prod_{i=1}^k p_i^{m_i},]
如果 primenu(n)
或者 (\nu(n)) 是:
[\nu(n) = k.]
例子
>>> from sympy.functions.combinatorial.numbers import primenu
>>> primenu(1)
0
>>> primenu(30)
3
参见
factorint
参考文献
[R670]
mathworld.wolfram.com/PrimeFactor.html
sympy.ntheory.factor_.primeomega(n)
计算一个正整数 n 的素因数的个数(包括重数)。
自版本 1.13 起已弃用:primeomega
函数已弃用。使用 sympy.functions.combinatorial.numbers.primeomega
替代。更多信息请参见其文档。详情请见 将符号函数从 ntheory 移至 functions。
如果 n 的素因数分解是:
[n = \prod_{i=1}^k p_i^{m_i},]
如果 primeomega(n)
或者 (\Omega(n)) 是:
[\Omega(n) = \sum_{i=1}^k m_i.]
例子
>>> from sympy.functions.combinatorial.numbers import primeomega
>>> primeomega(1)
0
>>> primeomega(20)
3
参见
factorint
参考文献
[R671]
mathworld.wolfram.com/PrimeFactor.html
sympy.ntheory.factor_.mersenne_prime_exponent(nth)
返回第 n 个 Mersenne prime 的指数 i
(形式为 (2^i - 1))。
例子
>>> from sympy.ntheory.factor_ import mersenne_prime_exponent
>>> mersenne_prime_exponent(1)
2
>>> mersenne_prime_exponent(20)
4423
sympy.ntheory.factor_.is_perfect(n)
返回 True 如果 n
是一个 perfect number,否则返回 False。
一个 perfect number 等于其正的真因数之和。
例子
>>> from sympy.functions.combinatorial.numbers import divisor_sigma
>>> from sympy.ntheory.factor_ import is_perfect, divisors
>>> is_perfect(20)
False
>>> is_perfect(6)
True
>>> 6 == divisor_sigma(6) - 6 == sum(divisors(6)[:-1])
True
参考文献
[R672]
mathworld.wolfram.com/PerfectNumber.html
[R673]
en.wikipedia.org/wiki/Perfect_number
sympy.ntheory.factor_.abundance(n)
返回一个数的正因数之和与该数之差。
例子
>>> from sympy.ntheory import abundance, is_perfect, is_abundant
>>> abundance(6)
0
>>> is_perfect(6)
True
>>> abundance(10)
-2
>>> is_abundant(10)
False
sympy.ntheory.factor_.is_abundant(n)
返回 True 如果 n
是一个 abundant number,否则返回 False。
一个 abundant number 是小于其正因数之和的数。
例子
>>> from sympy.ntheory.factor_ import is_abundant
>>> is_abundant(20)
True
>>> is_abundant(15)
False
参考文献
[R674]
mathworld.wolfram.com/AbundantNumber.html
sympy.ntheory.factor_.is_deficient(n)
返回 True 如果 n
是一个 deficient number,否则返回 False。
一个 deficient number 是大于其正因数之和的数。
例子
>>> from sympy.ntheory.factor_ import is_deficient
>>> is_deficient(20)
False
>>> is_deficient(15)
True
参考文献
[R675]
mathworld.wolfram.com/DeficientNumber.html
sympy.ntheory.factor_.is_amicable(m, n)
返回 True 如果数字 (m) 和 (n) 是“amicable”,否则返回 False。
Amicable numbers 是两个不同的相关数,其各自的真因数之和相等。
例子
>>> from sympy.functions.combinatorial.numbers import divisor_sigma
>>> from sympy.ntheory.factor_ import is_amicable
>>> is_amicable(220, 284)
True
>>> divisor_sigma(220) == divisor_sigma(284)
True
参考文献
[R676]
en.wikipedia.org/wiki/Amicable_numbers
sympy.ntheory.factor_.is_carmichael(n)
返回 True 如果数字 (n) 是 Carmichael number,否则返回 False。
参数:
n : 整数
参考文献
[R677]
en.wikipedia.org/wiki/Carmichael_number
[R678]
sympy.ntheory.factor_.find_carmichael_numbers_in_range(x, y)
返回在范围内的 Carmichael 数的列表。
另请参阅
is_carmichael
sympy.ntheory.factor_.find_first_n_carmichaels(n)
返回前 n
个 Carmichael 数。
参数:
n:整数
另请参阅
is_carmichael
sympy.ntheory.modular.symmetric_residue(a, m)
返回模 m
内的剩余部分,使其在模的一半范围内。
>>> from sympy.ntheory.modular import symmetric_residue
>>> symmetric_residue(1, 6)
1
>>> symmetric_residue(4, 6)
-2
sympy.ntheory.modular.crt(m, v, symmetric=False, check=True)
中国剩余定理。
在 m
中的模数被假设为两两互质。然后输出是一个整数 f
,使得对于 v
和 m
中的每对,都有 f = v_i mod m_i
。如果 symmetric
为 False,则返回一个正整数,否则 |f| 小于或等于模数的 LCM,因此 f
可能为负数。
如果模数不互质,则在发现结果不正确时将返回正确结果。如果没有解决方案,则结果将为 None
。
如果已知模数是互质的,则可以将关键字 check
设置为 False。
示例
例如,考虑一组剩余 U = [49, 76, 65]
和一组模数 M = [99, 97, 95]
。那么我们有:
>>> from sympy.ntheory.modular import crt
>>> crt([99, 97, 95], [49, 76, 65])
(639985, 912285)
这是正确的结果,因为:
>>> [639985 % m for m in [99, 97, 95]]
[49, 76, 65]
如果模数不互质,如果使用 check=False
,则可能会得到不正确的结果:
>>> crt([12, 6, 17], [3, 4, 2], check=False)
(954, 1224)
>>> [954 % m for m in [12, 6, 17]]
[6, 0, 2]
>>> crt([12, 6, 17], [3, 4, 2]) is None
True
>>> crt([3, 6], [2, 5])
(5, 6)
注意:相对于 crt
,gf_crt
的参数顺序被颠倒了,solve_congruence
接受剩余,模数对。
程序员注意:与检查所有模数对不共享 GCD(O(n**2) 测试)或者分解所有模数然后看是否有共同因子的方法不同,这里执行的操作是检查结果是否给出了指定的剩余部分 - 这是一个 O(n) 操作。
另请参阅
solve_congruence
sympy.polys.galoistools.gf_crt
低级 Chinese Remainder Theorem 运行在这个例程中使用。
sympy.ntheory.modular.crt1(m)
多个应用的 Chinese Remainder Theorem 的第一部分。
示例
>>> from sympy.ntheory.modular import crt, crt1, crt2
>>> m = [99, 97, 95]
>>> v = [49, 76, 65]
以下两个代码具有相同的结果。
>>> crt(m, v)
(639985, 912285)
>>> mm, e, s = crt1(m)
>>> crt2(m, v, mm, e, s)
(639985, 912285)
但是,当我们希望固定 m
并计算多个 v
时,速度更快,即以下情况:
>>> mm, e, s = crt1(m)
>>> vs = [[52, 21, 37], [19, 46, 76]]
>>> for v in vs:
... print(crt2(m, v, mm, e, s))
(397042, 912285)
(803206, 912285)
另请参阅
sympy.polys.galoistools.gf_crt1
低级 Chinese Remainder Theorem 运行在这个例程中使用。
sympy.ntheory.modular.crt
, sympy.ntheory.modular.crt2
sympy.ntheory.modular.crt2(m, v, mm, e, s, symmetric=False)
多个应用的 Chinese Remainder Theorem 的第二部分。
参见 crt1
以获取用法。
示例
>>> from sympy.ntheory.modular import crt1, crt2
>>> mm, e, s = crt1([18, 42, 6])
>>> crt2([18, 42, 6], [0, 0, 0], mm, e, s)
(0, 4536)
另请参阅
sympy.polys.galoistools.gf_crt2
低级 Chinese Remainder Theorem 运行在这个例程中使用。
sympy.ntheory.modular.crt
, sympy.ntheory.modular.crt1
sympy.ntheory.modular.solve_congruence(*remainder_modulus_pairs, **hint)
计算整数n
,当它被分别除以mi
时,余数为ai
,其中ai
和mi
作为此函数的一对给定:((a1, m1), (a2, m2), ...)。如果没有解决方案,则返回 None。否则返回n
及其模数。
mi
值不需要互质。如果已知模数不互质,则可以将提示check
设置为 False(默认为 True),并且将跳过通过 crt()(在模数互质时有效)寻找更快解决方案的检查。
如果提示symmetric
为 True(默认为 False),则n
的值将在模数的 1/2 范围内,可能为负。
示例
>>> from sympy.ntheory.modular import solve_congruence
数字是 2 模 3、3 模 5 和 2 模 7 的数字?
>>> solve_congruence((2, 3), (3, 5), (2, 7))
(23, 105)
>>> [23 % m for m in [3, 5, 7]]
[2, 3, 2]
如果您更喜欢将所有余数放在一个列表中,将所有模数放在另一个列表中,则可以像这样发送参数:
>>> solve_congruence(*zip((2, 3, 2), (3, 5, 7)))
(23, 105)
模数不需要互质;在这种情况下可能存在或不存在解决方案:
>>> solve_congruence((2, 3), (4, 6)) is None
True
>>> solve_congruence((2, 3), (5, 6))
(5, 6)
对称标志将使结果在模数的 1/2 范围内:
>>> solve_congruence((2, 3), (5, 6), symmetric=True)
(-1, 6)
另请参阅
crt
实现中国余数定理的高级程序
sympy.ntheory.multinomial.binomial_coefficients(n)
返回包含对{(k1,k2) : C_kn}
的字典,其中C_kn
为二项式系数且n=k1+k2
。
示例
>>> from sympy.ntheory import binomial_coefficients
>>> binomial_coefficients(9)
{(0, 9): 1, (1, 8): 9, (2, 7): 36, (3, 6): 84,
(4, 5): 126, (5, 4): 126, (6, 3): 84, (7, 2): 36, (8, 1): 9, (9, 0): 1}
另请参阅
二项式系数列表
,多项式系数
sympy.ntheory.multinomial.binomial_coefficients_list(n)
返回帕斯卡三角形的行作为二项式系数的列表。
示例
>>> from sympy.ntheory import binomial_coefficients_list
>>> binomial_coefficients_list(9)
[1, 9, 36, 84, 126, 126, 84, 36, 9, 1]
另请参阅
二项式系数
,多项式系数
sympy.ntheory.multinomial.multinomial_coefficients(m, n)
返回包含对{(k1,k2,..,km) : C_kn}
的字典,其中C_kn
为多项式系数,使得n=k1+k2+..+km
。
示例
>>> from sympy.ntheory import multinomial_coefficients
>>> multinomial_coefficients(2, 5) # indirect doctest
{(0, 5): 1, (1, 4): 5, (2, 3): 10, (3, 2): 10, (4, 1): 5, (5, 0): 1}
注释
该算法基于以下结果:
[\binom{n}{k_1, \ldots, k_m} = \frac{k_1 + 1}{n - k_1} \sum_{i=2}^m \binom{n}{k_1 + 1, \ldots, k_i - 1, \ldots}]
由 Yann Laigle-Chapuy 贡献的代码,经作者许可复制。
另请参阅
二项式系数列表
,二项式系数
sympy.ntheory.multinomial.multinomial_coefficients_iterator(m, n, _tuple=<class 'tuple'>)
多项式系数迭代器
通过利用从multinomial_coefficients(n, n)
获取的单项式元组t
去除零时的系数与后者的系数相同,因此,后者的系数预先计算以节省内存和时间。
>>> from sympy.ntheory.multinomial import multinomial_coefficients
>>> m53, m33 = multinomial_coefficients(5,3), multinomial_coefficients(3,3)
>>> m53[(0,0,0,1,2)] == m53[(0,0,1,0,2)] == m53[(1,0,2,0,0)] == m33[(0,1,2)]
True
示例
>>> from sympy.ntheory.multinomial import multinomial_coefficients_iterator
>>> it = multinomial_coefficients_iterator(20,3)
>>> next(it)
((3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), 1)
sympy.ntheory.partitions_.npartitions(n, verbose=False)
计算分割函数 P(n),即将 n 写为正整数和的方式数。
自版本 1.13 起弃用:npartitions
函数已弃用。请使用 sympy.functions.combinatorial.numbers.partition
。有关更多信息,请参阅其文档。详细信息请参见 将符号函数从 ntheory 移至 functions。
使用 Hardy-Ramanujan-Rademacher 公式计算 P(n) [R679]。
该实现的正确性已通过 (10^{10}) 的测试。
示例
>>> from sympy.functions.combinatorial.numbers import partition
>>> partition(25)
1958
参考文献
[R679] (1,2)
mathworld.wolfram.com/PartitionFunctionP.html
sympy.ntheory.primetest.is_fermat_pseudoprime(n, a)
如果 n
是素数或是与 a
互质且满足模算术同余关系的奇合数,则返回 True:
[a^{n-1} \equiv 1 \pmod{n}]
(其中 mod 指的是模运算)。
参数:
n:整数
n
是正整数。
a:整数
a
是正整数。a
和n
应该互质。
返回:
bool:如果 n
是素数,则总是返回 True
。
返回称为 Fermat 伪素数的合数。
示例
>>> from sympy.ntheory.primetest import is_fermat_pseudoprime
>>> from sympy.ntheory.factor_ import isprime
>>> for n in range(1, 1000):
... if is_fermat_pseudoprime(n, 2) and not isprime(n):
... print(n)
341
561
645
参考文献
[R680]
en.wikipedia.org/wiki/Fermat_pseudoprime
sympy.ntheory.primetest.is_euler_pseudoprime(n, a)
如果 n
是素数或是与 a
互质且满足模算术同余关系的奇合数,则返回 True:
[a^{(n-1)/2} \equiv \pm 1 \pmod{n}]
(其中 mod 指的是模运算)。
参数:
n:整数
n
是正整数。
a:整数
a
是正整数。a
和n
应该互质。
返回:
bool:如果 n
是素数,则总是返回 True
。
返回称为 Euler 伪素数的合数。
示例
>>> from sympy.ntheory.primetest import is_euler_pseudoprime
>>> from sympy.ntheory.factor_ import isprime
>>> for n in range(1, 1000):
... if is_euler_pseudoprime(n, 2) and not isprime(n):
... print(n)
341
561
参考文献
[R681]
en.wikipedia.org/wiki/Euler_pseudoprime
sympy.ntheory.primetest.is_euler_jacobi_pseudoprime(n, a)
如果 n
是素数或是与 a
互质且满足模算术同余关系的奇合数,则返回 True:
[a^{(n-1)/2} \equiv \left(\frac{a}{n}\right) \pmod{n}]
(其中 mod 指的是模运算)。
参数:
n:整数
n
是正整数。
a:整数
a
是正整数。a
和n
应该互质。
返回:
bool:如果 n
是素数,则总是返回 True
。
返回称为 Euler-Jacobi 伪素数的合数。
示例
>>> from sympy.ntheory.primetest import is_euler_jacobi_pseudoprime
>>> from sympy.ntheory.factor_ import isprime
>>> for n in range(1, 1000):
... if is_euler_jacobi_pseudoprime(n, 2) and not isprime(n):
... print(n)
561
参考文献
[R682]
en.wikipedia.org/wiki/Euler%E2%80%93Jacobi_pseudoprime
sympy.ntheory.primetest.is_square(n, prep=True)
如果 n == a * a 则返回 True,否则返回 False。如果怀疑 n 不 是一个平方数,则这是确认其不是的快速方法。
示例
>>> from sympy.ntheory.primetest import is_square
>>> is_square(25)
True
>>> is_square(2)
False
另请参阅
sympy.core.intfunc.isqrt
参考文献
[R683]
mersenneforum.org/showpost.php?p=110896
sympy.ntheory.primetest.mr(n, bases)
对 n 执行 Miller-Rabin 强伪素数测试,使用给定的基数/见证者列表。
示例
>>> from sympy.ntheory.primetest import mr
>>> mr(1373651, [2, 3])
False
>>> mr(479001599, [31, 73])
True
References
[R684]
Richard Crandall & Carl Pomerance (2005), “Prime Numbers: A Computational Perspective”, Springer, 2nd edition, 135-138
一个阈值列表及其所需的基数见这里:zh.wikipedia.org/wiki/Miller%E2%80%93Rabin%E8%B4%A8%E6%80%A7%E6%80%A7%E6%A3%80%E6%B5%8B#确定性变体
sympy.ntheory.primetest.is_lucas_prp(n)
标准 Lucas 复合性测试,使用 Selfridge 参数。如果 n 明确是复合数,则返回 False;如果 n 是 Lucas 可能质数,则返回 True。
这通常与米勒-拉宾测试结合使用。
示例
>>> from sympy.ntheory.primetest import isprime, is_lucas_prp
>>> for i in range(10000):
... if is_lucas_prp(i) and not isprime(i):
... print(i)
323
377
1159
1829
3827
5459
5777
9071
9179
References
[R685]
Robert Baillie, Samuel S. Wagstaff, Lucas Pseudoprimes, Math. Comp. Vol 35, Number 152 (1980), pp. 1391-1417, doi.org/10.1090%2FS0025-5718-1980-0583518-6
mpqs.free.fr/LucasPseudoprimes.pdf
[R686]
OEIS A217120: Lucas 伪素数 oeis.org/A217120
[R687]
zh.wikipedia.org/wiki/%E8%8E%B1%E5%8D%A1%E6%96%AF%E8%99%9A%E7%B4%A0%E5%88%86%E6%9E%90
sympy.ntheory.primetest.is_strong_lucas_prp(n)
强 Lucas 复合性测试,使用 Selfridge 参数。如果 n 明确是复合数,则返回 False;如果 n 是强 Lucas 可能质数,则返回 True。
这通常与米勒-拉宾测试结合使用,特别是与 M-R 基数 2 结合创建强 BPSW 测试时。
示例
>>> from sympy.ntheory.primetest import isprime, is_strong_lucas_prp
>>> for i in range(20000):
... if is_strong_lucas_prp(i) and not isprime(i):
... print(i)
5459
5777
10877
16109
18971
References
[R688]
Robert Baillie, Samuel S. Wagstaff, Lucas Pseudoprimes, Math. Comp. Vol 35, Number 152 (1980), pp. 1391-1417, doi.org/10.1090%2FS0025-5718-1980-0583518-6
mpqs.free.fr/LucasPseudoprimes.pdf
[R689]
OEIS A217255: 强 Lucas 伪素数 oeis.org/A217255
[R690]
zh.wikipedia.org/wiki/%E8%8E%B1%E5%8D%A1%E6%96%AF%E8%99%9A%E7%B4%A0%E5%88%86%E6%9E%90
[R691]
zh.wikipedia.org/wiki/Baillie-PSW%E8%B4%A8%E6%95%B0%E6%B5%8B%E8%AF%95
sympy.ntheory.primetest.is_extra_strong_lucas_prp(n)
额外强 Lucas 复合性测试。如果 n 明确是复合数,则返回 False;如果 n 是“额外强” Lucas 可能质数,则返回 True。
使用 P = 3,Q = 1 选择参数,然后递增 P 直到 (D|n) == -1. 测试本身如 [R692] 中所定义,出自 Mo 和 Jones 的预印本。参数选择和测试与 OEIS A217719、Perl 的 Math::Prime::Util 以及维基百科上的 Lucas 伪素数页面相同。
它比强测试快 20-50%。
由于选择了不同的参数,强 Lucas 伪素数和额外强 Lucas 伪素数之间没有关系。特别地,一个不是另一个的子集。
示例
>>> from sympy.ntheory.primetest import isprime, is_extra_strong_lucas_prp
>>> for i in range(20000):
... if is_extra_strong_lucas_prp(i) and not isprime(i):
... print(i)
989
3239
5777
10877
参考文献
[R692] (1,2)
Jon Grantham,Frobenius 伪素数,Math. Comp. Vol 70, Number 234 (2001),pp. 873-891,doi.org/10.1090%2FS0025-5718-00-01197-2
[R693]
OEIS A217719:额外强 Lucas 伪素数 oeis.org/A217719
[R694]
[zh.wikipedia.org/wiki/Lucas 伪素数
](https://zh.wikipedia.org/wiki/Lucas 伪素数)
sympy.ntheory.primetest.proth_test(n)
测试 Proth 数(n = k2^m + 1)是否为素数,其中 k 是正奇数且(2^m > k)。
参数:
n:整数
n
是 Proth 数
返回:
bool:如果为True
,则n
是 Proth 素数
抛出异常:
数值错误
如果
n
不是 Proth 数。
示例
>>> from sympy.ntheory.primetest import proth_test
>>> proth_test(41)
True
>>> proth_test(57)
False
参考文献
[R695]
[zh.wikipedia.org/wiki/Proth 素数
](https://zh.wikipedia.org/wiki/Proth 素数)
sympy.ntheory.primetest.is_mersenne_prime(n)
如果n
是梅森素数则返回 True,否则返回 False。
梅森素数是形如(2^i - 1)的素数。
示例
>>> from sympy.ntheory.factor_ import is_mersenne_prime
>>> is_mersenne_prime(6)
False
>>> is_mersenne_prime(127)
True
参考文献
[R696]
mathworld.wolfram.com/MersennePrime.html
sympy.ntheory.primetest.isprime(n)
测试 n 是否为素数(True)或不是素数(False)。对于 n < 2⁶⁴,答案是确定的;较大的 n 值可能实际上是伪素数。
负数(例如 -2)不被认为是素数。
第一步是寻找显著因子,如果找到则可以快速返回。接下来,如果筛选器足够大,则在筛选器上使用二分搜索。对于小数字,使用已知在其范围内没有反例的基数进行确定性 Miller-Rabin 测试。最后,如果数字大于 2⁶⁴,则执行强 BPSW 测试。虽然这是一个可能的素数测试,我们相信存在反例,但目前没有已知的反例。
示例
>>> from sympy.ntheory import isprime
>>> isprime(13)
True
>>> isprime(15)
False
注意
此程序仅适用于整数输入,不适用于可能代表数字的数值表达式。浮点数也会被拒绝,因为它们代表有限精度的数字。虽然允许 7.0 表示整数很诱人,但如果允许这样做可能会“悄悄地”传递错误:
>>> from sympy import Float, S
>>> int(1e3) == 1e3 == 10**3
True
>>> int(1e23) == 1e23
True
>>> int(1e23) == 10**23
False
>>> near_int = 1 + S(1)/10**19
>>> near_int == int(near_int)
False
>>> n = Float(near_int, 10) # truncated by precision
>>> n % 1 == 0
True
>>> n = Float(near_int, 20)
>>> n % 1 == 0
False
参见
sympy.ntheory.generate.primerange
在给定范围内生成所有素数
sympy.functions.combinatorial.numbers.primepi
返回小于或等于 n 的素数数量
sympy.ntheory.generate.prime
返回第 n 个素数
参考文献
[R697]
[R698]
Robert Baillie, Samuel S. Wagstaff, Lucas Pseudoprimes, Math. Comp. Vol 35, Number 152 (1980), pp. 1391-1417, doi.org/10.1090%2FS0025-5718-1980-0583518-6
mpqs.free.fr/LucasPseudoprimes.pdf
[R699]
en.wikipedia.org/wiki/Baillie-PSW_primality_test
sympy.ntheory.primetest.is_gaussian_prime(num)
测试 num 是否为高斯素数。
参考文献
[R700]
sympy.ntheory.residue_ntheory.n_order(a, n)
返回 a
模 n
的阶数。
参数:
a:整数
n:整数,n > 1。a
和 n
应该是相对素的。
返回:
整数:a
模 n
的阶数
引发:
数值错误
如果 (n \le 1) 或 (\gcd(a, n) \neq 1)。如果
a
或n
不是整数。
解释
a
模 n
的阶数是最小的整数 k
,使得 (a^k) 在 n
上留下余数。
示例
>>> from sympy.ntheory import n_order
>>> n_order(3, 7)
6
>>> n_order(4, 7)
3
另请参阅
is_primitive_root
当 a
模 n
的阶数等于 totient(n)
时,我们说 a
是 n
的原根。
sympy.ntheory.residue_ntheory.is_primitive_root(a, p)
如果 a
是 p
的原根,则返回 True。
参数:
a:整数
p:整数,p
> 1。a
和 p
应该是相对素的。
返回:
布尔:如果为 True,则 a
是 p
的原根。
引发:
数值错误
如果 (p \le 1) 或 (\gcd(a, p) \neq 1)。如果
a
或p
不是整数。
解释
如果 (\gcd(a, p) = 1) 并且 (\phi(p)) 是最小的正数,a
被称为 p
的原根。
(a^{\phi(p)} \equiv 1 \pmod{p})。
其中 (\phi(p)) 是欧拉函数。
只有对于 (p = 2, 4, q^e, 2q^e)(q
是奇素数),p
的原根才存在。因此,如果不是这样的 p
,则返回 False。要确定原根,我们需要知道 q-1
的素因子分解。这种确定的难度取决于这种复杂性。
示例
>>> from sympy.functions.combinatorial.numbers import totient
>>> from sympy.ntheory import is_primitive_root, n_order
>>> is_primitive_root(3, 10)
True
>>> is_primitive_root(9, 10)
False
>>> n_order(3, 10) == totient(10)
True
>>> n_order(9, 10) == totient(10)
False
另请参阅
primitive_root
sympy.ntheory.residue_ntheory.primitive_root(p, smallest=True)
返回 p
的一个原根或 None。
参数:
p:整数,p > 1
最小:如果为 True,则返回最小的原根,否则返回 None。
返回:
int | None:
如果存在原根,则返回
p
的原根。否则,返回 None。
引发:
数值错误
如果 (p \le 1) 或
p
不是整数。
解释
有关原根定义,请参阅 is_primitive_root
的解释。
对于p
的原始根仅适用于(p = 2, 4, q^e, 2qe)(`q`是奇质数)。现在,如果我们知道`q`的原始根,我们可以计算(qe)的原始根,如果我们知道(qe)的原始根,我们可以计算(2qe)的原始根。当不需要找到最小的原始根时,这个性质可以用来获取快速的原始根。另一方面,当我们想要最小的原始根时,我们会天真地确定是否是原始根。
示例
>>> from sympy.ntheory.residue_ntheory import primitive_root
>>> primitive_root(19)
2
>>> primitive_root(21) is None
True
>>> primitive_root(50, smallest=False)
27
另请参阅
is_primitive_root
引用
[R701]
- Stein“初等数论”(2011 年),第 44 页
[R702]
- Hackman“初等数论”(2009 年),第 C 章
sympy.ntheory.residue_ntheory.sqrt_mod(a, p, all_roots=False)
查找x**2 = a mod p
的根。
参数:
a :整数
p :正整数
all_roots :如果为 True 则返回根列表或 None
注:
如果没有根则返回 None;否则返回的根小于或等于p // 2
;一般来说不是最小的。仅当它是唯一的根时返回p // 2
。
仅在预期所有根都适合内存时使用all_roots
;否则使用sqrt_mod_iter
。
示例
>>> from sympy.ntheory import sqrt_mod
>>> sqrt_mod(11, 43)
21
>>> sqrt_mod(17, 32, True)
[7, 9, 23, 25]
sympy.ntheory.residue_ntheory.sqrt_mod_iter(a, p, domain=<class 'int'>)
迭代解决x**2 = a mod p
。
参数:
a :整数
p :正整数
domain :整数域,int
,ZZ
或 Integer
示例
>>> from sympy.ntheory.residue_ntheory import sqrt_mod_iter
>>> list(sqrt_mod_iter(11, 43))
[21, 22]
另请参阅
sqrt_mod
具有相同功能,但要求排序列表或仅有一个解。
sympy.ntheory.residue_ntheory.quadratic_residues(p) → list[int]
返回二次剩余的列表。
示例
>>> from sympy.ntheory.residue_ntheory import quadratic_residues
>>> quadratic_residues(7)
[0, 1, 2, 4]
sympy.ntheory.residue_ntheory.nthroot_mod(a, n, p, all_roots=False)
找到x**n = a mod p
的解决方案。
参数:
a :整数
n :正整数
p :正整数
all_roots :如果为 False,则返回最小的根,否则返回根列表
返回:
list[int] | int | None:
解决
x**n = a mod p
。输出类型表如下:
all_roots has roots Returns True Yes list[int] True No [] False Yes int False No None
引发:
ValueError
如果
a
,n
或p
不是整数。如果n
或p
不是正数。
示例
>>> from sympy.ntheory.residue_ntheory import nthroot_mod
>>> nthroot_mod(11, 4, 19)
8
>>> nthroot_mod(11, 4, 19, True)
[8, 11]
>>> nthroot_mod(68, 3, 109)
23
引用
[R703]
- Hackman“初等数论”(2009 年),第 76 页
sympy.ntheory.residue_ntheory.is_nthpow_residue(a, n, m)
如果x**n == a (mod m)
有解则返回 True。
引用
[R704]
- Hackman“初等数论”(2009 年),第 76 页
sympy.ntheory.residue_ntheory.is_quad_residue(a, p)
如果a
(mod p
)在p
的平方集中,则返回 True,即 a % p 在 set([i**2 % p for i in range(p)])中。
参数:
a :整数
p :正整数
返回:
bool :如果为 True,则x**2 == a (mod p)
有解。
引发:
ValueError
如果
a
,p
不是整数。如果p
不是正数。
示例
>>> from sympy.ntheory import is_quad_residue
>>> is_quad_residue(21, 100)
True
实际上,pow(39, 2, 100)
将是 21。
>>> is_quad_residue(21, 120)
False
即,对于任何整数x
,pow(x, 2, 120)
不是 21。
如果p
是奇质数,则使用迭代方法进行确定:
>>> from sympy.ntheory import is_quad_residue
>>> sorted(set([i**2 % 7 for i in range(7)]))
[0, 1, 2, 4]
>>> [j for j in range(7) if is_quad_residue(j, 7)]
[0, 1, 2, 4]
另请参阅
legendre_symbol
, jacobi_symbol
, sqrt_mod
sympy.ntheory.residue_ntheory.legendre_symbol(a, p)
返回勒让德符号 ((a / p))。
自版本 1.13 起已弃用:legendre_symbol
函数已弃用。请使用 sympy.functions.combinatorial.numbers.legendre_symbol
替代。更多信息请参阅其文档。参见 将符号函数从 ntheory 迁移到 functions 以获取详细信息。
对于整数 a
和奇素数 p
,勒让德符号定义为
[\begin{split}\genfrac(){}{}{a}{p} = \begin{cases} 0 & \text{如果 } p \text{ 整除 } a\ 1 & \text{如果 } a \text{ 是模 } p \text{ 的二次剩余}\ -1 & \text{如果 } a \text{ 是模 } p \text{ 的二次非剩余} \end{cases}\end{split}]
参数:
a : 整数
p : 奇素数
示例
>>> from sympy.functions.combinatorial.numbers import legendre_symbol
>>> [legendre_symbol(i, 7) for i in range(7)]
[0, 1, 1, -1, 1, -1, -1]
>>> sorted(set([i**2 % 7 for i in range(7)]))
[0, 1, 2, 4]
参见
is_quad_residue
, jacobi_symbol
sympy.ntheory.residue_ntheory.jacobi_symbol(m, n)
返回雅可比符号 ((m / n))。
自版本 1.13 起已弃用:jacobi_symbol
函数已弃用。请使用 sympy.functions.combinatorial.numbers.jacobi_symbol
替代。更多信息请参阅其文档。参见 将符号函数从 ntheory 迁移到 functions 以获取详细信息。
对于任意整数 m
和任意正奇整数 n
,雅可比符号定义为 n
的素因子对应的勒让德符号的乘积:
[\genfrac(){}{}{m}{n} = \genfrac(){}{}{m}{p{1}} \genfrac(){}{}{m}{p{2}} ... \genfrac(){}{}{m}{p{k}} \text{ 其中 } n = p_1^{\alpha_1} p_2^{\alpha_2} ... p_k^{\alpha_k}]
像勒让德符号一样,如果雅可比符号 (\genfrac(){}{}{m}{n} = -1),则 m
是模 n
的二次非剩余。
但与勒让德符号不同的是,如果雅可比符号 (\genfrac(){}{}{m}{n} = 1),则 m
可能是模 n
的二次剩余,也可能不是。
参数:
m : 整数
n : 奇正整数
示例
>>> from sympy.functions.combinatorial.numbers import jacobi_symbol, legendre_symbol
>>> from sympy import S
>>> jacobi_symbol(45, 77)
-1
>>> jacobi_symbol(60, 121)
1
jacobi_symbol
和 legendre_symbol
之间的关系可以如下所示:
>>> L = legendre_symbol
>>> S(45).factors()
{3: 2, 5: 1}
>>> jacobi_symbol(7, 45) == L(7, 3)**2 * L(7, 5)**1
True
参见
is_quad_residue
,legendre_symbol
sympy.ntheory.residue_ntheory.mobius(n)
Möbius 函数将自然数映射到{-1, 0, 1}
自版本 1.13 起不推荐使用mobius
函数。请使用sympy.functions.combinatorial.numbers.mobius
。有关详细信息,请参阅其文档。有关详细信息,请参阅从 ntheory 移动符号函数。
定义如下:
-
如果(n = 1),则为(1)。
-
如果(n)有一个平方素因子,则为(0)。
-
如果(n)是具有(k)个素因子的无平方正整数,则为((-1)^k)。
它是数论和组合数学中的重要乘性函数。在数学级数、代数数论以及物理学中(费米子算子与 Möbius 函数模型有着非常具体的实现)有应用。
参数:
n : 正整数
示例
>>> from sympy.functions.combinatorial.numbers import mobius
>>> mobius(13*7)
1
>>> mobius(1)
1
>>> mobius(13*7*5)
-1
>>> mobius(13**2)
0
参考文献
[R705]
en.wikipedia.org/wiki/M%C3%B6bius_function
[R706]
Thomas Koshy 的《Elementary Number Theory with Applications》
sympy.ntheory.residue_ntheory.discrete_log(n, a, b, order=None, prime_order=None)
计算a
以模n
为底的离散对数。
这是一个递归函数,用于将复合阶循环群中的离散对数问题简化为素数阶循环群中的问题。
它根据问题使用不同的算法(子群阶大小,是否为素数阶):
- 试乘法
- Baby-step giant-step 算法
- 波拉德 Rho 算法
- 指数演算法
- Pohlig-Hellman
示例
>>> from sympy.ntheory import discrete_log
>>> discrete_log(41, 15, 7)
3
参考文献
[R707]
mathworld.wolfram.com/DiscreteLogarithm.html
[R708]
“应用密码学手册”,Menezes, A. J., Van, O. P. C., & Vanstone, S. A. (1997).
sympy.ntheory.residue_ntheory.quadratic_congruence(a, b, c, n)
找到满足(a x² + b x + c \equiv 0 \pmod{n})的解。
参数:
a : int
b : int
c : int
n : int
正整数。
返回:
list[int]:
解的排序列表。如果没有解,则为
[]
。
示例
>>> from sympy.ntheory.residue_ntheory import quadratic_congruence
>>> quadratic_congruence(2, 5, 3, 7) # 2x² + 5x + 3 = 0 (mod 7)
[2, 6]
>>> quadratic_congruence(8, 6, 4, 15) # No solution
[]
参见
polynomial_congruence
解多项式同余方程
sympy.ntheory.residue_ntheory.polynomial_congruence(expr, m)
解多项式同余方程模 m。
参数:
expr : 整数系数多项式
m : 正整数
示例
>>> from sympy.ntheory import polynomial_congruence
>>> from sympy.abc import x
>>> expr = x**6 - 2*x**5 -35
>>> polynomial_congruence(expr, 6125)
[3257]
参见
sympy.polys.galoistools.gf_csolve
低级解决例程由该例程使用
sympy.ntheory.residue_ntheory.binomial_mod(n, m, k)
计算binomial(n, m) % k
。
参数:
n : 整数
m : 整数
k : 正整数
解释
使用 Lucas 定理的一般化与 Granville 的 R709(由中国剩余定理)提供的素数幂的余数计算时间为 O(log²(n) + q⁴log(n)log(p) + q⁴p*log³(p))。
例子
>>> from sympy.ntheory.residue_ntheory import binomial_mod
>>> binomial_mod(10, 2, 6) # binomial(10, 2) = 45
3
>>> binomial_mod(17, 9, 10) # binomial(17, 9) = 24310
0
参考文献
[R709] (1,2)
模素数幂的二项式系数,Andrew Granville,可用:web.archive.org/web/20170202003812/http://www.dms.umontreal.ca/~andrew/PDF/BinCoeff.pdf
sympy.ntheory.continued_fraction.continued_fraction(a) → list
返回有理数或二次无理数的连分数表示。
例子
>>> from sympy.ntheory.continued_fraction import continued_fraction
>>> from sympy import sqrt
>>> continued_fraction((1 + 2*sqrt(3))/5)
[0, 1, [8, 3, 34, 3]]
另见
continued_fraction_periodic
, continued_fraction_reduce
, continued_fraction_convergents
sympy.ntheory.continued_fraction.continued_fraction_convergents(cf)
返回一个连分数(CF)的收敛数迭代器。
参数应该是以下两种形式之一:- 部分商的列表,最后一个元素可能是重复的部分商列表,例如由 continued_fraction 和 continued_fraction_periodic 返回的形式。- 返回连分数的连续部分商的可迭代对象,例如由 continued_fraction_iterator 返回的形式。
在计算收敛分数时,连分数不必严格处于规范形式(全为整数,除了第一个为正数)。分数和负数元素可能存在于展开中。
例子
>>> from sympy.core import pi
>>> from sympy import S
>>> from sympy.ntheory.continued_fraction import continued_fraction_convergents, continued_fraction_iterator
>>> list(continued_fraction_convergents([0, 2, 1, 2]))
[0, 1/2, 1/3, 3/8]
>>> list(continued_fraction_convergents([1, S('1/2'), -7, S('1/4')]))
[1, 3, 19/5, 7]
>>> it = continued_fraction_convergents(continued_fraction_iterator(pi))
>>> for n in range(7):
... print(next(it))
3
22/7
333/106
355/113
103993/33102
104348/33215
208341/66317
>>> it = continued_fraction_convergents([1, [1, 2]]) # sqrt(3)
>>> for n in range(7):
... print(next(it))
1
2
5/3
7/4
19/11
26/15
71/41
另见
continued_fraction_iterator
, continued_fraction
, continued_fraction_periodic
sympy.ntheory.continued_fraction.continued_fraction_iterator(x)
返回 x 的连分数展开作为迭代器。
例子
>>> from sympy import Rational, pi
>>> from sympy.ntheory.continued_fraction import continued_fraction_iterator
>>> list(continued_fraction_iterator(Rational(3, 8)))
[0, 2, 1, 2]
>>> list(continued_fraction_iterator(Rational(-3, 8)))
[-1, 1, 1, 1, 2]
>>> for i, v in enumerate(continued_fraction_iterator(pi)):
... if i > 7:
... break
... print(v)
3
7
15
1
292
1
1
1
参考文献
[R710]
sympy.ntheory.continued_fraction.continued_fraction_periodic(p, q, d=0, s=1) → list
找到二次无理数的周期连分数展开。
计算有理数或二次无理数的连分数展开,即 (\frac{p + s\sqrt{d}}{q}),其中 (p),(q \ne 0),(d \ge 0) 为整数。
作为整数列表返回连分数表示(规范形式),对于二次无理数,可选地以表示重复数字的整数列表结束。
参数:
p : int
数字的有理数部分
q:整数
数字的分母
d:整数,可选
数字的分母
s:整数,可选
无理部分的系数
示例
>>> from sympy.ntheory.continued_fraction import continued_fraction_periodic
>>> continued_fraction_periodic(3, 2, 7)
[2, [1, 4, 1, 1]]
黄金分割率具有最简单的连分数展开:
>>> continued_fraction_periodic(1, 2, 5)
[[1]]
如果鉴别器为零或完全平方数,则该数字将是有理数:
>>> continued_fraction_periodic(4, 3, 0)
[1, 3]
>>> continued_fraction_periodic(4, 3, 49)
[3, 1, 2]
另见
continued_fraction_iterator
, continued_fraction_reduce
参考文献
[R711]
[R712]
K. Rosen. Elementary Number theory and its applications. Addison-Wesley, 3 Sub edition, pages 379-381, January 1992.
sympy.ntheory.continued_fraction.continued_fraction_reduce(cf)
将一个连分数化简为有理数或二次无理数。
从其终止或周期连分数展开计算有理或二次无理数。连分数展开(cf)应被视为提供展开项的终止迭代器。对于终止连分数,这相当于list(continued_fraction_convergents(cf))[-1]
,只是效率稍高。如果展开部分重复,迭代器应返回重复项的列表作为最后一个元素。这是由 continued_fraction_periodic 返回的格式。
对于二次无理数,如果分数处于规范形式(除第一个外所有项均为正),则返回找到的最大解,通常是所需解。
示例
>>> from sympy.ntheory.continued_fraction import continued_fraction_reduce
>>> continued_fraction_reduce([1, 2, 3, 4, 5])
225/157
>>> continued_fraction_reduce([-2, 1, 9, 7, 1, 2])
-256/233
>>> continued_fraction_reduce([2, 1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8]).n(10)
2.718281835
>>> continued_fraction_reduce([1, 4, 2, [3, 1]])
(sqrt(21) + 287)/238
>>> continued_fraction_reduce([[1]])
(1 + sqrt(5))/2
>>> from sympy.ntheory.continued_fraction import continued_fraction_periodic
>>> continued_fraction_reduce(continued_fraction_periodic(8, 5, 13))
(sqrt(13) + 8)/5
另见
continued_fraction_periodic
sympy.ntheory.digits.count_digits(n, b=10)
返回一个字典,其键是数字n
在给定基数b
中的数字,键表示数字出现在数字中的位置,值表示该数字出现的次数。
示例
>>> from sympy.ntheory import count_digits
>>> count_digits(1111339)
{1: 4, 3: 2, 9: 1}
返回的数字始终以十进制表示,但数字本身可以以 Python 理解的任何格式输入;如果数字的基数与 10 不同,也可以给出数字的基数:
>>> n = 0xFA; n
250
>>> count_digits(_)
{0: 1, 2: 1, 5: 1}
>>> count_digits(n, 16)
{10: 1, 15: 1}
默认情况下,字典将对任何未出现在数字中的数字返回 0。例如,77!
中出现 7 次的数字有哪些:
>>> from sympy import factorial
>>> c77 = count_digits(factorial(77))
>>> [i for i in range(10) if c77[i] == 7]
[1, 3, 7, 9]
另见
sympy.core.intfunc.num_digits
, digits
sympy.ntheory.digits.digits(n, b=10, digits=None)
返回n
在基数b
中的数字列表。列表中的第一个元素是b
(如果n
为负,则为-b
)。
参数:
n:整数
返回数字的位数。
b:整数
计算数字的基数。
digits:整数(或所有数字的 None)
要返回的数字位数(如有必要,用零填充)。
示例
>>> from sympy.ntheory.digits import digits
>>> digits(35)
[10, 3, 5]
如果数字为负,则负号将放在基数上(即返回列表中的第一个元素):
>>> digits(-35)
[-10, 3, 5]
可以选择除了 10(大于 1 的)以外的基数 b
:
>>> digits(27, b=2)
[2, 1, 1, 0, 1, 1]
如果需要某个特定数量的数字,请使用 digits
关键字:
>>> digits(35, digits=4)
[10, 0, 0, 3, 5]
另请参阅
sympy.core.intfunc.num_digits
, count_digits
sympy.ntheory.digits.is_palindromic(n, b=10)
如果在给定的基数 b
中从左到右或从右到左读取 n
时相同,则返回 True。
示例
>>> from sympy.ntheory import is_palindromic
>>> all(is_palindromic(i) for i in (-11, 1, 22, 121))
True
第二个参数允许您在其他基数中测试数字。例如,88 在十进制中是回文的,但在八进制中不是:
>>> is_palindromic(88, 8)
False
另一方面,一个数字在八进制中可能是回文的,但在十进制中不是:
>>> 0o121, is_palindromic(0o121)
(81, False)
或者在两个基础上是回文的:
>>> oct(121), is_palindromic(121, 8) and is_palindromic(121)
('0o171', True)
sympy.ntheory.egyptian_fraction.egyptian_fraction(r, algorithm='Greedy')
返回所述有理数 (r) 的埃及分数展开的分母列表[R713]。
参数:
r:有理数或(p,q)
正有理数
p/q
。
algorithm:{“贪婪”,“Graham Jewett”,“Takenouchi”,“Golomb”},可选
指定要使用的算法(默认为“贪婪”)。
示例
>>> from sympy import Rational
>>> from sympy.ntheory.egyptian_fraction import egyptian_fraction
>>> egyptian_fraction(Rational(3, 7))
[3, 11, 231]
>>> egyptian_fraction((3, 7), "Graham Jewett")
[7, 8, 9, 56, 57, 72, 3192]
>>> egyptian_fraction((3, 7), "Takenouchi")
[4, 7, 28]
>>> egyptian_fraction((3, 7), "Golomb")
[3, 15, 35]
>>> egyptian_fraction((11, 5), "Golomb")
[1, 2, 3, 4, 9, 234, 1118, 2580]
注意事项
目前支持以下算法:
-
贪婪算法
也称为 Fibonacci-Sylvester 算法[R714]。在每一步中,提取小于目标的最大单位分数,并用余数替换目标。
它具有一些独特的属性:
-
给定最简分数 (p/q),生成最大长度 (p) 的展开。即使分子变大,项数也很少超过一把。
-
使用最小内存。
-
术语可以爆炸(其标准示例为 5/121 和 31/311)。每一步中的分母最多是平方的(双指数增长),通常表现出单指数增长。
-
-
Graham Jewett 算法
根据 Graham 和 Jewett 的结果建议的算法。请注意,这倾向于爆炸:结果展开的长度始终为
2**(x/gcd(x, y)) - 1
。参见[R715]。 -
Takenouchi 算法
由 Takenouchi(1921)建议的算法。与 Graham-Jewett 算法仅在处理重复项时有所不同。参见[R715]。
-
Golomb 算法
由 Golumb(1962)提出的方法,使用模算术和逆元。它与 Bleicher(1972)提出的使用连分数的方法产生相同的结果。参见[R716]。
如果给定的有理数大于或等于 1,则使用贪婪算法对和谐序列 1/1 + 1/2 + 1/3 + …进行求和,直到添加一个以上的单位分数会比给定的数字大。这些分母的列表被添加到从剩余部分使用的请求算法的结果前面。例如,如果 r 是 8/3,使用贪婪算法,我们得到[1, 2, 3, 4, 5, 6, 7, 14, 420],其中序列的开头[1, 2, 3, 4, 5, 6, 7]是和谐序列的一部分,总和为 363/140,剩余 31/420,贪婪算法产生[14, 420]。埃及分数(Rational(8, 3),“Golomb”)的结果是[1, 2, 3, 4, 5, 6, 7, 14, 574, 2788, 6460, 11590, 33062, 113820],等等。
另见
sympy.core.numbers.Rational
参考文献
[R713] (1,2)
en.wikipedia.org/wiki/Egyptian_fraction
[R714] (1,2)
en.wikipedia.org/wiki/Greedy_algorithm_for_Egyptian_fractions
[R715] (1,2,3)
www.ics.uci.edu/~eppstein/numth/egypt/conflict.html
[R716] (1,2)
sympy.ntheory.bbp_pi.pi_hex_digits(n, prec=14)
返回一个包含prec
(默认为 14)个十六进制π的第 n 位开始的数字的字符串。数字的计数从 0 开始,小数点不计入,所以当 n = 0 时,返回的值以 3 开头;n = 1 对应于小数点后的第一个数字(在十六进制中为 2)。
参数:
n:非负整数
prec:非负整数,默认为 14
返回:
str:返回一个包含prec
个数字的字符串
从π的第 n 位开始的十六进制。如果
prec
= 0,则返回空字符串。
提升:
ValueError
如果
n
< 0 或prec
< 0。或n
或prec
不是整数。
示例
>>> from sympy.ntheory.bbp_pi import pi_hex_digits
>>> pi_hex_digits(0)
'3243f6a8885a30'
>>> pi_hex_digits(0, 3)
'324'
这些与以下结果一致
>>> import math
>>> hex(int(math.pi * 2**((14-1)*4)))
'0x3243f6a8885a30'
>>> hex(int(math.pi * 2**((3-1)*4)))
'0x324'
参考文献
[R717]
www.numberworld.org/digits/Pi/
ECM 函数
ecm
函数是一个次指数因数分解算法,能够在几秒钟内轻松地因式分解大约 35 位数的数字。ecm
的时间复杂度取决于该数的最小真因数。因此,即使数字很大,但其因子相对较小,ecm
也可以轻松地对其进行因式分解。例如,我们取带有 15 位数因子的 N,如 15154262241479,15423094826093,799333555511111,809709509409109,888888877777777,914148152112161。现在 N 是一个 87 位数。ECM
需要大约 47 秒来因式分解它。
sympy.ntheory.ecm.ecm(n, B1=10000, B2=100000, max_curve=200, seed=1234)
使用 Lenstra 的椭圆曲线方法进行因式分解。
此函数重复调用 _ecm_one_factor
来计算 n 的因子。首先通过试除法去除所有小因子。然后使用 _ecm_one_factor
逐个计算一个因子。
参数:
n : 待分解的数
B1 : 第 1 阶段界限。必须是偶数。
B2 : 第 2 阶段界限。必须是偶数。
max_curve : 生成的最大曲线数
seed : 初始化伪随机生成器
示例
>>> from sympy.ntheory import ecm
>>> ecm(25645121643901801)
{5394769, 4753701529}
>>> ecm(9804659461513846513)
{4641991, 2112166839943}
示例
>>> from sympy.ntheory import ecm
>>> ecm(7060005655815754299976961394452809, B1=100000, B2=1000000)
{6988699669998001, 1010203040506070809}
>>> ecm(122921448543883967430908091422761898618349713604256384403202282756086473494959648313841, B1=100000, B2=1000000)
{15154262241479,
15423094826093,
799333555511111,
809709509409109,
888888877777777,
914148152112161}
``` ## QS 函数
qs 函数是一种亚指数复杂度的因式分解算法,适用于 100 位数以内的最快因式分解算法。qs 的时间复杂度取决于数的大小,因此如果数包含大因子,则使用 qs。因此,在因式分解数时,首先使用 ecm 获得约 15 位数的较小因子,然后使用 qs 获取较大因子。
用于因式分解 \(2709077133180915240135586837960864768806330782747\) 的半素数,其具有两个 25 位数的因子。qs 能在大约 248 秒内完成对其的因式分解。
```py
sympy.ntheory.qs.qs(N, prime_bound, M, ERROR_TERM=25, seed=1234)
使用自初始化二次筛法进行因式分解。在 SIQS 中,令 N 为待分解的数,且 N 不应为完全幂。如果我们找到两个整数使得 X**2 = Y**2 modN
且 X != +-Y modN
,那么 (gcd(X + Y, N)) 将显示出 N 的一个正确因子。为了找到这些整数 X 和 Y,我们试图找到形如 t**2 = u modN 的关系,其中 u 是小质数的乘积。如果我们有足够多的这些关系,我们可以形成 (t1*t2...ti)**2 = u1*u2...ui modN
,使右边是一个平方数,因此我们找到了一个关系 X**2 = Y**2 modN
。
这里进行了几项优化,如使用多项式筛选、快速切换多项式和使用部分关系。部分关系的使用可以使因式分解加速 2 倍。
参数:
N : 待分解的数
prime_bound : 质因数基中的质数上限
M : 筛选区间
ERROR_TERM : 检查平滑性的误差项
threshold : 因式分解的额外平滑关系
seed : 生成伪质数
示例
>>> from sympy.ntheory import qs
>>> qs(25645121643901801, 2000, 10000)
{5394769, 4753701529}
>>> qs(9804659461513846513, 2000, 10000)
{4641991, 2112166839943}
参考文献
[R718]
pdfs.semanticscholar.org/5c52/8a975c1405bd35c65993abf5a4edb667c1db.pdf
[R719]
www.rieselprime.de/ziki/Self-initializing_quadratic_sieve
示例
>>> from sympy.ntheory import qs
>>> qs(5915587277*3267000013, 1000, 10000)
{3267000013, 5915587277}