ansys是融结构、流体、电场、磁场、声场分析于一体的大型通用有限元分析软件。在核工业、铁道、石油化工、航空航天、机械制造、能源、汽车交通、国防军工、电子、土木工程、造船、生物医学、轻工、地矿、水利、日用家电等领域有着广泛的应用。
最近需要解析一下ansys生成的文件，包括cdb和rst文件，但是没有学习过ansys方面的开发，意外中找到了Python关于解析ansys文件的包pyansys。国内目前没有关于此包的使用方法，借助此机会介绍Python通过pyansys解析ansys文件。
这是官方链接：https://akaszynski.github.io/pyansys/ansys_control.html
一、准备材料：
1. 安装pyansys
2. cdb文件，rst文件（有需要文末获取）
二、ANSYS APDL交互式控件
1. 非交互式绘图，需要安装ansys软件（目前来看需要15及其以上版本）。
在生成几何体和模数时绘制几何体和等点通常很有用，并且为了调试（或脚本）目的，在 pyansys 中绘制也很有用。要启用交互式绘图，在启动 ANSYS 时设置interactive_plotting=True。绘制命令（如 APLOT、EPLOT 和 KPLOT）将打开一个 matploblib。
import pyansys
# run ansys with interactive plotting enabledansys = pyansys.Mapdl(interactive_plotting=True)
# create a square area using keypointsansys.prep7()ansys.k(1, 0, 0, 0)ansys.k(2, 1, 0, 0)ansys.k(3, 1, 1, 0)ansys.k(4, 0, 1, 0)ansys.l(1, 2)ansys.l(2, 3)ansys.l(3, 4)ansys.l(4, 1)ansys.al(1, 2, 3, 4)
# sets the view to "isometric"ansys.view(1, 1, 1, 1)ansys.pnum('kp', 1) # enable keypoint numberingansys.pnum('line', 1) # enable line numbering
# each of these will create a matplotlib figure and pause executionansys.aplot()ansys.lplot()ansys.kplot()
运行完成后
2. 互动断点（需安装ansys软件）
在大多数情况下，尤其是在生成几何图形时，不打开APDL GUI就是不可能的。使用内联绘图很难识别几何图形项，因此pyansys有一种open_gui方法可以无缝打开GUI，而又不会丢失工作或重新启动会话，可以使用以下代码。
import pyansys
# run ansys with interactive plotting enabledansys = pyansys.Mapdl()
# create a square area using keypointsansys.prep7()ansys.k(1, 0, 0, 0)ansys.k(2, 1, 0, 0)ansys.k(3, 1, 1, 0)ansys.k(4, 0, 1, 0)ansys.l(1, 2)ansys.l(2, 3)ansys.l(3, 4)ansys.l(4, 1)ansys.al(1, 2, 3, 4)
# open up the guiansys.open_gui()
# it resumes where you left off...ansys.et(1, 'MESH200', 6)ansys.amesh('all')ansys.eplot()
这种方法避免了必须在交互式会话和脚本会话之间来回切换的麻烦。相反，可以进行一个脚本会话，并从脚本会话中打开GUI，而不会丢失工作或进度。此外，在GUI中所做的任何更改都不会影响脚本。可以在GUI中进行实验，python脚本不会受到影响。
3. 运行批处理（需要安装ansys软件）
可以定义运行ansys的函数，而不是通过使用输入文件调用ANSYS来运行ANSYS批处理。本示例基于具有扭转载荷的圆柱体的最大应力进行网格收敛研究。
import numpy as npimport pyansys
def cylinder_batch(elemsize, plot=False): """ Report the maximum von Mises stress of a Cantilever supported cylinder"""
# clear ansys.finish() ansys.clear()
# cylinder parameters radius = 2 h_tip = 2 height = 20 force = 100/radius pressure = force/(h_tip*2*np.pi*radius)
ansys.prep7() ansys.et(1, 186) ansys.et(2, 154) ansys.r(1) ansys.r(2)
# Aluminum properties (or something) ansys.mp('ex', 1, 10e6) ansys.mp('nuxy', 1, 0.3) ansys.mp('dens', 1, 0.1/386.1) ansys.mp('dens', 2, 0)
# Simple cylinder for i in range(4): ansys.cylind(radius, '', '', height, 90*(i-1), 90*i)
ansys.nummrg('kp')
# mesh cylinder ansys.lsel('s', 'loc', 'x', 0) ansys.lsel('r', 'loc', 'y', 0) ansys.lsel('r', 'loc', 'z', 0, height - h_tip) # ansys.lesize('all', elemsize*2) ansys.mshape(0) ansys.mshkey(1) ansys.esize(elemsize) ansys.allsel('all') ansys.vsweep('ALL') ansys.csys(1) ansys.asel('s', 'loc', 'z', '', height - h_tip + 0.0001) ansys.asel('r', 'loc', 'x', radius) ansys.local(11, 1) ansys.csys(0) ansys.aatt(2, 2, 2, 11) ansys.amesh('all') ansys.finish()
if plot: ansys.view(1, 1, 1, 1) ansys.eplot()
# new solution ansys.slashsolu() ansys.antype('static', 'new') ansys.eqslv('pcg', 1e-8)
# Apply tangential pressure ansys.esel('s', 'type', '', 2) ansys.sfe('all', 2, 'pres', '', pressure)
# Constrain bottom of cylinder/rod ansys.asel('s', 'loc', 'z', 0) ansys.nsla('s', 1)
ansys.d('all', 'all') ansys.allsel() ansys.psf('pres', '', 2) ansys.pbc('u', 1) ansys.solve() ansys.finish()
# access results using ANSYS object result = ansys.result
# to access the results you could have run: # resultfile = os.path.join(ansys.path, '%s.rst' % ansys.jobname) # result = pyansys.read_binary(result file)
# Get maximum von Mises stress at result 1 nodenum, stress = result.principal_nodal_stress(0) # 0 as it's zero based indexing
# von Mises stress is the last column # must be nanmax as the shell element stress is not recorded maxstress = np.nanmax(stress[:, -1])
# return number of nodes and max stress return nodenum.size, maxstress
# initialize ANSYSansys = pyansys.Mapdl(override=True, loglevel='error')
result_summ = []for elemsize in np.linspace(0.6, 0.15, 15): # run the batch and report the results nnode, maxstress = cylinder_batch(elemsize, plot=False) result_summ.append([nnode, maxstress]) print('Element size %f: %6d nodes and maximum vom Mises stress %f' % (elemsize, nnode, maxstress))
# Exit ANSYSansys.exit()
运行结果
Element size 0.600000: 9657 nodes and maximum vom Mises stress 142.623505Element size 0.567857: 10213 nodes and maximum vom Mises stress 142.697800Element size 0.535714: 10769 nodes and maximum vom Mises stress 142.766510Element size 0.503571: 14177 nodes and maximum vom Mises stress 142.585388Element size 0.471429: 18371 nodes and maximum vom Mises stress 142.825684Element size 0.439286: 19724 nodes and maximum vom Mises stress 142.841202Element size 0.407143: 21412 nodes and maximum vom Mises stress 142.945984Element size 0.375000: 33502 nodes and maximum vom Mises stress 142.913437Element size 0.342857: 37877 nodes and maximum vom Mises stress 143.033401Element size 0.310714: 59432 nodes and maximum vom Mises stress 143.328842Element size 0.278571: 69106 nodes and maximum vom Mises stress 143.176086Element size 0.246429: 110547 nodes and maximum vom Mises stress 143.499329Element size 0.214286: 142496 nodes and maximum vom Mises stress 143.559128Element size 0.182143: 211966 nodes and maximum vom Mises stress 143.953430Element size 0.150000: 412324 nodes and maximum vom Mises stress 144.275406
4. ANSYS对象方法
classpyansys.Mapdl(exec_file=None, run_location=None, jobname='file', nproc=2, override=False, loglevel='INFO', additional_switches='', start_timeout=120, interactive_plotting=False, log_broadcast=False, check_version=True, prefer_pexpect=True, log_apdl='w')
该类在后台打开ANSYS，并允许将命令传递给持久会话。
exec_file (str, optional)– ANSYS可执行文件的位置。保留为默认值None时将使用缓存的位置。run_location (str, optional) – ANSYS工作目录。默认为临时工作目录。jobname (str, optional) – ANSYS作业名。默认为'file'。nproc（int ，optional）–处理器数量。默认为2。override (bool, optional)–尝试删除run_location处的锁定文件。在先前的ANSYS会话过早退出并且未删除锁定文件时很有用。wait（bool ，optional）–为True时，等待初始化ANSYS之后再初始化python ansys对象。将此设置为False进行调试。loglevel（str ，optional）–设置将哪些消息打印到控制台。默认的“ INFO”打印所有ANSYS消息，“ WARNING”仅打印包含ANSYS警告的消息，“ ERROR”仅打印错误消息。Additional_switches（str ，optional）–用于ANSYS的其他开关，例如aa_r和学术研究许可证，将添加：Additional_switches =”-aa_r”避免添加-i -o或-b之类的开关，因为它们已经包含在启动ANSYS MAPDL服务器中。start_timeout（float ，optional）–引发ANSYS无法启动的错误之前要等待的时间。interactive_plotting (bool, optional)–使用启用交互式绘图matplotlib。默认为False。log_broadcast（bool ，optional）– ansys解决方案进度的其他日志记录。默认值True，在日志级别“ INFO”下可见。check_version（bool ，optional）–检查二进制文件的版本，并在无效时引发异常。preferred_pexpect（bool ，optional）–启用后，将避免在CORBA服务器模式下使用ansys APDL，并将生成一个进程并使用pexpect控制它。默认为False。log_apdl（str ，optional）–在当前ANSYS工作目录中打开一个APDL日志文件。默认为“ w”。设置为“ a”以追加到现有日志。
三、读写ANSYS档案（不需要安装ansys软件）
1. 阅读ANSYS档案
可以使用Archive加载包含实体元素（传统元素和现代元素）的ANSYS存档文件，然后将其转换为vtk对象：
import pyansysfrom pyansys import examples
# Sample *.cdbfilename = 'examples/sector.cdb'
# Read ansys archive filearchive = pyansys.Archive(filename)
# Print raw data from cdbfor key in archive.raw: print("%s : %s" % (key, archive.raw[key]))
# Create a vtk unstructured grid from the raw data and plot itgrid = archive.parse_vtk(force_linear=True)grid.plot(color='w', show_edges=True)
输出运行结果（三维可旋转）
此外还可以通过启用read_parameters参数来有选择地读入存档文件中的任何已存储参数。
import pyansysarchive = pyansys.Archive('examples/HexBeam.cdb')
# parameters are stored as a dictionaryprint(archive.raw.keys())
所有的参数以字典形式存在，上述代码中输出了所有的key值，对应可以获取value。
2. 编写ANSYS档案
使用VTK生成的非结构化网格可以转换为ANSYS APDL存档文件，并使用加载到任何版本的
import pyvista as pvfrom pyvista import examplesimport pyansys
grid = pv.UnstructuredGrid(examples.hexbeamfile)script_filename = '/tmp/grid.cdb'pyansys.save_as_archive(script_filename, grid)
# read in archive in ANSYS and generate cell shape quality reportansys = pyansys.ANSYS()ansys.cdread('db', script_filename)ansys.prep7()ansys.shpp('SUMM')
3. 使用ANSYS结果文件（rst）
ANSYS结果文件是FORTRAN格式的二进制文件，其中包含从ANSYS分析写入的结果。结果至少包含所分析模型的几何形状以及节点和单元结果。根据分析，这些结果可能是从模态位移到节点温度的任何结果。目前，仅支持以下结果：
# Load the reader from pyansysimport pyansysfrom pyansys import examples
# Sample result filerstfile = 'examples/file.rst'
# Create result object by loading the result fileresult = pyansys.read_binary(rstfile)
# Beam natural frequenciesprint(result.resultheader)result.plot_nodal_solution(0, 'x', label='Displacement')
其中result可以获取以下属性。
result.resultheader #结果标头
result.time_values #分析的时间和频率值
result.nnum #排序节点
result.enum #元素编号
result.geometry #几何内容
其中result.geometry包含以下键：
本期内容和BIM没有太大关系，但是也需要思考怎样使用Python将解析出来的结果附着到BIM模型中去，比如IFC或者常用的Revit软件。
微信公众号【BIM技术应用交流】回复“35”获取本期软件安装包和代码程序。
知乎专栏【BKM — 换个角度看BIM】，获取更多脑洞应用。
提升卡
置顶卡
沉默卡
喧嚣卡
变色卡
千斤顶
照妖镜