向量积分的应用#

要在区域上积分标量场或向量场，我们必须首先定义一个区域。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^2\) . 求总质量。

>>> triangle = ParametricRegion((x, y), (x, 0, 5), (y, 0, 5 - x))
>>> vector_integrate(C.x*C.y, triangle)
625/24

求一个圆柱体的密度，它的半径是a轴的中心 \(\rho = x^2 + y^2\:kg/m^2\) .

>>> 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^2 \mathbf{{\hat{{k}}}}\) 在x-y平面上，还有一个场 \(E(x, y, z) = y^2 \mathbf{{\hat{{k}}}}\) 在x-y平面以下。这个向量场通过一个边长为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

斯托克定理的验证#

看到了吗https://en.wikipedia.org/wiki/Stokes%27_定理

例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

验证散度定理#

看到了吗https://en.wikipedia.org/wiki/difference_定理

例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)