在 LS-DYNA 中将应力转换为用户定义的坐标系

介绍

通常，使用 LS-DYNA 或 Ansys Mechanical 等仿真工具解决工程问题需要将张量结果与解析解进行比较。一个这样的例子是加压圆柱体，其中圆周应力或环状应力是感兴趣的主要应力度量。例如，如果对具有复杂端部处理的几何结构进行此类仿真，则希望使用环向应力（例如在容器中心）来验证仿真是否在不存在解析结果的更复杂几何结构区域中提供准确的结果。

在大多数现代后处理器中，应力分量可能在笛卡尔坐标、主方向，甚至圆柱或球坐标中绘制条纹。但是，可能需要将结果将张量结果旋转到所述后处理器中“开箱即用”的方向。在本文中，我们将提供一个可验证的张量转换示例，该示例使用 Python 3.11 脚本从 LS-PrePost 中以笛卡尔坐标读取加压圆柱体的应力条纹图结果，然后将应力张量旋转为圆柱坐标以计算环向应力。

个案研究

我们将用来说明这个概念的模型是一个 d = 100 mm 平均直径 x l = 150 mm 长 x t = 2mm 厚的圆柱体，内部压力为 p = 4 MPa，使用 LS-DYNA 求解。教科书般的图像1如下所示：

 

 

有限元模型

圆柱体的材料是 Workbench Engineering Data 中的标准线性各向同性弹性结构钢，网格由 2.5 毫米大小的壳单元组成，使用默认截面属性。网格尺寸提供均匀的网格。

 

 

载荷和边界条件

根据众所周知的公式，内部压力p选定为4 MPa，这样环向应力将等于100 MPa

 

其中r = 平均圆柱半径，t = 壁厚。此外，均衡约束（又称 3-2-1 约束）应用于黄色显示的顶点，以提供无应力约束。

 

分析设置

分析类型设置为隐式，伪结束时间为 0.01 秒。选择隐式分析是为了加快分析运行时间。

模拟结果

利用 Ansys Mechanical 中的圆柱坐标系，我们可以轻松地对（未平均的）环向应力进行后处理，得到 100 MPa 的结果，从而与理论相匹配。

 

这个结果表明我们的模拟与理论相符。

使用 Python 重建环向应力

现在已经使用 Ansys Workbench LS-DYNA 找到了环向应力结果，我们将通过处理从 LS-PrePost 输出的笛卡尔应力分量来重建相同的环向应力结果。

导出笛卡尔应力分量

为了获得旋转应力张量所需的数据，必须从 LS-PrePost 导出所有六个分量。导出每个应力分量的过程将针对 X 分量进行说明，并且必须针对 Y、Z、XY、YZ 和 ZX 分量重复此过程：

1. 像往常一样加载二进制 D3Plot 文件并将结果时间设置为结束状态。
2. 在边缘图中绘制 X 应力。

   

3. 使用下面显示的设置输出结果，例如，将文件命名为“xstress”。

   

4. 点击“完成”并对其他 5 个压力分量重复此操作。
5. 确保所有六个文件都在同一个目录和python脚本中。

输出文件格式

输出文件的格式为 LS-DYNA 关键字格式，这是一种固定宽度的文件格式。Python 脚本利用固定宽度的结构来解析文件的每个部分。

*ELEMENT_SHELL 部分

元素连接表位于文件的 *ELEMENT_SHELL 部分中，其中的一小部分如下所示：

列定义和宽度如下：

柱子	数据	宽度
1	元素 ID	8
2	零件编号	8
3	节点 1	8
4	节点 2	8
5	节点 3	8
6	节点 4	8

*NODE 部分

类似地，节点坐标在 *NODE 部分中给出，其中的摘录如下所示：

列定义和宽度如下：

柱子	数据	宽度
1	节点 ID	8
2	X 坐标	16
3	Y 坐标	16
4	Z 坐标	16

 

$SHELL_ELEMENT_RESULTS 部分

$SHELL_ELEMENT_RESULTS 部分提供元素每个积分点的输出结果。在本例中，我们使用具有单个积分点的壳元素。本节的摘录如下所示：

 

列定义和宽度如下： 

柱子	数据	宽度
1	元素 ID	8
2	元素结果	16

 

$SHELL_ELEMENT_CENTROID 部分

$SHELL_ELEMENT_CENTROID 部分提供变形状态下每个元素的质心和体积。本节摘录如下：

列定义和宽度如下：

柱子	数据	宽度
1	元素 ID	10
2	X 坐标	16
3	Y 坐标	16
4	Z 坐标	16
5	体积	16

 

Python 算法

下面的伪代码解释了用于转换压力的算法：

1. 读取六个应力文件中的每一个并解析各个部分以获得元素表、节点坐标、元素结果（非平均应力）和元素质心坐标（在变形配置中）。
2. 将元素表与节点坐标和元素结果合并。
3. 形成元素结果（应力）张量S。
4. 根据右手定则，使用两个元素边缘的叉积计算每个元素质心处的法向单位向量。我们将其称为元素Z方向。 
   1. 对于本例中使用的特定壳单元，两条边从节点 1 到 4 ( 1->4 ) 形成，从节点 1 到 2 ( 1->2 ) 形成。因此， 1->4 x 1->2提供向外的单元法线，然后进行归一化以使其成为单位向量。
5. 由于我们使用的是柱坐标系，因此我们通过将全局 Z 轴与元素Z方向向量相交来形成元素X方向单位向量。
   1. 该元素X方向向量与元素法向量不正交，因此我们利用 Gram-Schmidt 正交化对其进行正交归一化。
6. 有了元素X和Z单位向量，我们取Z x X来获得元素Y方向单位向量，从而完成我们的正交集。
7. 使用X - 、Y - 和Z - 元素单位向量作为列，构成变换矩阵Q。
8. 使用张量变换规则S loc  = Q ' S Q旋转应力张量，其中 ()' 表示转置。
9. 将结果写入输出 .csv 文件。

Python 计算环向应力

此处显示了 Python 脚本结果的摘录，其中左图显示笛卡尔坐标中的元素应力；右图显示了相应的环向应力S11_loc ，它是由变换和变换矩阵Q的分量引起的：

 

   

 

从右图中我们可以看到，计算出的环向应力S11_loc与 Ansys Mechanical 的输出非常匹配。

结论

虽然可以轻松从 Workbench Ansys LS-DYNA 或 LS-PrePost 中获得环向应力（作为第一主应力），但本文中的步骤展示了如何通过将结果转换为笛卡尔坐标来实现张量结果的旋转以获得更一般的方向。上述过程也可以使用 PyAnsys 的工具在更大的自动化框架中使用。

进一步

1. 修改 Python 脚本来查看其他方向的应力（或应变）或应用于您自己的模型来获得感兴趣的应力（或应变）结果。
2. 修改代码以用于具有多个积分点的实体元素或壳元素。
3. 使用此脚本和 PyAnsys 来自动化该过程。

可下载项目

Ansys 2024 R2 Workbench 存档，python 代码和结果

参考

1. 压力容器图像取自Stress and Strain - Cylindrical bar of cross-sectional area - Engineersfield