概括
在本系列讨论的第一部分中，我将探索并验证手工计算与物体掉落到刚性表面上的时间相关响应之间的关系，以及在撞击期间发生的偏转和应力。下面讨论的概念对于理解如何使用 ANSYS Mechanical 设置和解决瞬态结构分析至关重要。
以下是本次讨论中涉及的主题列表（按出现顺序）：
- 势能
- 弹性能量
- 方向刚度
- 静态结构分析
- 平均压缩力
- 动能
- 冲击速度
- 影响期
- 固有频率
- 模态有效质量
- 瞬态结构分析
- 分析持续时间
- 时间步长频率
- 复杂瞬态位移和应力结果
- 平均压缩应力
细节
我开始探索手工计算与动态事件（例如一个物体撞击另一个物体）产生的最大挠度和应力结果之间的关系。我使用 ANSYS Mechanical 来验证假设并提供这些动态事件的更多细节。我对有限元方法充满信心，并期待与您分享我的发现。
我的探索始于一个简单的例子：一根圆柱形杆掉落在刚性表面上。在这个场景中，以及与此讨论相关的所有场景中，我假设材料具有弹性行为，并且任何载荷都不会导致塑性或损坏。这意味着几何结构中存在能量守恒。
在此过程中，我做出了一系列假设，并进行了有限元分析以验证或反驳这些假设。这个过程让我提出了更多问题，进行了更多分析，并得到了更多答案。最后，我能够自信地描述这一动态事件的许多方面，我将尽可能简洁地介绍所有这些方面。
让我们从最简单的想法开始。
第一个例子探讨了圆柱形棒材跌落到坚硬表面上的情况。该棒材的直径为 25.4 毫米，长度为 254 毫米，质量密度为 7.85e-06 [kg/mm³]。该棒材将从距底面 1 米处跌落。
图 2：圆柱形杆
在这种情况下，势能必定等于弹性能。
势能：
弹性能量：
让我们将下落圆柱体的势能与其弹性能相等，以了解我们是否能够准确预测挠度。
我们假设杆将以负 Y 方向移动的方式下落。因此，第一步是确定几何体在 Y 方向上的刚度。
对于我们的模型，该刚度可以从理论上推导，也可以通过使用有限元方法（施加已知载荷并将该载荷除以计算出的 Y 方向挠度）来推导。这样我们就可以将有限元估算的结果与理论推导值进行比较，并对我们的方法建立信心。
让我们推导出理论刚度。对于 Y 方向，将应用以下公式：
现在，让我们通过两种不同的方法通过有限元分析来计算几何刚度，以便我们了解不同类型的负载所涉及的不同考虑因素。
这里我们考虑一个固定底座，同时对圆柱体的远端施加一个力，使得底座处的反作用力与施加在圆柱体顶部的力大小相等，符号相反。只要我们的材料属性是线性的，并且我们不考虑几何非线性（大变形），这种压缩载荷的大小并不重要。
图 3：从顶部加载并在底部支撑的杆的静态结构分析。
在这种情况下，顶部表面施加了 9.9079 N 的载荷，而底部表面因 Y 位移而固定。此载荷的大小稍后将变得更加清晰。
图 4：静态结构挠度结果
解决静态结构分析显示，我们的压缩挠度最大值为 2.4833e-5 毫米。基于 9.9079 N 负载，我们可以计算出以下平均刚度。
可以看出，有限元方法与理论方法的差异仅为 0.00125%。这让我们感觉很好😊
接下来，让我们考虑不同类型的载荷，例如重力载荷。在这里，力从顶面移除，但底面仍然受到约束，重力以 9806.6 mm/sec^2 施加。底部的反作用力为 9.9079 N，但作用在顶面上的力肯定为零。挠度结果显示不同的挠度大小为 1.2421e-5 mm。
图 5：静态结构分析挠度结果
如果我们尝试使用与上述公式相同的形式来计算刚度，则会得出错误的值，因为我们不会考虑作用于几何体的平均压缩力。在这种情况下，我们可以轻松估算平均压缩力，因为我们的几何体具有恒定的横截面，并且我们估计的平均刚度可以按如下方式计算。
该估算值与我们理论计算值相差0.0348%。因此，通过两种不同的加载方法，我们可以自信地估算几何方向刚度。
利用我们估计的方向刚度以及模型和材料属性，我们可以估算出下面的最大挠度。
约为 ~0.223 毫米。现在我们需要一个瞬态结构模型来模拟圆柱杆与地面碰撞时及之后的挠度。
图 6：瞬态结构分析模型和约束
我们的模型假设杆以 4428 毫米/秒的速度沿负 Y 方向移动。这是杆从 1 米高处掉落时的速度。我们可以通过将杆的势能与撞击瞬间的动能相等来计算。
动能：
所以：
添加约束以保持几何体的垂直高度，同时在底座上添加仅受压的支撑。定义分析持续时间（稍后会详细介绍），以及收集结果的频繁捕获率。
分析结果很复杂，因为当杆的顶部在撞击后向下偏转时，杆的底部可能会因与仅压缩支撑相关的刚度而略微压缩。这两个表面之间的偏转差异描述了杆的整体压缩。
图 7：瞬态结构挠度结果
蓝线的平坦时期“撞击区杆变形”表示底部发生接触的时期，或与撞击相关的停留时间。如果我们重新绘制数据以关注这一时期，我们可以更轻松地进行多项观察。
图 8：瞬态结构压缩挠度结果
我添加了一条绿线，表示钢筋顶部和底部表面挠度之间的差异，从而得出钢筋的总压缩量随时间的变化。两条垂直虚线表示与橙色和绿色线的最小结果相关的时间间隔。可以看出，钢筋的最小压缩量与钢筋顶部的最小压缩量同时发生。由此，我们可以估计钢筋的最大压缩量为 0.216 毫米，比我们估计的最大挠度 0.223 毫米小不到 3%。
然而，从这个例子及其结果中，还有更多值得探索和学习的东西。
现在让我们看一下设置瞬态结构分析的方面，并专门讨论与我们的事件相关的时间段，并确定如何为任何几何图形估计这个时间段。
如果我们观察蓝线从负位移变为正位移的时间（0.00010026 秒），这等于与轴压缩相关的周期。
我们的轴可以看作是一根弹跳的弹簧。冲击持续时间反映了其在压缩方向上的固有频率周期的 ½。因此，以下应该是正确的。
现在让我们探讨一下是否可以支持这一理论，以及如何估计这一理论的固有频率，以及我们希望在未来分析的任何几何形状。
我们可以改变应用于上述测试模型的“仅压缩”支撑，并在相同表面上施加 Y 方向约束，并直接在 ANSYS 中求解纵向振动频率。
图 9：模态分析基频
这里我们看到纵向振动频率为 4967.4 Hz。这表明我们根据瞬态结构位移结果估算的固有频率准确度在 0.4% 以内。这一发现证实了我们的理论，即对于两个撞击物体中较软的物体，撞击持续时间等于基频周期的 ½。
图 10：影响期
在我们的示例分析中，我们只有一个部件，我们假设仅受压的垂直支撑代表另一个非常坚硬的物体。因此，我们的轴是这两个撞击伙伴中更灵活的一个。我们的瞬态结构分析被定义为分析持续时间等于该固有频率的周期，我们预计撞击将发生在该周期的前 1/2 秒内。
但是如果我们要建立瞬态结构分析，如何计算这个固有频率呢？
我们可以像我刚才演示的那样，对几何结构进行模态分析，同时使用约束来保持冲击姿势并直接求解冲击压缩的振型。但是有没有一种分析方法可以完成同样的事情呢？
为了探索这一点，我们需要计算几何体在冲击方向上的固有频率，计算方法如下。
我们已经探讨了如何估计方向刚度，但现在让我们考虑几个例子来更好地理解模态有效质量的考虑。
我们的系统的总质量是：
由于几何形状的一端受到约束，因此系统的有效质量将小于总质量。根据几何形状的复杂性和用于固定几何形状的约束，可能很难估计参与的质量量。在我们的例子中，几何形状很简单，对于固定在一端的挤压横截面，有效质量通常等于系统总质量的 1/3。如果是这种情况，那么我们的基频如下：
更明确地可以表述为：
这个值（5478 Hz）比我们用有限元法计算的固有频率（4967.4 Hz）大10％以上……为什么？
答案与模态有效质量的估计有关。如果我们对纵向刚度的计算有信心（……我们确实有信心……），那么我们可以更改模型设置，以便更轻松地验证基频的计算。
我们将模型材料的质量密度设置为零，然后在几何体的远端添加 1.01031 千克的质量，这将产生以下固有频率。
运行有限元分析后，我们发现基本纵向频率为 3154.6 Hz。
图 11：模态分析基频
这与我们手工计算的结果有 0.27% 的差异。因此，我们对有限元方法更加有信心，并确定了我们在利用整个几何结构中的分布质量进行手工计算时使用的有效质量估计存在缺陷。
通过重新排列我们的方程，我们可以通过考虑 ANSYS 计算的固有频率并考虑分布质量来求解该有效质量，如下所示：
我们发现该有效质量为系统总质量的 0.405。
因此，对于其他更复杂的几何形状，我们可以预期这个质量分数是唯一的，并且很难从理论上推导，并且我们可以使用有限元方法可靠地求解固有频率。
现在我们了解了如何确定给定几何形状的固有频率，我们需要确定在瞬态分析期间应以什么频率收集分析结果。
为此，我们应该从两个不同的角度考虑分析步长如何影响分析。一个角度是使用非线性响应（例如正弦波）的分段线性表征可以预期什么。另一个角度涉及实际分析结果如何随捕获率频率的变化而变化。因为我们正在考虑冲击问题，并且我们假设冲击周期与柔性部件的基本振动有关，所以探索步长与它如何影响表示正弦曲线的精度之间的关系是有意义的。
为了评估这种关系，我们将正弦曲线前半部分的幅度加权平均值与同一跨度的理论平均值进行比较，其等于：
现在我们来探索一下如何使用不同的步数来准确地表示这条正弦曲线，然后使用梯形法则计算前半部分曲线下的面积，然后将该面积除以该跨度的周期。
图 12：4 步；估计加权平均值 = 0.5；差异为 21.46%
图 13：8 步；估计加权平均值 = 0.60536；差异 5.19%
图 14：16 步；估计加权平均值 = 0.6284；差异 1.29%
图 15：32 步；估计加权平均值 = 0.6346；差异 0.32%
图 16：512 步；估计加权平均值 = 0.6366；差异 0.001%
事实上，我们从未实现与理论值的完全一致，我们可以绘制出随着时间步长的增加，这种准确度如何变化，如下所示。
图 17：估计平均震级收敛图
从这个检查中，我们可以看到，随着近似步骤数的增加，准确度会迅速大幅提高。但这并没有考虑到随着时间步骤数的增加，结果会如何发展。
为了探索对分析结果的影响，我们将绘制最大杆压缩结果随时间的变化图，同时考虑解决方案中定义的不同时间步长数。
图 18：考虑 4 个时间步长的杆位移。
图 19：考虑 8 个时间步长的杆位移。
图 20：通过 16 个时间步长，我们可以看到更大的杆压缩位移。
图 21：采用 32 个时间步长的瞬态最大杆偏转。
根据我们之前对计算平均正弦幅度的分析，我们发现，通过考虑 32 个步骤，我们应该预计与理论结果的差异高达 0.32%。但是，我们可以看到，与之前的时间步长探索相比，我们的杆压缩位移响应仍在形成中。如果我们考虑额外的时间步长会发生什么？
图 22：采用 64 个时间步长的瞬态最大杆偏转。
64 个时间步长产生的最大压缩杆偏转稍大一些，但大致相似，但是，有迹象表明在杆回弹过程中会产生更高频率的振荡。
图 23：采用 128 个时间步长的瞬态最大杆偏转。
考虑 128 个时间步长时，轴的峰值压缩比时间步长较少时大。最大压缩挠度比考虑 32 个时间步长时大 6%。
图 24：使用最多 512 个时间步长的瞬态最大杆偏转。
将 256 个时间步长的分析结果与考虑 512 个时间步长的结果进行比较，可以发现通过考虑更多时间步长，可以计算出动态响应的进一步分辨率。下图说明了在瞬态结构分析中，我们的最大压缩杆挠度将如何随时间步长的变化而变化。
图 25：最大杆压缩位移与分析时间步数
最终，对于给定的网格尺寸，我们达到收益递减点......其中延长分析时间和增大结果文件所涉及的努力无法通过结果幅度的显著变化来证明。
基于这两个考虑，我们可以看到，通过使用 32 个时间步长，我们可以期望表征出差异高达 0.32% 的正弦曲线，但与使用 512 个时间步长的相同分析相比，可能只能产生准确度为 94% 的位移结果，与我们的理论最大挠度相比，准确度仅在 91% 以内。我们在下图中看到我们的最大杆压缩结果与我们理论估计的最大值的比较情况。
| 子步骤数 | 最大压缩 | 理论最大值 | ％ 不同之处 |
| 4 | 0.152139 | 0.222859 | 31.73％ |
| 8 | 0.174023 | 0.222859 | 21.91％ |
| 16 | 0.191486 | 0.222859 | 14.08％ |
| 三十二 | 0.202746 | 0.222859 | 9.03% |
| 64 | 0.210712 | 0.222859 | 5.45％ |
| 128 | 0.214455 | 0.222859 | 3.77％ |
| 256 | 0.215752 | 0.222859 | 3.19％ |
| 512 | 0.216232 | 0.222859 | 2.97％ |
图 26：最大杆压缩位移与分析时间步骤
随着时间步长的增加，我们的有限元分析结果与理论推导值之间的一致性增加。
图 27：与理论最大挠度的百分比差异
我们从该图中看到，当时间步数在 64 到 128 之间时，收敛速度有显著的改善，但是随着时间步数的增加，变化很小。
但对于压力等其他结果而言，这一切意味着什么？
下面，我将绘制几种分析场景下每个时间的最大冯·米塞斯应力。
总结一下，这是我们在观察压力时所看到的：
| 子步骤数 | 最大平均应力 | 最大峰值应力 | 平均应力差异百分比 |
| 4 | 120.04 | 301.98 | 29.62% |
| 8 | 136.4 | 316.61 | 20.03% |
| 16 | 150.57 | 312.42 | 11.72％ |
| 三十二 | 159.15 | 301.68 | 6.69％ |
| 64 | 165.4 | 283.51 | 3.03% |
| 128 | 168.84 | 303.99 | 1.01％ |
| 256 | 170.05 | 354.35 | 0.30% |
| 512 | 170.56 | 442.08 | 0.00% |
图 33：最大杆压缩位移与分析时间步骤
我们可以看到，随着时间步长的增加，峰值冲击应力会变得更大，在冲击开始时显示最大值，并在整个事件过程中变小。杆中的平均应力确实从零上升到最大压缩时的峰值，然后在冲击结束后回落到零。我们可以将这些最大平均应力绘制在下面。
图 34：最大杆压缩位移与分析时间步骤
为了正确理解这些最大平均应力，我们可以按如下方式计算杆中的理论最大平均应力：
再次绘制我们的最大平均应力与理论计算的最大平均值，我们可以总结出百分比差异如下：
| 子步骤数 | 最大平均应力 | 理论最大平均应力 | ％ 不同之处 |
| 4 | 120.04 | 175.5 | 31.60％ |
| 8 | 136.4 | 175.5 | 22.28％ |
| 16 | 150.57 | 175.5 | 14.21％ |
| 三十二 | 159.15 | 175.5 | 9.32% |
| 64 | 165.4 | 175.5 | 5.75％ |
| 128 | 168.84 | 175.5 | 3.79％ |
| 256 | 170.05 | 175.5 | 3.11％ |
| 512 | 170.56 | 175.5 | 2.81％ |
图 35：最大杆压缩位移与分析时间步骤
再次，就像偏转一样，我们承认，在考虑事件的 64 到 128 个时间步长时，我们的最大平均应力与理论计算的最大平均应力相比的百分比差异有显著改善。然而，使用更多的时间步长并不能提高我们结果的准确性。
结论
我们验证了用于估计动态事件（例如弹性物体掉落到坚硬表面上）的最大挠度和平均应力的能量方法可以产生合理的结果，并且当考虑 64 到 128 个时间步长时，精度会得到最显著的提高，该时间步长的周期等于我们撞击物体的固有频率（针对撞击方向计算）。我们还发现，近似该固有频率的最佳方法不是手工计算，而是使用有限元方法执行的模态分析，特别是因为参与模态有效质量的数量可能难以通过任何其他方法估算。
