ANSYS计算节点内力的宏

今天根据ansys.net上的一个宏稍微修改了一下，做了一个计算节点内力的宏，贴在下面。这个宏的原始作者是Shen-Yeh Chen。这里是原始描述和原始宏。

这个宏可以计算指定节点上的内力。通常情况下，ANSYS只计算被限定自由度的节点的反应力。如果模型里有两个材料，并且这两个材料用连续网格模拟(而不是接触)，有些时候我们想知道这两个材料间的相互作用力，而在ANSYS里面并没有直接的命令。使用这个宏可以做到这一点。

这个宏的使用法方法是：0.当然是把这个宏命名并放在相应的目录下；1.选择要计算内力的节点；2.选择这些节点一边的单元，因为这些节点代表一个截面，而这个截面会有两边有单元。选择不同的边，内力的方向是不同，如果两边都选的话，内力就会平衡为零了。因此必须要先选择一边，但是不必要只选和节点相邻的一层，因为宏里面已经有命令去选择和节点相连的单元；3.运行宏；4.在工作目录下，会有一个x_nfor2f.ans文件，这个文件的最后一行，就是在x,y,z三个方向的内力分量了。

这个宏必须要在通用后处理器中运行。

! ========================================

! Shen-Yeh Chen

! 05/01/1998

! 12/10/2001

! =========================================

!

! This macro will transfer the element force on the selected

! node into cross-section force on each node, and write

! all the force on each node into the "F" command. Note that

! the macro go through each SELECTED node, sum up all SELECTED

! elements that are attached to the SELECTED node. Then the

! force on each nodes is added a negative sign and printed

! into the 'F' command.

! ---------------------------------------------------------

! Syntax :

! x_nfor2f

!

! Input :

! None

!

! Output :

! file : x_nfor2f.ans : nodal force for the section load

! x_nfor2f.cdb : nodal force & nodal coordinates/ID in

! CDB format

! nodal components : x_nfor2f : defined in x_nfor2f.cdb

! X_NFOR2F_FX : summation of FX

! X_NFOR2F_FY : summation of FY

! X_NFOR2F_FZ : summation of FZ

! Temporary components

! X_NFOR2F

! X_NFOR2F_E

!

! Temporary components

! x_nfor2f.db

!

! Temporary variables

! AR30,AR31,AR32,AR33

!

! External Programs called

! None

!

!

! ========================================

/POST1

*GET,ZNNODE,NODE,,NUM,MAX

CM,X_NFOR2F,NODE

CM,X_NFOR2F_E,ELEM

ESLN,R

!

AR30=0

*CFOPEN,x_nfor2f,ans

*VWRITE

('/PREP7')

X_FOR2F_FX=0.0

X_FOR2F_FY=0.0

X_FOR2F_FZ=0.0

!

! LOOP THROUGH ALL SELECTED NODES

:LINE01

AR30=NDNEXT(AR30)

*IF,AR30,EQ,0,THEN

*GO,:LINE02

*ENDIF

NSEL,R,NODE,,AR30

FSUM

!

!

! GET THE FORCE

*GET,AR31,FSUM,0,ITEM,FX

*GET,AR32,FSUM,0,ITEM,FY

*GET,AR33,FSUM,0,ITEM,FZ

!

AR31=-AR31

AR32=-AR32

AR33=-AR33

X_FOR2F_FX=X_FOR2F_FX+AR31

X_FOR2F_FY=X_FOR2F_FY+AR32

X_FOR2F_FZ=X_FOR2F_FZ+AR33

!

! WRITE TO FILE

*VWRITE,AR30,AR31

('F,',F6.0,',FX,',E12.5)

*VWRITE,AR30,AR32

('F,',F6.0,',FY,',E12.5)

*VWRITE,AR30,AR33

('F,',F6.0,',FZ,',E12.5)

!

! SELECT THE NODES

CMSEL,,X_NFOR2F

*GO,:LINE01

!

!

! SAVE THE STATE TO x_nfor2f.db

:LINE02

*VWRITE,X_FOR2F_FX,X_FOR2F_FY,X_FOR2F_FZ

('FX = ',E12.5,', FY = ',E12.5,', FZ = ',E12.5)

*CFCLOS

STAR-CD和STAR-CCM+ ：冠军的选择

最近CD-Adapco网站上登出了一则新闻。新闻的题目是：STAR-CD and STAR-CCM+: Still the choice of champions (STAR-CD和STAR-CCM+：冠军的选择)。这则新闻被很多网站转载。新闻说的是2006赛季F1冠军车队雷诺车队使用的是CD-Adapco提供的CFD软件。我每天上班的路上看的一份小报几乎每天都要介绍F1大赛，我对当中一些人得的名字还是比较熟悉的。

新闻稿中说，自2001年雷诺车队成立以来，CD-Adapco就是雷诺车队的CFD软件提供商。RS26引擎 和R26底盘 的设计都和CD-Adapco的STAR-CD或STARE-CCM+软件有关。CD-Adapco是F1大赛中最成功的CFD软件提供商，因为，在2006赛季有得分的车队中，有70%的车队使用STAR-CD或者STAR-CCM+作气动或者引擎流动模拟。看来STAR-CD在汽车行业的优势是当之无愧的。关于F1大赛的更多信息可以参见搜狐体育的F1频道 。

多面体网格的优点

多面体网格(polyhedral mesh)是最近CFD求解器开发的热点，它已经被使用在CD-Adapco的STAR-CCM+中，STAR-CD V4也支持多面体网格，Fluent V6.3也将支持多面体网格。其实在CFX中，由于采用基于节点的控制体积(Vertex-based control volume)，用户提供给求解器的网站是正常的四面体或六面体网格，而在求解器内部会生成多面体网格，并且声称已经使用这个技术超过15年。下面的文章翻译自CD-Adapco上面的一个新闻发布(链接)，作者是Milovan Peric和Stephen Ferguson。

自CD-Adapco推出STAR-CCM+求解器和自动多面体网格生成器以来， 我们经常被问到一个问题：为什么我们会采用多面体网格，而不是四面体网格？ 在这篇文章里，我们会就这个问题提供答案和实例。

四面体单元是最简单的体积单元；他们的表面是平面碎片 (plane segment)， 因此面和体的中心点很容易定义。四面体网格也是相对容易自动生成的。它们几乎是自动网格生成技术中的标准应用，主流CFD求解器都支持四面体网格。但是，在另一方面，四面体网格不能被拉伸太多，因此为了保证在边界层、长的管道、或小的缝隙的精确度， 需要使用比结构化网格(六面体)多的多的四面体网格。棱形网格层可以部分解决这个问题。

四面体控制体只有四个相邻网格，在使用标准近似(线性形函数)来计算单元中心的梯度时，会比较困难。一个困难是相邻网格节点的空间位置，因为它们可能会几乎处在一个平面上，使得几乎不可能计算垂直于该平面的梯度。当单元邻近边界时会有另外一个可能碰到的问题：即使只有一个面是边界面，其它三个面的分布可能会很不好；尤其是在邻近线或者角落时，一个单元可能只有两个或者一个相邻单元，这不光会降低精确度，而且会引起严重的数值问题。

为了在四面体网格上实现精确求解和好的收敛性，需要特别的离散技术和大量的单元。这些解决方案并不是最优化的；第一种方法使得程序更加复杂、更难扩展和维护，而第二个方法增加了内存和计算时间需求。多面体网格提供了和四面体网格一样的自动网格生成功能，并且能解决这些缺点。

多面体网格的主要优点是它有多个相邻单元(通常10数量级) ，和四面体网格相比，使用多面体网格能更好的近似梯度(使用线性形函数以及最近单元的信息)。即使靠近线和角，一个多面体单元也可能有几个相邻单元，从而可以合理的预测梯度和当地流动分布。更多的相邻单元意味着更多的内存使用量和更多的计算操作，这可以从更高的精度中得到补偿。下面将要讨论这个问题。

多面体网格也不如四面体网格多拉伸敏感。智能化的网格生成和优化技术也提供了无限的可能性：通过加入新的节点、单元或面，单元可以自动结合、分裂、或修改。实际上，网格质量在将来会有显著提高，并提高求解器效率和精确度。而且，求解器中在使用四面体网格时需要的特别处理将不再必要。在使用多面体网格时，要模拟当地网格细化、滑移网格界面、周期性边界条件等，只需要生成特别的多面体网格，对求解器来说它们都是一样的。

多面体网格在处理回流时更有优势。比如测试显示在模拟正方体空腔驱动流时，要得到指定的精度所需的多面体单元数量要比笛卡尔六面体(通常认为对长方体计算域，六面体是最理想的)少得多。这可以解释为：对六面体单元，有三个最佳流动方向可以达到最大精度(和三对平行面垂直的方向) ；对有十二个面的多面体单元，有六个最佳流动方向，并且有更多的相邻单元，因此可以用较少的单元数得到较高的精确度。对于不同网格类型特点和一个测试结果可以参考Peric (2004) 发表的文章。

很多实际应用的比较分析显示，要达到和四面体网格同样的精确度，使用多面体网格只需要四分之一的单元数、一半的内存、十分之一到五分之一的计算时间。另外，多面体网格有更好的收敛性，通常不需要调整默认的求解器参数，如图1至3所示。

图1. 模拟发动机水套内流动所用的多面体网格

图2. 最细的四面体网格(左图)和最细的多面体网格(右图)所预测的水套壁面压力分布

图3. 不同四面体网格和多面体网格所预测的水套入口出口压力差

这个测试模拟的是发动机水套；监测的是进口和出口的压力差。计算分别在6个多面体网格(单元数21872至593888)和6个四面体网格(单元数39587至2322106)上进行。在所有情形中近壁面处均使用了棱形网格层。图1是一个多面体网格。图2是在最细的多面体网格和最细的四面体网格上计算所得的压力分布。因为网格很小，这两个计算结果非常相似。为了分析结果对网格的依赖性，图3显示了在不同网格上计算得到的压差。该图显示这两种网格类型都可以收敛得到网格无关解。所有的情形使用了同样的离散和求解方法(对流项采用二阶迎风差分格式)。

重要的是，在多面体网格上得到的解总是比有相近网格数量的四面体网格要精确。实际上，从有65513单元的多面体网格上得到的解比从有393273个单元(大约6倍多)的四面体网格上得到的解还要稍微精确。多面体网格所需的计算时间只是得到同等精度的四面体网格的十分之一。到目前为止所进行的比较得到了相似的结果，这验证了CD-Adapco发展这一技术的决定是正确的。随着CD-Adapco的全自动多面体网格生成器的面世，多面体网格将会成为业界标准。

参考文献：

1 M. Peric: Flow simulation using control volumes of arbitrary polyhedral shape, ERCOFTAC Bulletin, No. 62, September 2004.

一周没有发言了

原因是参加了一个关于CAE的conference。

Abaqus 6.6学生版发布

达索系统(Dassault systems，全球网站，中国网站)在本月10号发布了Abaqus
6.6的学生版(Abaqus全球网站，中文网站，新闻发布)。Abaqus学生版本是研究版本的教学版本的补充，对学习Abaqus进行有限元分析的学生和专业人士有很大的帮助。



在Abaqus网站上对产品的特色有详细介绍。主要的限制是用户不能求解和后处理超过1000个节点的模型。另外，用户不能编程扩展。Abaqus也原则上不提供技术支持。对软件的功能和模型没有特别的限制，也就是说用户可以接触练习所有Abaqus模型。价格是99美元一套。这样的价格应该是非常合理的。不知道ANSYS有没有提供这样的学生版本，如果有的话，我一定会买一套，但是在ANSYS的官方网站上并没有相关信息。看来ANSYS对教育市场还是比较保守的。

压力和应力的区别

今天有个朋友问我压力和应力的区别。压力(pressure)和应力(stress)有相同的单位帕斯卡(Pa, N/m^2)，但是它们的区别是什么呢？

通常来说，固体内部只有应力，只有流体内部才既有应力又有压力。应力是物体对外力的方应，而压力和温度、压力一样是流体的一种状态。应力是矢量，有方向，而压力是标量没有方向。应力是由于材料的弹性或粘性而产生的，压力是由于流体内部分子的无规则热运动而产生的。

虽然在固体力学中，一些外力可以定义成压力。比如，当我们用手指摁桌面时，会给桌面施加一个压力。这里的压力是指和应力相对应物体内部的压力。

综合来说，应力是在外力作用下固体或流体由于材料的弹性或粘性对外力的反应。压力是由于流体内部的无规则分子热运动而产生的相互作用力，和材料的粘性无关。

压力的定义是相对简单的(wikipedia ，NASA Glenn Research Center )。因为压力是由于无规则分子热运动产生的，因此各向同性，是个标量。而应力的定义就要复杂得多。

1. 应力是表面力

当一个物体受到外力作用时，物体内部会相应产生内力来平衡外力。而内力的分布和该物体的强度以及断裂等直接相关。内应力(stress)就是是描述物体内部力的分布的矢量。值得注意的是，应力是表面力(surface force，单位面积上的力)，而不是体积力(volume force，单位体积上的力，比如重力)。如果我们把一个材料划分成无限小的连续微粒，微粒之间的相互作用都是表面力(因为微粒在面上接触)。对于单个微粒而言，在静力学中，由于受力平衡，体积力为零；在动力学中，体积力正比于加速度。显然，为了表示材料的强度，应该使用表面力。而体积力可以应用在动力学上，比如在流体动力学中，动量方程的左边是体积力。

2. 应力是对外力的反应

如上图所示，一个外力(traction force)作用在任意物体表面上。要完整表示这个外力，必须既要描述这个力(大小和方向)和该力所作用的表面(面积和方向)。比如，只给定外力的大小和方向，而没有给定作用表面的面积和方向，该力的作用效果是不确定的(或者说定义不完整)。

作用在一个无限小表面上的外力，T，可以表示为:

外力T是一个二维张量，表示无限小表面上表面力的大小，单位是帕斯卡。

应力是对外力的反应。我们可以用同样的方法来定义内应力。假设如上图所示，一个物体被任意曲面剖开。在该剖面上任意无限小表面上的应力可以表示为

因此，应力可以理解为在物体内部平面上的内力的分布。要描述应力，必须要先定义一个平面。

3. 物体内一点的应力

应力是物体内部平面上的内力分布，通常用二阶张量来表示物体内一点(无限小的连续微粒)的应力。

在笛卡尔坐标系中，物体内一个无限小微粒的应力分布可以用如上图所示的无限小正方体来表示(在笛卡尔坐标系中，用正方体可以得到最简洁的形式)。在正方体的六个表面中，每个表面有三个相互正交的应力分量(一个正应力分量，，两个切应力分量，)。显然，和其他矢量一样，各分量的大小是和坐标系的选取有关的。

在该正方体表面上，实际有6×3=18个分量，因此这种用正方体表示应力的方法容易引起混淆。当这个无限小的正方体收缩成一个点时，该点的应力可以用下图来表示。通过该点有三个正交平面，每个平面上有三个分量，共有3×3=9个分量。这样，我们就得到了小正方体图中的矩阵：

在该矩阵中，代表正应力(垂直于作用平面)，代表切应力(平行于作用平面)，第一个下标代表应力作用平面的法线方向，第二个下标代表应力分量的方向。比如代表作用在x平面，y方向的切应力。

值得注意的是，上面的矩阵总是相对于对角线对称的。也就是说，，，。这是因为微粒表面上的力矩是平衡的，否则该微粒将会有无限大的旋转加速度。如上图所示的二维情况，如果，该微粒表面上的力相对于A点的力矩为(假设正方形边长为1，力矩的定义为)。要满足力矩平衡，必须有。

因此，应力可以由六个独立的分量表示：

至此，我们已经得到了最常见的应力的定义。

飞昂软件技术（上海）有限公司是美国ANSYS公司在中国的唯一全资子公司

“飞昂软件技术（上海）有限公司是美国ANSYS公司在中国的唯一全资子公司”，这是Fluent China网站最下方的一行小字。

同时网站上还有一则新闻“ANSYS中国独资子公司在上海揭牌”。看来Fluent中国已经加速投入ANSYS的怀抱了，但是另外一边还有代表ANSYS中国的安世亚太。

现在很难看出谁是有名无实，事件的发展还有待观察。