在NumPy中确定静态平衡#

在分析物理结构时,了解使其保持稳定的力学原理至关重要。施加在地板、横梁或任何其他结构上的力会产生反作用力和力矩。这些反作用力是结构在不破坏的情况下抵抗运动的表现。在结构尽管受到力作用但仍不移动的情况下,牛顿第二定律指出,系统中的加速度和所有方向上的力之和都必须为零。您可以使用NumPy数组来表示和解决这个概念。

您将做什么:#

  • 在本教程中,您将使用NumPy数组创建向量和力矩。

  • 解决涉及电缆和地板支撑结构的问题。

  • 编写NumPy矩阵来隔离未知数。

  • 使用NumPy函数执行线性代数运算。

您将学到什么:#

  • 如何用NumPy表示点、向量和力矩。

  • 如何找到向量的法线

  • 使用NumPy进行矩阵计算。

您将需要什么:#

用以下命令导入

import numpy as np
import matplotlib.pyplot as plt

在本教程中,您将使用以下NumPy工具

  • np.linalg.norm:此函数确定向量的量度大小。

  • np.cross:此函数接受两个矩阵并产生叉积。

用牛顿第二定律求解平衡#

您的模型由一根受到一系列力和力矩作用的横梁组成。您可以从牛顿第二定律开始分析这个系统。

\[\sum{\text{force}} = \text{mass} \times \text{acceleration}.\]

为了简化所考察的例子,假设它们是静止的,加速度\(=0\)。由于我们的系统存在于三维空间中,因此考虑在每个维度上施加的力。这意味着您可以将这些力表示为向量。对于力矩,您得出相同的结论,力矩是由于力作用于离物体质心一定距离处而产生的。

假设力\(F\)表示为三维向量

\[F = (F_x, F_y, F_z)\]

其中三个分量分别代表在相应方向上施加的力的大小。还假设向量中的每个分量

\[r = (r_x, r_y, r_z)\]

是每个力分量施加点与系统质心之间的距离。那么,力矩可以计算为

\[r \times F = (r_x, r_y, r_z) \times (F_x, F_y, F_z).\]

从一些简单的力向量示例开始

forceA = np.array([1, 0, 0])
forceB = np.array([0, 1, 0])
print("Force A =", forceA)
print("Force B =", forceB)
Force A = [1 0 0]
Force B = [0 1 0]

这定义了forceA为在\(x\)方向上大小为1的向量,forceB为在\(y\)方向上大小为1的向量。

可视化这些力可能有助于更好地理解它们之间的相互作用。Matplotlib是一个包含可视化工具的库,可用于此目的。Quiver图将用于演示三维向量,但该库也可用于二维演示

fig = plt.figure()

d3 = fig.add_subplot(projection="3d")

d3.set_xlim(-1, 1)
d3.set_ylim(-1, 1)
d3.set_zlim(-1, 1)

x, y, z = np.array([0, 0, 0])  # defining the point of application.  Make it the origin

u, v, w = forceA  # breaking the force vector into individual components
d3.quiver(x, y, z, u, v, w, color="r", label="forceA")

u, v, w = forceB
d3.quiver(x, y, z, u, v, w, color="b", label="forceB")

plt.legend()
plt.show()
../_images/554b0d786ab48f9201684a1150327fd10ec8192b96061247607cb2196bc80240.png

有两个力从一个点发出。为了简化这个问题,您可以将它们相加,以找到力的总和。请注意,forceAforceB都是三维向量,由NumPy表示为具有三个分量的数组。由于NumPy旨在简化和优化向量之间的操作,您可以轻松地计算这两个向量的总和,如下所示

forceC = forceA + forceB
print("Force C =", forceC)
Force C = [1 1 0]

力C现在充当一个单一的力,代表A和B。您可以绘制它来查看结果。

fig = plt.figure()

d3 = fig.add_subplot(projection="3d")

d3.set_xlim(-1, 1)
d3.set_ylim(-1, 1)
d3.set_zlim(-1, 1)

x, y, z = np.array([0, 0, 0])

u, v, w = forceA
d3.quiver(x, y, z, u, v, w, color="r", label="forceA")
u, v, w = forceB
d3.quiver(x, y, z, u, v, w, color="b", label="forceB")
u, v, w = forceC
d3.quiver(x, y, z, u, v, w, color="g", label="forceC")

plt.legend()
plt.show()
../_images/57f25a6b4b08d491594ade0e955f5119fd3a7f44b1bb828bd944c7707597ba87.png

然而,目标是平衡。这意味着您希望您的力的总和为\((0, 0, 0)\),否则您的物体将经历加速度。因此,需要有另一个力来抵消之前的力。

您可以将这个问题写为\(A+B+R=0\),其中\(R\)是解决问题的反作用力。

在这个例子中,这意味着

\[(1, 0, 0) + (0, 1, 0) + (R_x, R_y, R_z) = (0, 0, 0)\]

分解为\(x\)\(y\)\(z\)分量,这给您

\[\begin{split}\begin{cases} 1+0+R_x=0\\ 0+1+R_y=0\\ 0+0+R_z=0 \end{cases}\end{split}\]

求解\(R_x\)\(R_y\)\(R_z\),得到向量\(R\)\((-1, -1, 0)\)

如果绘制出来,先前示例中看到的力应该被抵消。只有当没有力剩余时,系统才被认为是处于平衡状态。

R = np.array([-1, -1, 0])

fig = plt.figure()

d3.set_xlim(-1, 1)
d3.set_ylim(-1, 1)
d3.set_zlim(-1, 1)

d3 = fig.add_subplot(projection="3d")

x, y, z = np.array([0, 0, 0])

u, v, w = forceA + forceB + R  # add them all together for sum of forces
d3.quiver(x, y, z, u, v, w)

plt.show()
../_images/6cbc62722ac4f8219c2fac75d2d38048d171b778c939f4ce5dbe0537bc307770.png

空白图表示没有外部力。这表示系统处于平衡状态。

将力矩和求解平衡#

接下来,让我们转向一个更复杂的应用。当力不都作用在同一点时,会产生力矩。

与力类似,这些力矩的总和也必须为零,否则将发生旋转加速度。与力的总和类似,这为空间中的三个坐标方向中的每一个都创建了一个线性方程。

一个简单的例子是力作用在一个固定在地面的静止杆上。杆不移动,所以它必须施加一个反作用力。杆也不旋转,所以它也必须产生一个反作用力矩。求解反作用力和力矩。

假设一个5N的力垂直作用在杆基座上方2m处。

f = 5  # Force in newtons
L = 2  # Length of the pole

R = 0 - f
M = 0 - f * L
print("Reaction force =", R)
print("Reaction moment =", M)
Reaction force = -5
Reaction moment = -10

用物理性质寻找值#

假设不是力垂直作用于横梁,而是通过一根也连接到地面的电线施加到我们的杆上。给定这根电线的张力,您只需要这些物体的物理位置就可以解决这个问题。

Image representing the problem

为了响应作用在杆上的力,基座在x和y方向上产生了反作用力,以及一个反作用力矩。

将杆的基座标记为原点。现在,假设电线在x方向上连接到地面3m处,并在z方向上连接到杆2m处。

将这些空间点定义为NumPy数组,然后使用这些数组来找到方向向量。

poleBase = np.array([0, 0, 0])
cordBase = np.array([3, 0, 0])
cordConnection = np.array([0, 0, 2])

poleDirection = cordConnection - poleBase
print("Pole direction =", poleDirection)
cordDirection = cordBase - cordConnection
print("Cord direction =", cordDirection)
Pole direction = [0 0 2]
Cord direction = [ 3  0 -2]

为了在与力相关的方面使用这些向量,您需要将它们转换为单位向量。单位向量的大小为一,只表示力的方向。

cordUnit = cordDirection / np.linalg.norm(cordDirection)
print("Cord unit vector =", cordUnit)
Cord unit vector = [ 0.83205029  0.         -0.5547002 ]

然后,您可以将此方向与力的大小相乘,以找到力向量。

假设电线的张力为5N

cordTension = 5
forceCord = cordUnit * cordTension
print("Force from the cord =", forceCord)
Force from the cord = [ 4.16025147  0.         -2.77350098]

为了找到力矩,您需要力向量和半径的叉积。

momentCord = np.cross(forceCord, poleDirection)
print("Moment from the cord =", momentCord)
Moment from the cord = [ 0.         -8.32050294  0.        ]

现在您只需找到反作用力和力矩。

equilibrium = np.array([0, 0, 0])
R = equilibrium - forceCord
M = equilibrium - momentCord
print("Reaction force =", R)
print("Reaction moment =", M)
Reaction force = [-4.16025147  0.          2.77350098]
Reaction moment = [0.         8.32050294 0.        ]

另一个例子#

让我们看一个稍微复杂的模型。在这个例子中,您将观察一根带有两根电缆和一个作用力的横梁。这次您需要找出电缆中的张力以及横梁的反作用力。(来源:工程师矢量力学:静力学与动力学,问题4.106)

image.png

定义距离a为3米

和以前一样,首先将每个相关点的位置定义为一个数组。

A = np.array([0, 0, 0])
B = np.array([0, 3, 0])
C = np.array([0, 6, 0])
D = np.array([1.5, 0, -3])
E = np.array([1.5, 0, 3])
F = np.array([-3, 0, 2])

从这些方程中,您首先使用单位向量确定向量方向。

AB = B - C
AC = C - A
BD = D - B
BE = E - B
CF = F - C

UnitBD = BD / np.linalg.norm(BD)
UnitBE = BE / np.linalg.norm(BE)
UnitCF = CF / np.linalg.norm(CF)

RadBD = np.cross(AB, UnitBD)
RadBE = np.cross(AB, UnitBE)
RadCF = np.cross(AC, UnitCF)

这让您将作用在系统上的张力(T)和反作用力(R)表示为

\[\begin{split}\left[ \begin{array} ~1/3 & 1/3 & 1 & 0 & 0\\ -2/3 & -2/3 & 0 & 1 & 0\\ -2/3 & 2/3 & 0 & 0 & 1\\ \end{array} \right] \left[ \begin{array} ~T_{BD}\\ T_{BE}\\ R_{x}\\ R_{y}\\ R_{z}\\ \end{array} \right] = \left[ \begin{array} ~195\\ 390\\ -130\\ \end{array} \right]\end{split}\]

力矩表示为

\[\begin{split}\left[ \begin{array} ~2 & -2\\ 1 & 1\\ \end{array} \right] \left[ \begin{array} ~T_{BD}\\ T_{BE}\\ \end{array} \right] = \left[ \begin{array} ~780\\ 1170\\ \end{array} \right]\end{split}\]

其中\(T\)是相应电缆中的张力,\(R\)是相应方向上的反作用力。那么您就有六个方程

\(\sum F_{x} = 0 = T_{BE}/3+T_{BD}/3-195+R_{x}\)

\(\sum F_{y} = 0 = (-\frac{2}{3})T_{BE}-\frac{2}{3}T_{BD}-390+R_{y}\)

\(\sum F_{z} = 0 = (-\frac{2}{3})T_{BE}+\frac{2}{3}T_{BD}+130+R_{z}\)

\(\sum M_{x} = 0 = 780+2T_{BE}-2T_{BD}\)

\(\sum M_{z} = 0 = 1170-T_{BE}-T_{BD}\)

您现在有五个未知数和五个方程,可以求解

\(\ T_{BD} = 780N\)

\(\ T_{BE} = 390N\)

\(\ R_{x} = -195N\)

\(\ R_{y} = 1170N\)

\(\ R_{z} = 130N\)

总结#

您已经学会了如何使用数组来表示三维空间中的点、力和力矩。数组中的每个条目都可以用来表示分解为方向分量的物理属性。然后,这些属性可以轻松地使用NumPy函数进行操作。

其他应用#

同样的过程可以应用于动力学问题或任何维数的问题。本教程中的例子假设的是静力平衡中的三维问题。这些方法可以很容易地用于更多样化的问题。更多或更少的维度需要更大或更小的数组来表示。在经历加速度的系统中,速度和加速度也可以类似地表示为向量。

参考文献#

  1. 工程师矢量力学:静力学与动力学 (Beer & Johnston & Mazurek & et al.)

  2. NumPy参考