基于PANDA框架的非线性静力学有限元
- 格式:docx
- 大小:16.35 KB
- 文档页数:6
面向对象有限元并行计算框架PANDA的并行机制李于锋;张亚林【摘要】为实现面向对象有限元并行计算框架PANDA对高性能计算的支持,分别从并行计算流程、区域分解、分区信息和通信封装等部分设计PANDA框架在并行计算方面的数据结构.在计算流程中建立区域分解和并行求解器的配合协作方式,进而描述进行区域分割的3种网格剖分方法;对分区边界单元和节点信息的组织以及对并行通信操作的封装使复杂的并行通信调用简单、易行.用PANDA框架开发有限元程序时,只需以串行方式编写代码就可得到并行计算服务.对某夹具的算例表明PANDA框架并行计算实现的正确性和可靠性.【期刊名称】《计算机辅助工程》【年(卷),期】2011(020)001【总页数】6页(P24-28,41)【关键词】并行计算框架;PANDA;有限元;软件开发;并行机制【作者】李于锋;张亚林【作者单位】中国工程物理研究院计算机应用研究所,四川绵阳,621900;中国工程物理研究院计算机应用研究所,四川绵阳,621900【正文语种】中文【中图分类】TP319;O2460 引言有限元法是当今主流的数值计算方法,被广泛应用于工程领域的各个方面.大量高性能计算机的出现,特别是集群系统的广泛应用,使大规模复杂数值模拟得以实现,但必须采用并行技术进行程序开发.有限元并行计算框架开发的目的是使在复杂网格上进行的大型数值模拟程序开发变得更加快捷、高效.近年来,面向对象有限元技术[1-4]发展迅速,而开发框架却为数不多,且功能模块有较强耦合[5].为此,中国工程物理研究院总体工程研究所和计算机应用研究所联合开发有限元并行计算框架PANDA[6].框架中各组件松散耦合,可适应各种线性和非线性有限元法,不仅提供单元公式和本构模型的添加接口,也提供各种复杂边界条件的实现方法.基于消息传递在分布存储模型的并行机上开发并行程序,能显著控制数据分配和通信,这种控制又促进在大规模分布存储并行机上的高性能编程.[7]针对目前高性能集群的广泛应用,采用 MPI(Message Passing Interface)以分布存储方式实现有限元并行是正确的选择.面向对象技术对并行编程过程有重要影响(如对数据和函数的封装、继承和多态性),许多基于类等的编程思想都来源于串行编程领域,但也开始在并行编程中使用,尤其是恰当的抽象可隐藏那些使编程复杂的并行性.[7]以面向对象技术为基础设计的PANDA框架正体现该优点.1 框架的数据结构框架的主要数据结构有网格、场变量、单元组和通信数据管理等.网格是整个模拟计算区域的离散,支持线段、三角形、四边形、四面体和六面体等几何形状.一个计算区域可包含不同形状的网格子集.所有网格单元被划分为若干单元集,一个单元属于且仅属于一个单元集.场变量即自由度及其微分,框架以场的形式组织变量.任意一个自由度变量都归属于某个场,场通常定义在网格节点或单元上.根据时间积分算法,变量可能存储的数量不一,即场可能含有d阶微分,均要开设存储.所有的场均由节点管理器管理. 单元组是PANDA框架的算法核心类,实现一类有限元公式,在特定单元形状、形函数(含单元空间积分方式)、计算参考系(Euler或 Lagrange或ALE)和物理固有非线性(如接触)等情况下形成一种特殊的计算模式,即如何形成系统方程的单个单元贡献.每个单元组都会关联一个场,同时关联若干个单元集,负责遍历所施加的单元,进行单元矩阵和右端项的计算,并提供组装方程组的接口.通信数据管理基于网格剖分,所有单元集作为一个整体进行剖分.在采用镜像节点后,通信只发生在节点信息上.2 并行的实现策略PANDA框架基于MPI实现并行,具体策略为:首先,基于几何网格的区域分解,将子域映射到各进程;然后,每个进程在其子域上计算单元矩阵和单元右端项并提供组装接口,借助并行线性代数求解器实现大规模矩阵的分布式存储、组装接口以及方程求解,求解后,交换场数据边界信息;最后,进程在局部进行时间积分更新.所有的有限元相关计算都是分布式并行的,而大规模方程求解的主要并行任务交由线性求解器库负责.目前,已在PANDA框架中集成Aztec,SPOOLES,SuperLU,HYPRE[8],PETSc,Trilinos和PHG(由中国科学院科学与工程计算国家重点实验室开发)等著名求解器库.3 框架的并行计算流程框架并行设计的关键是数据变量的分布式存储和场数据时间更新算法的并行化.框架默认对所有的场变量基于网格区域分解进行并行通信.框架并行计算流程见图1. 图1 框架并行计算流程Fig.1 Parallel computation process of framework由图1可知,并行计算首先需进行区域分解,然后将局部模型信息读入框架,再进行场变量的编号(局部编号在进程内部确定,全局编号则通过全局通信确定).基于编号的局部和全局映射关系,每个子域将遍历其所存单元,进行单元矩阵和右端残力的计算和组装,并组装到在初始阶段就设置好的并行求解器中.有限元计算对每个自由度都形成一个代数方程,方程的编号很重要,框架采用的编号方式为:局部编号按节点进行,节点内按场顺序进行,场内按自由度顺序进行;全局编号按进程号排序,并保持局部编号的连续性.框架的并行是基于节点的剖分,所以按节点编号的方式对于并行处理十分方便.得到自由度编号和活动方程的编号(包括全局编号和局部编号对应关系),就可将从单元算法得到的单元矩阵和残量组装得到整体方程的矩阵和残量.通过并行求解器得解的本地部分.并行求解器主要完成A(u)=f的求解.作为一个独立的外部接口模块,其输入包含向量f(或者还有初始值u)的分布式存储、算子A的操作对u的影响方式.其中,可能涉及稀疏矩阵的分布式存储、A(u)的分布式计算和规约;输出为解向量u的分布式最终结果,其在误差允许范围内近似满足方程或报告不收敛状态.求得分布式存储后,下一步对场的时间积分算法更新.为计算内力,需从相邻分区获得外部影像节点场信息,需用到分区管理和通信管理.PANDA框架将分区和通信进行封装,使上层只需1~2个函数的调用即可完成外部信息的交换,设计灵活、巧妙.在进行场信息的更新以及内力的计算后,就可推进时间步进行下一次的循环迭代:根据算法选择进入局部单元遍历计算和组装,或只进行右端残量的遍历计算和组装,再到并行求解器求解,最后到场的更新和内力计算,循环往复.3.1 区域分解框架采用网格区域剖分的方法实现计算任务的分解.目前,框架网格分解以串行实现,且以节点分解为基础,进而对单元和通信作出安排.PANDA框架提供的分解方法有3种:图分解方式、索引分解方式和空间分解方式.图分解方式先将网格组织成图的形式,然后按照图的分解方法进行分解,最后将结果转变为网格单元和节点的分解;索引分解方式是按照模型读入的节点顺序依次分配给p个进程;空间分解方式是按照空间坐标进行箱子分解,每个箱子的尺寸都相同,先形成很多小箱子,然后按照箱子的邻接关系,借助图剖分将箱子分成p个组.由于有限元计算大多是以单元为单位进行的计算密集型作业,对于复杂几何构型的离散化,图剖分方式似乎是最好的选择.将网格组成图的方式由单元连接关系和节点场关联信息共同确定.场关联信息存在于一些复杂相依的约束条件(如边界条件或接触条件)中.目前,PANDA框架调用METIS[9]开源软件包进行多水平图剖分,也可采用其他方式分解.由图剖分得到网格节点的划分情况,然后需确定如何对有限元计算形成分区信息,用这种信息指导数据通信.3.2 分区信息框架采用类PartitionT管理所有的分区信息.首先明确的是哪些节点或单元信息需进行通信,哪些在本地局部范围不需要与其他分区通信.从区域分解得到节点-分区号映射后,结合图的表示进行节点的类别划分,再结合单元连接关系进行单元的类别划分.根据图和单元连接关系可将本地存储的节点和单元进行分类,形成3类节点和2类单元,见图2.图中网格被虚线分为2个子域,节点和单元的分类基于右边的子域.内部节点指没有与其他分区节点有连接关系的本分区所拥有的节点;边界节点指至少与某个其他分区节点有连接关系的本分区所拥有的节点;外部节点指与本分区边界节点有连接关系的其他分区的节点.内部单元特征是其所有节点均为内部节点,边界单元特征是至少一个节点为边界节点或外部节点.图2 并行计算中分区的单元和节点分类Fig.2 Node and element classificationin a partition of parallel computation在节点和单元分类后,就可以明确消息的来源和目的.为方便管理,采用分区对象处理这些信息.调用分区函数,得到分区信息的初始化,其中主要函数有4个:(1)Set()设置分区总数、ID号和本地3类节点,接收信息的分区号和节点号等;(2)SetOutgoing()设置发送信息的节点号;(3)InitElement Blocks()设置2类单元数组,单元映射和逆映射数组的主维数为单元块数目;(4)SetElements()按照单元连接关系将单元划分为2类,设置各单元数组及映射数组.通过这些函数可对分区对象进行必要的设置,另一种方法是读取分区文件.在设置完成后就可使用PartitionT对象,得到的重要信息包含本分区所有节点号、所有单元号、编号局部到全局的映射及其逆映射等.得到分区信息后,可根据需要写几何输出文件.设总分区数目为M,本分区号为m(0≤m<M),则几何分区文件分别为file.nM.pm.geom和file.nM.partm.在geom文件中指定的其他文件(如节点坐标文件、单元集子文件、节点集子文件)中的编号都是从1开始的局部编号.节点和单元的映射信息存储于partm文件中,其中存储的内部、边界、外部节点,进、出节点,内、外单元均为从0开始的局部编号,映射数组中的编号为全局编号.在框架中有限元管理器读入分区文件,构造本地分区PartitionT对象.通信管理器设置其成员fPartition.在并行模拟中,可通过PartitionT实现对通信消息的掌握,通过层次封装,上层只能得到通信管理器和有限元管理器掌握的节点的处理器位置信息.3.3 通信封装边界节点信息的传递是整个通信的核心任务.在图剖分的区域分解中,通过分区类PartitionT的NodesOut和NodesIn函数可获得本分区需互换信息的节点号.节点数据的分布在框架中影响以下重要的计算步骤:大型稀疏矩阵的分布式存储和场信息的更新.求解器的初始化函数设置方程的全局方程总数、局部方程总数和局部方程起始号;在单元贡献集成时需判断方程号是否在本地范围内,不在则不集成;场信息的更新由节点管理器负责,在其成员函数的实现中,调用通信管理器的AllGather函数(已经过封装,在图剖分情况下属于点到点的非阻塞通信).图3 通信的封装层次Fig.3 Hierarchy of communication encapsulation分区信息通过通信的封装影响系统方程的并行求解和场信息的更新.通信在框架中的重点类有CommunicatorT,MessageT和 CommManagerT,封装层次见图3,上层的类使用下层的类提供服务,应用程序可使用以下各层,但最方便的是使用通信管理类,其与分区之间的协作见图4.图4 通信管理类与分区之间的协作Fig.4 Relationship between communication classes and partition通信器类CommunicatorT是对标准MPI函数的封装.涉及的全局通信有全规约、收集、全收集、广播和同步等;局部通信有非阻塞和阻塞的发送与接收.开发者在需要其他数据通信时,可使用这种封装的通信函数进行.该方法体现模块性,提高错误处理的能力,但用户使用时需理解具体的通信细节.MessageT是对各种并行消息通信模式的封装,有 PointToPointT,AllGatherT和 CartesianShiftT等 3个派生类,分别完成点对点、全收集和笛卡儿网格移位通信,对应于图剖分、索引分解和空间分解等3种区域分解方式的通信. PointToPointT提供点对点的阻塞、非阻塞通信.fPartition成员提供剖分信息,其他成员包括非阻塞的接收缓冲区、发送缓冲区、接收和发送请求.成员函数Initialize()按分区的进出节点信息检查或分配接收和发送缓冲区,接收发送请求.AllGatherT()依据分区进行非阻塞的通信,先投递接收信息,然后对要发送信息的节点信息进行发送投递,接着等待处理接收到的消息并组装到欲通信的数组里,最后等待发送的完成.AllGatherT类提供全收集的通信功能,经过通信,每个进程都得到相同的全局数据场,主要服务于基于索引的区域分解方式.成员函数Initialize()设置本进程发送数据的长度以及确定本进程数据在目标数组中的位移.AllGatherT()执行全收集操作.应用于索引分解的全收集通信只是计算功能的并行,在数据空间上并没有分布存储,每个进程都存储全局数据.CartesianShiftT提供笛卡儿网格上的移位通信,主要成员为构造函数和AllGatherT().针对空间分解方式,在空间的方向和位移上交换数据.该类的主要成员有周围进程号,通信移位模式以及发送、接收信息的节点号.在实际应用中,可利用框架提供的3种消息模式进行通信,只需创建对象,初始化数据分布,然后调用AllGatherT()就可实现.CommManagerT是框架中对并行操作的顶层封装.它维护分区信息、管理消息列表,作为框架的重要组成部分负责节点场信息的通信.使用CommManagerT进行数值模拟物理量的通信非常方便:其将所有的通信消息管理起来,开发者只需提供欲通信的数组,在进行注册后就可通过简单的调用完成通信. CommManagerT的主要接口函数有Init_AllGatherT()和AllGatherT().前者根据区域分解方式初始化相应的消息通信方式,将分配好的消息追加存储到成员fCommunications数组;后者调用相应消息的通信函数(对于图剖分则调用点对点消息的AllGatherT()函数).典型的CommManagerT使用方式为nArray2DT <double> values(num_nodes,num_dofs);int id=fCommManagerT.Init_AllGatherT(values);…/* values assingment locally*/fCommManagerT.AllGatherT(id,values);fCommManagerT.Clear_AllGatherT(id);其中:nArray2DT为框架工具箱中二维数值数组模板类.使用CommManagerT很方便,隐藏具体的分区通信方式,易于上层应用,并行通信在框架中逐层抽象封装,应用程序可选择使用不同封装层次的类.4 算例用基于PANDA框架开发的有限元模拟程序对某夹具进行计算,采用开源并行预条件库HYPRE求解代数方程组.模型自由度约为10万,得到的位移云图见图5. 图5 某夹具位移云图Fig.5 Displacement clound picture of a clamp由图5可知,并行计算的结果与串行计算结果一致,并行计算中单计算节点的时间和内存消耗见表1,测试平台为某Linux集群,单计算节点处理器为Intel(R)Xeon(TM)CPU 3.40 GHz,内存4 GB.可知,基本呈现出线性加速比.表1 并行计算中单计算节点的时间和内存消耗Tab.1 Time and memory cost of single computation node in parallel computation计算节点个数1 2 4 8 t /s 327 157 79 37内存/MB 170 95 52 305 结论PANDA框架可实现区域分解、分区管理、通信封装、场的通信等的无缝集成、相互协作,完成工程计算的并行通信任务,屏蔽具体通信底层实现的复杂性,与外部线性代数求解器的接口设计巧妙,可配置多种求解器.基于PANDA框架开发有限元程序,开发者只需以串行的方式编写与物理相关的代码,就可得到并行的计算服务.在PANDA框架中利用封装好的CommManagerT进行自定义变量的通信也十分方便,甚至只需1~2个函数的调用就能实现多计算节点间复杂拓扑关系的通信. 参考文献:【相关文献】[1]BOSE A,CAREY G F.A class of data structures and object-oriented implementation for finite element methods on distributed memory systems[J].Comput Methods Appl Mech Eng,1999,171(1-2):109-121.[2]PATZAK B,BITTNAR Z.Design of object-oriented finite element code[J].Adv Eng Software,2001,32(10-11):759-767.[3]PANTALÉ O,CAPERAA S,RAKOTOMALALA R.Development of an object-oriented finite element program:application to metal-forming and impact simulations[J].J Comput Appl Math,2004,168(1-2):341-351.[4]COMMEND S,ZIMMERMANN T.Object-oriented nonlinear finite element programming:a Primer[J].Adv Eng Software,2001,32(8):611-628.[5]YU Lichao,KUMMAR A V.An object-oriented modular framework for implementing the finite element method[J].Computers & Structures,2001,79(9):919-928.[6]史光梅,何颖波,吴瑞安,等.面向对象有限元并行计算框架PANDA[J].计算机辅助工程,2010,19(4):8-14.SHI Guangmei,HE Yinbo,WU Ruian,et al.Object-oriented finite element parallel computation framework PANDA[J].Comput Aided Eng,2010,19(4):8-14.[7]DONGARRA J,FOSTER I,FOX G,et al.并行计算综论[M].莫则尧,陈军,曹小林,等,译.北京:电子工业出版社,2005:198-199.[8]FALGOUT R D,YANG U M.Hypre:a library of high performance preconditioners[C]//SLOOT P M A,TAN C J K,DONGARRA J J,et al.Proc ICCS 2002.Berlin:Spinger-Verlag,2002:632-641.[9]KARYPIS G,KUMAR V.METIS 4.0:unstructured graph partitioning and sparse matrix ordering system[R].Minneapolis:Univ Minnesota,1998.。
非线性有限元分析1 概述在科学技术领域内,对于许多力学问题和物理问题,人们已经得到了它们所应遵循的基本方程(常微分方程或偏微分方程)和相应的定解条件(边界条件)。
但能够用解析方法求出精确解的只是少数方程性质比较简单,并且几何形状相当规则的问题。
对于大多数工程实际问题,由于方程的某些特征的非线性性质,或由于求解区域的几何形状比较复杂,则不能得到解析的答案。
这类问题的解决通常有两种途径。
一是引入简化假设,将方程和几何边界简化为能够处理的情况,从而得到问题在简化状态下的解答。
但是这种方法只是在有限的情况下是可行的,因为过多的简化可能导致误差很大甚至是错误的解答。
因此人们多年来一直在致力于寻找和发展另一种求解途径和方法——数值解法。
特别是五十多年来,随着电子计算机的飞速发展和广泛应用,数值分析方法已成为求解科学技术问题的主要工具。
已经发展的数值分析方法可以分为两大类。
一类以有限差分法为代表,主要特点是直接求解基本方程和相应定解条件的近似解。
其具体解法是将求解区域划分为网格,然后在网格的结点上用差分方程来近似微分方程,当采用较多结点时,近似解的精度可以得到改善。
但是当用于求解几何形状复杂的问题时,有限差分法的精度将降低,甚至发生困难。
另一类数值分析方法是首先建立和原问题基本方程及相应定解条件相等效的积分提法,然后再建立近似解法并求解。
如果原问题的方程具有某些特定的性质,则它的等效积分提法可以归结为某个泛函的变分,相应的近似解法实际上就是求解泛函的驻值问题。
诸如里兹法,配点法,最小二乘法,伽辽金法,力矩法等都属于这一类方法。
但此类方法也只能局限于几何形状规则的问题,原因在于它们都是在整个求解区域上假设近似函数,因此,对于几何形状复杂的问题,不可能建立合乎要求的近似函数。
1960年,R.W.CLOUGH发表了有限单元法的第一篇文献“The Finite Element Method in Plane Stress Analysis”,这同时也标志着有限单元法(FEM)的问世。
基于PANDA框架的非线性静力学有限元论文导读:基于PANDA框架。
能够分析千万自由度规模的弹塑性静力学问题。
非线性求解策略。
形成了面向对象有限元并行计算框架PANDA。
并行计算,基于PANDA框架的非线性静力学有限元。
关键词:PANDA,静力学,非线性,有限元,并行计算1 引言特种武器结构复杂,在整个库存到靶序列(Stockpile to TargetSequence,STS)全寿命周期内要经历复杂严酷的载荷和环境条件,结构响应呈现出高度的材料非线性、边界非线性和几何非线性。
为提高特种武器的设计、试验和库存维护水平,对武器结构在各种条件下响应的精细建模和分析至关重要,需要充分考虑结构的几何细节和物理内涵,所建立的有限元模型可达上千万自由度规模乃至更高,而传统的商用有限元程序由于国外对我国的出口限制,非线性有限元模型的分析规模被限制在几百万自由度以下,且计算周期较长,无法快速响应设计和维护的需要。
为了提升特种武器的工程数值模拟能力,适应不断提高的武器工程数值模拟需求,迎接和加速由现阶段小规模低效率计算向大规模高效并行计算的转变,2007年中国工程物理研究院启动了院预研重大项目“武器工程大规模并行计算框架研究及基础平台开发”。
该项目在已有源码程序的基础上,通过在有限元并行计算方法方面开展研究与软件开发,初步形成了面向对象有限元并行计算框架PANDA,并基于PANDA框架初步开发了可应用于部分静力、振动、冲击和传热武器工程问题求解的大规模有限元并行计算模拟程序。
针对特种武器研制中的非线性静力学有限元大规模精细分析需求,充分消化吸收开放源代码的程序设计思想和技巧,基于PANDA框架,开发非线性静力学有限元分析所需的单元类型、材料模型、非线性并行求解策略,集成大规模线性方程组并行求解算法,初步形成了可求解小应变、有限应变线弹性和弹塑性静力学问题的非线性静力学程序。
悬臂梁弹塑性有限元分析模型达到了千万自由度规模,并行求解时间低于一小时。
基于PANDA框架的并行有限元模态分析程序开发和应用李健;郝志明;宁佐贵【摘要】为提高大型结构振动分析的规模、精度和效率,基于面向对象有限元并行计算框架PANDA和高性能矩阵特征问题并行求解算法,开发出适用于大规模结构振动问题计算的并行有限元模态分析程序;在超级计算机银河YH和曙光5000A上,通过不同算例验证该程序的正确性和可靠性.以某靶室结构为研究对象演示该程序的应用,指出实际应用时需注意加速收敛和谱变换等关键问题.该程序可作为大型振动工程计算的有效工具.【期刊名称】《计算机辅助工程》【年(卷),期】2011(020)001【总页数】4页(P29-32)【关键词】并行计算框架;有限元;软件开发;PANDA;模态分析;特征问题【作者】李健;郝志明;宁佐贵【作者单位】中国工程物理研究院总体工程研究所,四川绵阳,621900;中国工程物理研究院总体工程研究所,四川绵阳,621900;中国工程物理研究院总体工程研究所,四川绵阳,621900【正文语种】中文【中图分类】TB123;TP311;O2420 引言有限元模态分析技术的发展对结构振动分析至关重要.随着现代工程结构的日益复杂、工作环境愈加多样化,人们对结构振动分析的要求也越来越高.其关键和核心是实现大规模、高精度和高效率计算分析,以满足复杂结构设计中越来越高的可靠性和安全性要求.传统的串行计算机、串行版有限元软件以及节点数很少的并行有限元商用软件很难满足这种对大型复杂结构振动问题的分析要求,因此发展大规模有限元并行模态分析软件具有十分重要的意义.结构模态分析的数学基础是大型稀疏矩阵的广义特征值求解,属于当前国内外一个活跃的研究领域.[1-4]近年来,人们对矩阵特征值问题的并行计算研究比较深入,也提出许多行之有效的并行算法,并开发出一系列高性能特征值问题计算软件包.[5-6]然而,要实现大型复杂工程问题的并行有限元模态分析,仅有矩阵特征值问题并行求解算法还不够,必须将其与有限元分析方法结合起来.过去有较多的工作是基于串行版有限元分析软件的单元处理、整体矩阵组装及前后处理功能,利用并行特征值求解算法替代原串行特征值求解程序,实现模态分析的有限元并行计算.[7-9]其优点是工作量小、可充分借助成熟软件的分析能力,但采用串行方法进行有限元分析难以大幅提高模态分析问题的规模和效率.本文基于并行有限元程序框架PANDA[10]和高性能矩阵特征问题并行求解算法,开发并行有限元模态分析程序,通过不同算例验证程序的正确性和可靠性.1 矩阵广义特征值问题的并行求解算法本文重点分析科莱诺夫-舒尔算法,它是当前大型稀疏矩阵特征问题求解的最新、最具效率的算法之一,在Trilinos和SLEPc等著名软件包中均有集成.其特点是将正迭代投影算法、Arnoldi分解过程与重启动技术结合起来,而主要的并行工作为Arnoldi分解在各处理器上的并行执行.科莱诺夫-舒尔算法首先循环执行Arnoldi 过程,将稀疏矩阵化简为准三角矩阵;当Krylov子空间维数达到设定值时,如果Ritz对仍不完全收敛,则将收敛的Ritz向量保存起来,利用未收敛的Ritz向量中较好的一组线性组合构建新的初始迭代向量,重新执行Arnoldi过程,即将矩阵限制在已收敛特征向量张成空间的正交补空间上进行重新迭代;依次循环,直至全部希望的特征对收敛.有限元模态分析的广义特征值问题为式中:K和M分别为结构的总体刚度矩阵和质量矩阵;λ和x分别为特征值和特征向量.由于科莱诺夫-舒尔算法易于收敛至按模最大的高阶特征对,而在实际工程中却对低阶特征对的运用最为广泛.可对式(1)进行平移-逆谱变换,将原计算低阶特征对的特征问题转换为易于收敛的计算高阶特征对的新特征问题,即式中:μ为平移量;θ=1/(λ-μ).令 A=(K-μM)-1M,则式(2)转化为求矩阵A的标准特征值问题Ax=θx.假设求解n维矩阵A的k个特征对,Krylov 子空间维数限定为 m=k+p,m≪n,k<p,则科莱诺夫-舒尔算法的求解过程为(1)计算m步Arnoldi分解.选择范数值为1的初始向量v1(该处为2-范数),计算得式中:Vm满足 VTmVm=Im,是 Krylov子空间的基,Vm=Span(v1,Av1,…,Am-lv1);Hm 为上 H 矩阵,其下三角部分除次对角元素外均为0;f为残差向量.(2)对Hm作QR分解,得m×m阶的正交相似变换矩阵Q1.利用Ql对m阶的方阵Hm进行正交相似变换,将其转化为准三角矩阵,即舒尔型矩阵,即式中:Tm为准三角矩阵.式(3)可转化为(3)利用式(5),Ritz值可从分块对角矩阵Tm的主对角元素中获取.Ritz值分为Ωw 和Ωu这2部分,Ωw含有p个(p<m)中意的Ritz值,而Ωu含有m-p个不中意的Ritz值.对式(5)采用正交变换Q2可将中意的Ritz值移动到Tm主子矩阵中,得式中:m=VmQlQ2;Tw的全部 Ritz值λ(Tw)=Ωw,同理,λ(Tu)=Ωu;b*w为向量eTmQ1Q2前p个分量构成的向量.(4)判断p个(p<m)中意的Ritz对中是否存在k个满足精度要求的Ritz值,如果没有,则构建新的初始迭代向量,并对式(6)的分解进行截断,在第p+1步重新启动Arnoldi分解.(5)计算原问题的特征对和特征向量.2 PANDA框架及模态分析程序开发PANDA框架[10]是由中国工程物理研究院总体工程研究所和计算机应用研究所共同开发的一个并行有限元程序框架,主要瞄准并行自适应非结构网格上的多物理场数值模拟.PANDA框架采用区域分解算法和MPI技术实现有限元并行计算. 在有限元求解过程中,PANDA框架涉及到2个关键数据类:GlobalMatrixT类和SolverT类.其中,SolverT类为实际问题的求解提供具体的计算流程,并组集具有特定存储格式的描述系统方程总体矩阵;GlobalMatrixT类提供系统方程总体矩阵的存储格式.利用PANDA框架和科莱诺夫-舒尔算法开发并行模态分析程序的技术途径如下:(1)在GlobalMatrixT类基础上建立一种新的稀疏行压缩总体矩阵类型EpetraCRSMatrix,其4个一维成员数组 factive,frowptr,fcolind和 fnzval用于具体保存有限元分析中的总体矩阵.成员函数包括获取并行计算的通信域、将PANDA框架生成的总体有限元矩阵转换为特征问题求解器(如SLEPc)的矩阵格式等.(2)修改总体矩阵组集控制代码,生成类型为EpetraCRSMatrix的总体质量和刚度矩阵,并将其转换为特征问题求解器所需的并行矩阵格式.(3)编写特征问题求解的驱动代码.创建求解对象,构造特征问题方程,并提取特征算法所需参数,调用特征求解函数,最后将计算结果返回为PANDA框架数据格式.驱动代码的编写入口之一在SolverT类的CloseStep()函数体中.为利用外部求解器设置的各种运行时命令行参数,需对SolverT类的其他函数进行修改.(4)对外部求解器与PANDA 框架程序进行统编链接.利用上述方法,基于C++的继承、多态实现机制和Epetra对象模型,本文实现科莱诺夫-舒尔算法在PANDA框架中的集成,从而实现有限元并行模态分析程序的开发.3 并行有限元模态分析程序验证和工程实例分析3.1 有限元模态分析程序验证分别在万亿次银河YH(含64个双处理器节点)和曙光5000A(含16个八处理器节点)超级计算机上对开发的并行有限元模态分析程序进行统编链接.通过悬臂梁、锥壳结构算例考核程序计算结果的正确性,并与ANSYS计算结果对比,见表1和2.其中,表1对应的梁模型为402个单元、869个节点、2 466个自由度,采用10节点四面体单元;表2的锥壳模型为1 800个单元、12 900个节点、37 800个自由度,采用20节点等参数实体元.表1 PANDA框架与ANSYS对悬臂梁模态分析结果对比Tab.1 Comparison of modal analysis for cantilever beam using PANDA framework and ANSYS固有频率阶次/Hz PANDA框架 ANSYS 精确解(弯曲)1 58.26 58.59 57.49 2 69.56 69.47 3 351.48 355.95 360.34 4 412.58 416.86 5 542.06 543.74 6 893.40 896.85 7 939.42 960.18 8 1 080.60 1 097.30 1 009.10表2 PANDA框架与ANSYS对锥壳结构模态分析结果对比Tab.2 Comparison ofmodal analysis for cone-shaped shell structure using PANDA framework and ANSYS阶次固有频率/Hz PANDA框架ANSYS 1 127.82 128.33 2 167.38 169.36 3 297.28 297.78 4 299.24 304.09 5 388.55 390.37 6 414.11 418.70 7 473.20 483.32 8 580.42 589.52 9 599.74 600.93 10 671.48 676.40用梁模型测试该程序的并行分析能力.该模型共有200 000个六面体单元、220 941个节点、66万自由度.在曙光5000A上用不同数目的处理器计算前50个最小特征对所需的时间见表3.用单处理器计算时,由于内存限制,计算失败.可知,随着处理器个数的增加,计算的效率显著提高.表3 PANDA框架对66万自由度模型振动模态分析结果Tab.3 Results of modal analysis for the model with 660 000 degree of freedom using PANDA framework并行计算处理器个数计算时间/min 2 110.5 8 80.7 16 55.3 32 45.4 64 39.1 128 35.33.2 复杂结构模态并行分析工程应用实例某靶室结构的有限元模型共有39 014个节点、18 604个八节点实体单元,总自由度为116 934个,见图1.分析时,计算前200阶模态参数.在银河YH机上分析该靶室结构的振动模态,并与ANSYS分析结果比较,见表4.可知,本文开发的并行模态分析程序和ANSYS的分析结果非常接近.图1 某靶室结构有限元模型Fig.1 Finite element model of the target drone structure表4 靶室结构模态分析结果Tab.4 Modal analysis result for target drone structure固有频率阶次/Hz PANDA框架ANSYS 1 2.501 5 2.501 4 2 2.572 1 2.572 0 3 3.562 1 3.561 9 4 5.408 7 5.408 8 5 5.477 9 5.478 0 6 5.499 25.499 2 7 5.663 8 5.663 8 8 5.721 7 5.721 8 96.263 8 6.263 6 107.124 0 7.124 1另外,在实际应用时发现:(1)对具有密集模态的特征问题,计算的收敛速度较慢,甚至难以收敛;(2)当特征问题的矩阵形态较差(如矩阵元素存在数量级的差异时),低阶特征对迭代计算的收敛速度往往很慢,且计算结果的精度不高.为此,有必要采用谱变换方法将原计算低阶特征对的特征问题转换为易于收敛的计算高阶特征对的新特征问题.(3)新特征问题在求解时,需计算一个大型矩阵的逆矩阵与向量相乘的问题,该问题归结为求解一个大型线性系统方程,用迭代法求解该线性系统可能会失效,需采用直接法求解,这对内存容量提出较高的要求.4 结论(1)分析并行有限元程序框架PANDA和高性能矩阵特征问题并行求解算法原理;(2)开发一种并行有限元模态分析程序,并在超级计算机上验证程序的正确性和可靠性;(3)对复杂工程实例的具体应用表明该并行模态分析程序具有较高的性能.参考文献:【相关文献】[1]MITIN A V.New methods for calculations of the lowest eigenvalues of the real symmetric generalized eigenvalue problem[J].J Comput Phys,2000,161(2):653-667. [2]STEWART G W.A Krylov-Schur algorithm for large eigenproblems[J].SIAM J Matrix Anal& Applications,2001,23(3):601-614.[3]ZHOU Yunkai,SAAD Y.Block Krylov-Schur method for large symmetric eigenvalue problems[J].Numer Algorithms,2008,47(4):341-359.[4]魏立峰,李晓梅.计算实对称矩阵广义特征值问题的并行算法[J].计算机工程与应用,2001,37(11):4-5.WEI Lifeng,LI Xiaomei.Parallel algorithm for generalized eigenvalue problemof symmetric matrix pencils[J].Comput Eng & Application,2001,37(11):4-5.[5]张妮,曹建文.面向对象数值软件Trilinos及其线性代数包剖析[J].计算机工程与设计,2007,28(5):993-998.ZHANG Ni,CAO Jianwen.Research and analysis on object-orientednumerical software Trilinos and its linear algebra package[J].Comput Eng& Des,2007,28(5):993-998.[6]HERNANDEZ V,ROMAN J E,TOMAS A,et al.SLEPc:a scalable and flexible toolkit for the solution of eigenvalue problems[J].ACM Trans Math Software,2005,31(3):351-362.[7]张友良,冯夏庭.岩土工程百万以上自由度有限元并行计算[J].岩土力学,2007,28(4):684-688.ZHANG Youliang,FENG Xiating.Parallel computation of finite element model with millions of DOF in geotechnical engineering[J].Rock & Soil mech,2007,28(4):684-688.[8]李丽君,金先龙.有限元软件结构分析模块的并行开发及应用[J].上海交通大学学报,2004,38(8):1354-1357.LI Lijun,JIN Xianlong.Parallel development and application of FEA software structural analysis module[J].J Shanghai Jiaotong Univ,2004,38(8):1354-1357.[9]李丽君,金先龙.一种新型并行化有限元结构模态分析集成系统[J].计算力学学报,2004,21(5):546-550.LI Lijun,JIN Xianlong.A new parallel integration system of finite element modal analysis[J].Chin J Comput Mech,2004,21(5):546-550.[10]史光梅,何颖波,吴瑞安,等.面向对象有限元并行计算框架PANDA[J].计算机辅助工程,2010,19(4):8-14.SHI Guangmei,HE Yingbo,WU Ruian,et al.Object-oriented finite element parallel computation framework PANDA[J].Comput Aided Eng,2010,19(4):8-14.。
PANDA的功能可扩展性设计与应用孙乐;何颖波;莫军【摘要】针对有限元分析软件功能的可扩展性需求,介绍面向对象有限元并行计算框架PAN DA的单元设计和单元类的扩充方法.在PANDA中,单元分析被抽象成单元算法和材料本构模式,并提供单元类型、材料模型的接口服务,用户可以根据不同的应用需求构成并实例化不同的C++抽象基类来解决实际问题.单元类的扩展步骤为:在ElementListT类的单元对象列表中添加单元对象;为该单元类设置相应的属性和操作,以实现程序功能的扩充;修改编译文件以保证程序能够顺利地通过编译.以轴对称导热分析为例验证单元类的扩充方法在PANDA中的实现.PANDA的单元分析设计有力地保证其功能的可扩展性.【期刊名称】《计算机辅助工程》【年(卷),期】2011(020)003【总页数】5页(P55-59)【关键词】有限元并行计算框架;PANDA;软件开发;功能可扩展性;单元分析【作者】孙乐;何颖波;莫军【作者单位】中国工程物理研究院总体工程研究所,四川绵阳621900;中国工程物理研究院,四川绵阳621900;中国工程物理研究院总体工程研究所,四川绵阳621900【正文语种】中文【中图分类】TP311;0242在大型结构工程设计中,数值模拟是研究复杂系统综合性能的主要方法之一,而有限元分析软件在其中发挥着重要作用.近年来,随着航空、航天以及武器工程等领域的不断发展,出现很多大型和超大型结构,同时结构力学行为分析涉及的物理领域也不断扩展,使有限元分析软件逐渐面向大规模、多物理场的应用[1-3].将客观世界准确、简单而又自然地表达出来是计算机软件设计所追求的目标.良好的软件设计需要软件开发方法的配合,传统的结构化程序设计方法具有模块化特点,执行流程是各程序段根据给定的算法依次对数据进行处理,任何数据的添加或形式改变都会导致一系列相关过程的变化.随着分析问题类型的急剧扩大,有限元分析软件模块越来越多,模块之间的数据传递也越来越频繁,结构化程序设计导致程序的扩展能力有限,不易维护、调试复杂且代码的重用性低[4-6].从20世纪80年代后期兴起的面向对象的程序设计方法克服上述结构化设计的不足,其封装、继承和多态性特点为功能模块的集成化、扩充的灵活性提供有利条件,在一定程度上提高软件的开发效率,逐渐成为大型软件的主流设计方法[7-8].基于框架的现代软件开发理念和模式,贯彻软件模块化的思想,更注重沉淀和复用应用软件的共性技术,不仅可以实现代码级的重用,还可以重用设计与分析,在更大的粒度上实现系统级的复用,极大地降低软件研制难度、缩短软件研制周期,同时领域内的软件结构一致性好,确保软件的集成性以及可升级、可扩展性.[9]为面向复杂工程结构的大规模、多物理场计算分析应用,中国工程物理研究院于2007年采用面向对象、层次化、模块化的设计模式,形成面向对象有限元并行计算框架PANDA[10]的基本版本,见图1.该框架通过建立网格、场等底层数据结构,集成高性能数值求解器和工具包,如线性方程组求解器、特征值求解器、应用参数解析、数据输入输出、积分器、区域分割软件和并行通信等,以服务的方式提供一个快速建立高性能并行计算程序的开发环境;通过将有限元分析的数据和过程逐级封装,形成一个构架良好的有限元分析系统.在该系统中,有限元管理器以访问者的身份协调不同功能组件按照约定的接口协同工作,传递数据,完成有限元分析;系统设计采用继承和多态机制,易于集成或添加工程结构有限元计算所需的各种服务组件,如求解器、单元类型库和材料模型库及其接口服务,为软件开发提供一个扩展性良好的基础平台,支撑工程结构有限元高性能并行计算应用开发能力的不断扩充和完善.目前,已具备初步的结构静力学、振动力学、冲击力学、传热学以及热力耦合等并行计算程序开发能力.本文通过分析PANDA框架的单元分析设计,探讨其功能的可扩展性;并在现有单元库及单元接口服务的基础上给出添加二维轴对称导热单元类的实例,实现不修改原模块并在相同操作界面功能的有规则扩充,展示PANDA框架设计的优越性.单元分析是整个有限元分析过程的基础,尽管不同物理分析类型的单元分析算法各不相同,但其分析过程均可归结为依据单元的几何属性和材料属性等完成单元特性矩阵的计算以及系统特性矩阵的组装.[11]据此,PANDA框架将单元分析抽象成单元算法和材料本构模式两方面,利用面向对象的设计模式分别提供单元类型和材料模型的接口服务.单元(库)、材料(库)及其相关的接口服务处于PANDA框架顶层的应用个性层,与实际应用紧密相联;根据不同的应用需求,可以构成不同的C++抽象基类,通过对这些基类的实例化,可在PANDA框架上求解实际问题.PANDA框架的单元接口服务由单元基类ElementBaseT和单元对象构造器ElementListT共同完成(见图2).ElementBaseT主要负责向应用程序提供单元算法服务,它抽象出各单元类的共性(属性和操作),为单元算法提供统一界面,根据具体的单元算法可逐级派生得各种具体的单元类.ElementListT为单元类型的添加和识别设置统一接口,其主要功能为定义单元对象的列表;通过单元对象列表中的单元名称创建具体的单元对象;并通过一个单元类的指针指向已创建的单元对象,执行具体的单元算法.它以不改变程序主体结构的方式,实现向应用程序添加和识别新单元类的服务.PANDA框架的材料接口设计与此类似,不再赘述.面向多物理场应用的大型有限元软件涉及多种物理分析类型,因此在设计上需要考虑多物理场分析类型的可扩展性.该设计模式在面向多物理场应用时,可通过继承和多态机制在相同操作界面下方便地添加新的单元和材料类,并以不修改程序原模块的方式实现新的单元算法和材料本构模式,从而不断扩充已有程序的功能,有力地保障物理场分析功能的多样性、可扩展性以及代码的可重用性.从上文的分析可知,PANDA框架的单元分析设计有力地保证系统分析功能的可扩展性.在程序中添加新的单元类(添加材料类与其类似)需要完成以下3个步骤:在ElementListT类的单元对象列表中添加该单元对象,以便程序能识别该新增的单元类;为该单元类设计相应的属性和操作,实现程序功能的扩充;修改编译制导文件,保证程序能顺利通过编译.在实际工程的传热分析中,常需要采用轴对称算法简化整个分析过程,本文以轴对称热传导单元的添加为例,对上述3个步骤作详细说明.为最大限度地实现对代码的重用,分析PANDA框架的单元组织结构和已有的导热单元代码,分别以现有单元库中的线性傅里叶导热单元类DiffusionElementT,非线性傅里叶导热单元类NLDiffusionElementT和非傅里叶导热单元类HyperbolicDiffusionElementT作为父类,继承它们的属性和方法,再根据目标通过重定义的方式修改或添加某些方法,派生出各自对应的轴对称导热单元类,实现轴对称导热问题的有限元分析.详细的扩充过程如下:(1)在单元对象构造器中添加轴对称导热单元名称并定义新的单元类.①修改ElementListT.cpp文件,添加预定义宏命令,以便识别轴对称导热单元类DiffusionAxiElementT, NLDiffusionAxiElementT 和HyperbolicDiffusionAxiElementT的头文件.②修改函数 ElementListT::DefineInlineSub,在单元列表中添加轴对称导热单元的名称.③修改函数ElementListT::NewElement,按照单元列表中添加轴对称导热单元的名称创建新单元类的实例.(2)为轴对称导热单元类设计相应的属性和操作,实现程序功能的扩充.其中重要的成员函数如下:InternalEnergy计算内能;FormStiffness计算热传导矩阵;FormKd计算内力对载荷的贡献;LHSDriver计算热容及热传导矩阵;RHSDriver 计算内热源及热容项对载荷项的贡献;TractionBC_RHS计算热交换对载荷的贡献;TractionBC_LHS计算热交换对热传导矩阵的贡献;Axisymmetric判断是否为轴对称分析;B_axi计算轴对称情况下形函数对坐标的导数.(3)修改编译制导文件makefile.PANDA采用编译制导文件makefile进行分级管理.将添加的单元类 DiffusionAxiElementT,NLDiffusionAxiElementT和HyperbolicDiffusionAxiElementT的头文件与实现文件放入已创建好的单元目录panda\src\elements\continuum\diffusion下.修改该目录下的编译制导文件makefile,在其中加入新单元类的目标文件、源文件和头文件等信息.至此,完成PANDA框架单元库中轴对称导热单元的添加.在执行相应的应用程序时,有限元管理器将协调新的单元类及其他计算功能组件,使各功能组件通过约定的接口协同工作、传递数据,并最终完成轴对称热传导的有限元分析.为检验程序功能的正确性,采用如下算例进行测试:某双层空心圆筒的高度为8 m,内筒为厚度1 m的燃料层,外筒为厚度3 m的铝.燃料层有1 500 W/m3的内热源,热导率为35 W/(m·℃);铝层的热导率为100 W/(m·℃);铝筒外表面受到温度为150℃的高压水冷却,表面传热系数为3 500 W/(m2·℃).算例测试表明PANDA轴对称导热分析计算结果与商用有限元程序计算结果完全一致(见图3),验证程序的正确性.与传统的面向过程结构化程序设计方法相比,基于框架的先进开发模式和面向对象的主流设计思想可更好地适应程序对可重用性、可扩展性及复杂性的需要.PANDA框架采用面向对象、层次化、模块化的设计模式对有限元法单元分析过程的共性问题进行抽象,将其概括为单元算法和本构模型2个方面分别封装数据及操作;通过构建基类,为单元及本构计算分别建立统一界面;通过构建单元、材料对象构造器,为单元类型、材料模型的添加、识别设置统一接口.实例表明,该设计模式在面向多物理场应用时可在相同操作界面下方便地添加新的单元和材料类,并采用继承和多态机制以不修改程序原模块的方式实现新的单元算法和材料本构模式,从而不断扩充已有程序的功能,为多物理场的有限元分析程序开发提供一个扩展性良好的基础平台,大大简化软件开发过程,提高软件的开发效率.【相关文献】[1]Sandia National Laboratories.Computational activities in engineering sciences at Sandia National Laboratories,SAND2002-0392p[R].2002:5-6.[2]MCGLAUN M.Sandia’s engineering code development under ASCI,SAND2003-1238p[R].2003:9.[3]Sandia National Laboratories.Functional requirements for SIERRA version 1.0 beta,SAND1999-2587[R].1999:3.[4]李海江,杨刚,易南概.面向对象的串并行有限元分析系统[J].计算力学学报,2006,23(5):588-593.LI Haijiang,YANG Gang,YI Nangai.Object-oriented serial/parallel finite element analysis system[J].Chin J Comput Mech,2006,23(5):588-593.[5]曹骥,袁勇.面向对象有限元方法研究进展[J].力学季刊,2002,23(2):241-248.CAO Ji,YUAN Yong.Development of object-oriented finite element method[J].Chin Q Mech,2002,23(2):241-248.[6]马永其,冯伟.面向对象有限元程序研究综述[J].计算机应用研究,2001,18(10):7-9.MA Yongqi,FENG Wei.Advance in object-oriented finite element programming [J].Application Res Computers,2001,18(10):7-9.[7]魏泳涛,于建华,陈君楷.面向对象有限元程序设计——程序构架[J].四川大学学报:工程科学版,2001,33(4):21-25.WEI Yongtao,YU Jianhua,CHEN Junkai.Object-oriented approach to the finite elementprogramming:the application architecture[J].J Sichuan Univ:Eng Sci,2001,33(4):21-25.[8]张妮,曹建文.面向对象数值软件Trilinos及其线性代数包Epetra剖析[J].超级计算通讯,2006,4(1):36-45.ZHANG Ni,CAO Jianwen.Research and analysis on object-oriented numerical software Trilinos and its linear algebra package Epetra[J].Supercomputing Newsletter,2006,4(1):36-45.[9]Sandia National Laboratories.SIERRA framework version 3:core services theory and design,SAND2002-3616[R].2002:11-14.[10]史光梅,何颖波,吴瑞安,等.面向对象有限元并行计算框架PANDA[J].计算机辅助工程,2010,19(4):8-14.SHI Guangmei,HE Yingbo,WU Ruian,et al.Object-oriented finite element parallel computation framework PANDA[J].Comput Aided Eng,2010,19(4):8-14.[11]魏泳涛,于建华,陈君楷.面向对象有限元程序设计——单元过程设计[J].四川大学学报:工程科学版,2001,33(3):9-12.WEI Yongtao,YU Jianhua,CHEN Junkai.Object-oriented approach to the finite element programming:design of element procedure[J].J Sichuan Univ:Eng Sci,2001,33(3):9-12.。
基于PANDA的并行显式有限元程序开发
陈成军;柳阳;张元章;白小勇;何颖波
【期刊名称】《计算力学学报》
【年(卷),期】2011(028)0z1
【摘要】冲击响应数值模拟在军用与民用领域有广泛的应用需求,但现有串行分析程序难以满足求解规模和效率的需求,因此需要开发并行显式有限元程序.本文在分析了显式有限元串行基本算法的基础上,设计了相应的并行算法,采用模块化程序开发模式基于PANDA程序框架开发了并行显式有限元分析程序DynPack,并进行了算例验证.算例表明,DynPack实现了结构冲击响应的大规模并行计算,并具有较高的并行效率.
【总页数】5页(P204-207,214)
【作者】陈成军;柳阳;张元章;白小勇;何颖波
【作者单位】中国工程物理研究院,总体工程研究所,绵阳,621900;中国工程物理研究院,总体工程研究所,绵阳,621900;中国工程物理研究院,总体工程研究所,绵阳,621900;中国工程物理研究院,总体工程研究所,绵阳,621900;中国工程物理研究院,总体工程研究所,绵阳,621900
【正文语种】中文
【中图分类】TP311.O242
【相关文献】
1.基于显式有限元方法的二维楔形刚体入水砰击载荷并行计算预报 [J], 骆寒冰;吴景健;王珊;徐慧
2.基于PANDA的并行显式有限元程序开发 [J], 陈成军;柳阳;张元章;白小勇;何颖波
3.基于PANDA框架的并行有限元模态分析程序开发和应用 [J], 李健;郝志明;宁佐贵
4.基于有效并行求解策略的显式有限元分析并行算法 [J], 付朝江;王天奇;林悦荣
5.基于 AFEPack 软件包的并行显式有限元程序开发 [J], 郭永辉;浦锡锋;姚成宝;董楠
因版权原因,仅展示原文概要,查看原文内容请购买。
如何利用非线性有限元法进行力学分析非线性有限元法是一种用于数值分析问题的计算方法,其主要应用于力学分析领域。
这种方法在于其对于复杂结构的建模能力和高精度数值计算能力而备受推崇。
在本文中,将介绍如何对力学问题进行分析,以及如何应用非线性有限元法对力学分析进行模拟。
1. 引言力学分析整体上分为两种类型:静力学分析和动力学分析。
静力学分析研究对于物体的力和静止条件进行研究,其中力一般会造成物体的运动。
而动力学分析则研究运动物体的变化,特别是再一定条件下物体的振动问题等。
因为力学分析问题具有很高的复杂性,很多时候需要使用非线性有限元法来得到更准确的结果。
下面我们将详细介绍使用非线性有限元法进行力学分析的方法和流程。
2. 有限元法简介有限元法是一种现代数值计算方法,它将大工程结构分割为小的有限元。
在每个有限元内,结构的物理性质可以被认为是常量。
(具体内容可以自己百度)3. 如何利用非线性有限元法进行力学分析使用非线性有限元法进行力学分析的核心是将宏观问题转变为微观问题来进行模拟计算。
其中需要注意下面几点:3.1 确定力学分析的类型根据要进行分析的结构本身的性质和应用场景,可能涉及到静力学分析或者动力学分析。
其中静力学分析的计算主要涉及到结构在平衡状态下的情况,而动力学分析主要涉及到结构在某种条件下的运动和振动情况。
因此,在进行力学分析之前需要确定其类型,以便进行后续的计算。
3.2 建立结构模型根据具体情况,需要对结构进行建模。
建模可以通过一定的工具软件实现,或者手工建立结构模型。
模型的建立需要考虑到其复杂性和具体的应用场景。
构建好结构模型之后,需要对其进行精细化剖分得到单元网格,并进行编号。
3.3 确定边界条件在进行力学分析时,还需要考虑结构的边界条件。
边界条件可以通过指定某些点的坐标或者某些角度的变化来确定。
因此,在进行计算时需要根据具体情况设定边界条件,以便进行后续的计算。
3.4 进行数值模拟计算运用有限元法的基本原理,将每个单元的机械性质进行计算,根据力学分析的情况,可以得到结构节点的位移、应变和应力等参数。
基于PANDA框架的非线性静力学有限元
论文导读:基于PANDA框架。
能够分析千万自由度规模的弹塑性静力学问题。
非线性求解策略。
形成了面向对象有限元并行计算框架PANDA。
并行计算,基于PANDA框架的非线性静力学有限元。
关键词:PANDA,静力学,非线性,有限元,并行计算
1 引言
特种武器结构复杂,在整个库存到靶序列(Stockpile to TargetSequence,STS)全寿命周期内要经历复杂严酷的载荷和环境条件,结构响应呈现出高度的材料非线性、边界非线性和几何非线性。
为提高特种武器的设计、试验和库存维护水平,对武器结构在各种条件下响应的精细建模和分析至关重要,需要充分考虑结构的几何细节和物理内涵,所建立的有限元模型可达上千万自由度规模乃至更高,而传统的商用有限元程序由于国外对我国的出口限制,非线性有限元模型的分析规模被限制在几百万自由度以下,且计算周期较长,无法快速响应设计和维护的需要。
为了提升特种武器的工程数值模拟能力,适应不断提高的武器工程数值模拟需求,迎接和加速由现阶段小规模低效率计算向大规模高效并行计算的转变,2007年中国工程物理研究院启动了院预研重大项目“武器工程大规模并行计算框架研究及基础平台开发”。
该项目在已有源码程序的基础上,通过在有限元并行计算方法方面开展研究与软件开发,初步形成了面向对象有限元并行计算框架PANDA,并基于PANDA框架初步开发了可应用于部分静力、振动、冲击和传热武器工
程问题求解的大规模有限元并行计算模拟程序。
针对特种武器研制中的非线性静力学有限元大规模精细分析需求,充分消化吸收开放源代码的程序设计思想和技巧,基于PANDA框架,开发非线性静力学有限元分析所需的单元类型、材料模型、非线性并行求解策略,集成大规模线性方程组并行求解算法,初步形成了可求解小应变、有限应变线弹性和弹塑性静力学问题的非线性静力学程序。
悬臂梁弹塑性有限元分析模型达到了千万自由度规模,并行求解时间低于一小时。
本文介绍了基于PANDA框架的单元类型、材料模型、非线性求解策略设计,并初步验证了非线性静力学有限元并行计算程序的计算精度和千万自由度规模分析能力。
2 基于PANDA框架的非线性静力学有限元并行计算程序设计
通过中国工程物理研究院的预研重大项目,采用面向对象、层次化、组件化的设计思想,对工程结构非结构网格有限元分析程序的基本数据结构、并行通信、求解控制等方面的共性和可重用部分进行抽象和程序实现,并集成了区域分割、解法器等服务组件,形成了面向对象有限元并行计算框架PANDA,提供经过系统规划设计的应用程序开发接口,以提供服务的形式引导应用程序的设计和实现,初步建立了结构分析有限元并行计算应用程序的集成开发环境。
科技论文,并行计算。
基于PANDA框架,结构分析有限元并行计算应用程序的开发工作变得较为简单和高效,程序开发工作量大为减少。
在PANDA框架既设的应用软件架构下,应用程序开发者可以将精力集中到本应用程序独
有的个性部分,并充分利用框架中集成的经过充分验证的高效解法器等服务组件,场和网格数据的组织、存储和管理由框架负责,应用程序开发者无需关心其底层数据结构等实现细节。
对非线性静力学有限元并行计算程序的开发而言,通过继承PANDA框架中场、节点、单元类型、材料模型、空间积分器、求解控制等共性部分,开发适用于非线性静力学有限元计算的场、节点、单元类型、材料模型和非线性求解策略,使用由框架提供的输入参数解析、并行通信、区域分割、线性方程组解法器等组件,就可较快速形成可求解大规模非线性静力学问题的高性能有限元并行计算程序。
2.1 单元类型和材料模型
PANDA框架中的ElementBaseT基类抽象出了各种单元类型的共性(属性和操作),在应用程序中,具体单元类型都可由该基类逐级派生并添加自身特有的属性和操作得到,添加新单元类型较为方便。
目前,在PANDA非线性静力学有限元并行计算程序中已基于PANDA框架开发了结构分析中常用的八节点六面体小应变单元和更新拉格朗日有限应变单元,后者在每一载荷增量步对几何构形进行更新。
基于PANDA框架,应用程序中具体材料模型在各个类继承层次上定义,在程序结构设计时便于添加各种具体的材料模型。
目前已实现了适用于小应变单元和更新拉格朗日单元的二维、三维各向同性线弹性和弹塑性材料模型。
其中的弹塑性材料模型是线弹性模型与J2各向同性屈服条件的组合,塑性变形阶段的当前屈服应力由下述四种各向同性塑性硬化函数之一进行描述。
(1)线性各向同性硬化函数
(1)
上式中为初始屈服应力,为硬化模量
K=(EET)/(E-ET) (2)
其中,E为弹性模量,ET为割线模量。
科技论文,并行计算。
当硬化模量K为0时,采用线性各向同性硬化函数的小应变各向同性弹塑性材料模型就成为小应变理想弹塑性材料模型。
(2)linear-exponentialsaturation各向同性硬化函数
(3)
上式中为初始屈服应力,为线性硬化模量。
为saturation硬化模量,称为saturation应变,是一特征等效塑性应变量,用于描述等效塑性应变对saturation硬化的影响,的值越大,相同等效塑性应变时的当前屈服应力越小。
随着等效塑性应变值增大,linear-exponentialsaturation各向同性硬化函数逐步趋近于由表示的线性各向同性硬化函数(其初始屈服应力为)。
在ANSYS中,Nonlinear Isotropic Hardening (NLISO)材料模型对塑性硬化规律的描述与此相似。
(3)general cubic spline(立方样条)各向同性硬化函数
(4)
上式中,n为用于定义样条的点数(对应于样条的n-1个子区间和2个端点外区域)。
在每个子区间上,样条具有下面的形式:
(5)
上式中,。
从样条点计算系数a(i)的条件为:样条通过这些数据,且样条的一阶导数在子区间之间连续。
此外,样条的端条件被规定为:,以满足样条点数据范围之外的线性延伸。
科技论文,并行计算。
尽管不是必须和强制性的,人们一般选择,以便用定义初始屈服应力。
(4)power law(幂函数)各向同性硬化函数
(6)
由上式可知,初始屈服应力。
2.2 非线性求解策略
PANDA非线性静力学有限元并行计算程序实现了标准牛顿法、带线性搜索的非线性牛顿法、带线性搜索的非线性预处理共轭梯度法等三种策略。
它将整个载荷划分为多个载荷增量(可根据迭代求解过程的情况进行自动增减),在每个载荷增量内进行反复迭代求解直至满足某个收敛准则,而每一迭代过程中通过调用PANDA框架中集成的线性方程组解法器组件实现线性方程组的求解(线性方程组的求解可选用直接解法或迭代解法)。
三种求解策略使用最多17个参数对非线性求解的迭代过程、精度、自动载荷增量步、线性搜索、预处理过程等进行控制,都使用基于残余力范数的两个收敛评估准则(绝对容差和相对容差,后者是两次迭代的残余力范数的比值),只要满足其中一个准则,一次非线性迭代结束。
3 测试算例
为了初步验证基于PANDA框架开发的非线性静力学有限元并行计算
程序的计算精度、并行计算规模和计算速度,我们用该程序分别计算了杆、悬臂梁和及来自于工程实际问题的某离心机的静力学响应,并与理论解或商用有限元程序计算结果进行了对比。
对比结果表明,该并行程序具有很好的计算精度和计算速度,能够分析千万自由度规模的弹塑性静力学问题。