Lev Brovchenkov user avatar Lev Brovchenkov · Updated May. 26, 22 · Tutorial
翻译/Adso2004
我们将在本文中介绍模拟N维空间中矢量场（例如电磁场）的方法；
本文主要分为四个部分：
- 定义理论基础（如numpy中的arrays）
- 定义粒子间通过场相互作用的机制
- 可视化电场
- 可视化例子在电磁场中的运动
理论基础
向量
描述任何物理场的基本元素都是向量。所以，我们需要了解什么？向量的距离、模的运算方法。我们将list参数传入向量，这是它的初始化方法：
def __init__(self, *el):
for e in el:
self.append(e)
现在，我们可以创建一个向量：
v = Vector(1, 2, 3)
让我们定义向量的加法运算：
def __add__(self, other):
if type(other) is Vector:
assert len(self) == len(other), "Error 0"
r = Vector()
for i in range(len(self)):
r.append(self[i]) + other[i])
return r
else:
other = Vector.emptyvec(lens=len(self), n=other)
return self + other
结果：
v2 = [2, 57, 23.2]
v1 + v2
>>> [3, 59, 26.2]
类比于上面的方法，我们可以定义向量所有的运算（文末附上完整的代码）
现在，我们需要一个距离函数。我可以简单地定义一个dist(v1, v2)——但这并不美观，我们将重新定义一个%操作：
...
def __mod__(self, other):
return sum((self - other)**2) ** 0.5
结果：
v2 = Vector(2, 57, 23.2)
v1 % v2
>>> 58.60068258988115
为了在提高运算效率的同时获得更加美观的输出，我们需要用到一系列方法，这里没有什么疑难之处，以上就是Vector类的全部代码。
粒子
这儿，从理论上来说，一切都简单易懂——粒子的坐标、速度以及加速度。此外，它还有质量和其他自定义参数（比如说，对于电磁场，我们可以设置它的电量）
粒子的初始化如下：
def __init__(self, coords, mass=1.0, q=1.0 speed=None, **properties):
self.coords = coords
if speed is None:
self.speed = Vector(*[0 for i in range(len(coords))])
else:
self.speed = speed
self.acc = Vector(*[0 for i in range(len(coords))])
self.mass = mass
self.__params__ = ["coords", "speed", "acc", "q"] + list(properties.keys())
self.q = q
for prop in properties:
setattr(self, prop, properties[prop])
下面的代码可以实现粒子的移动、固定和加速：
...
def move(self, dt):
self.coords = self.coords + self.speed * dt
def accelerate(self, dt):
self.speed = self.speed + self.acc * dt
def accinc(self, force): # Considering exerted force the point obtains acceleration
self.acc = self.acc + force / self.mass
def clean_acc(self):
self.acc = self.acc * 0
这样粒子的设置基本完成了。
交互物理场
场是一种包含了空间中所有粒子及其相互间作用力的物质。我们将其视为一种特例，因此我们将自定义一个场（当然，它很容易推广）。 首先声明构造函数并添加一个点：
def __init__(self, F): # F - is a custom force, F(p1, p2, r), p1, p2 - points, r - distance inbetween
self.points = []
self.F = F
def append(self, *args, **kwargs):
self.points.append(Point(*args, **kwargs))
现在，我们需要声明一个返回粒子“张力”的函数。尽管这个概念是从电磁场中引申出的，不过在我们的例子当中，场被抽象为向量，我们将沿着这些向量移动粒子。在这种情况下，我们将得到C粒子的属性，而在特殊情境下可以是点电荷的电量（实际上无论是什么都行，甚至是向量）。 那么，什么是张力（tension）？是形如：
$$
\vec{E}(\vec{C})=\sum\vec{F_i}(\vec{C})
$$
C点的张力
C粒子处的场强等于作用于该点的所有力的矢量和。
...
def intensity(self, coord):
proj = Vector(*[0 for i in range(coord.dim())])
single_point = Point(Vector(), mass=1.0, q=1.0) # That's our "Single point"
for p in self.points:
if coord % p.coords < 10 ** (-10): # Check whether we compare coord with a point P where P.coords = coord
continue
d = p.coords % coord
fmod = self.F(single_point, p, d) * (-1)
proj = proj + (coord - p.coords) / d * fmod
return proj
这样，我们已经可以将向量场可视化了，但我们还需要一点收尾工作：
...
def step(self, dt):
self.clean_acc()
for p in self.points:
p.accinc(self.intensity(p.coords) * p.q)
p.accelerate(dt)
p.move(dt)
对于每个点，我们先确定其坐标，再确定其处场强及作用力：
$$
\vec{F}=\vec{E}*q
$$
粒子运动及矢量场的可视化
终于到了最有趣的部分，让我们开始吧…
模拟粒子在电磁场中的运动
for i in range(3):
u.append(Vector.randvec(2) * 10, q=random.random() - 0.5)
事实上，系数k的值应当是$9*10^{-9}$，然而由于它会在$t \rightarrow 0$时衰灭，从而更容易使它们成为正态数，因此在物理学中定义k=30000。
接下来，我们在每条轴0到10的位置上添加十个点（二维空间）。同时，我们在每个点上设置一个电量从-0.25到0.25的电荷。然后，我们设计一个循环来通过它们的坐标（和轨迹）绘制图像：
for i in range(130):
u.step(0.0006)
xd, yd = zip(*u.gather_coords())
X.extend(xd)
Y.extend(yd)
plt.figure(figsize=[8, 8])
plt.scatter(X, Y)
plt.scatter(*zip(*u.gather_coords()), color="orange")
plt.show()
译者注
- init函数
在init()这个初始化实例对象的方法中，self就是实例对象；而所谓初始化，例如self.name = name和self.age = age，就是给self绑定name, age两个属性，便于我们之后给这两个属性赋值。
- assert语句
- 原作者虽说给出完整代码，但不难发现这些代码已经无迹可寻，请读者自行完成。
0 条评论