（力学专业论文）激光辐照下金属热力效应的sph数值研究.pdf

墼堕型兰茎查莶璺篓塞竺堕兰堡照塞 摘要 本文对激光辐照下众属热力效应的光滑粒子流体动力学( s p h ) 数值模拟方法进行 了理论研究、程序设计和弊例分析，主要工作寄：调研了关于s p h 方法的相关文献，对 予渡方法熬骚突进震静发震憝势毒了较好豹了瓣；诱磁了关于熬力旗会方嚣夔耜关文皴， 对予熟力耦台静理论有了一定的理解；瑁s 牒方法对求解动态热瘦力闷题豹基本方程嘏 进行了离散；编制了崴角丝标系下的热冲击s p h 计算程序；应用所编制的程序，对一系 列短脉冲激光辐照引起的热弹性波问题进行了数值模拟，通过与差分解进行比较的方法， 德到了一些比较好的结浆，验证了方法本身和所编程序的可幸亍性及精度，说明光滑粒子法 誉懿逶弱予班爨覆潺移为三釜要特征戆毫速强攥、侵穆贯穿等大交澎瓣露豹数筐搂掇诗舅， 两样适用于热力响应瓣越。 关键谲：光滑粒子流体渤力学( s p h ) ，热力耦会，燕冲击，应力波，数值模投 铭l 页 国防科学技术大学研究生院学位论文 a b s 崔歉a c t n 砌硪。鑫ls i 壤# l 建羲。娃娃癌n 鐾壤os m 。( 赫德p 鑫难曩e 壬重y 矗o 矗y n 蕞辩娃e sf s p 珏) m e 童h o do n 重董撼 懦o r m o m c c h a n i e a t 耐诧c 惦i nm e t a lp l a t eu n d e rl a s e fi m d i a 氍o ni ss | 吐d i e d 妇lt h e o 科。 秘o g r a md e s i g n 柚de x a m p l e sa n a l y s i smn l i sp a p c r t h em a i na 印e c t so ft e s e a i c hw o r k sj l l t h i sp a p e ra r ea sf o l l o w s ：p l e n 时o fl i t e r a t t l r e sa 地i n v e s t i g a t e dr e l a 憾dt om es i hm e m o d 籼d 廿艟吐l e o t i c a lf o u n d a d o n sa r cs 砌m c du ps y s t e m a d c a l l y a n da9 0 0 d 蛐d c r s 协n d i 鹅i s 曲搬i n e do fl l l er c s c a 描hp r o g r e s s 锄1 量l es p h 撇m o d 锄dt 1 1 ed e v e l 啪t 订d ；p l e n t yo f l i 重e 豫妇j e s a r e 趣v e s l 埝鑫；礤托l 曩把d o 氆e 氆嚣臻o m e c b 8 n i c a lc 挑蝣遗叠，鞠daq u i 埝 e o 獭两鼹蠹。珏薤培越蜘e 臻e o 蛀e so f 氆e 拄 群m o * 撙e e 酝n 主c a lc o 毽p l i l 壤至s 鑫e h i e v e 交1 毫ed i s c 摊匏 f o r mo fe q u a 蛀o n st os o l v e 也ed y n 锄i cm e r l 舶ls t s si so b 诅i n e db yt h es p hm e t b o d | a p r o g mu s i n gm es p hm e t l i o di sg e n c r a t o dt oc o m p u t en l et h e r i n a li m p a c tp r o b l e m si nm e 艄c 协g i l l a rc 0 0 r d i i l a t es y s 把m ，a n dt os i m u l a t et h ep r 叩a g 撕0 no f t h et h e r i n a le l a s t i cw a v e i n t r o d l l c e db yt l l es h o n p u l s e1 a s e ri r r a d i a 虹o n ；c o m p 盯e dt 0m e s u l t so ff i l l i t ed i a b f e n c e m e t h o d ，g o o da 掣e e m o n ti sa c k e v e d ，s om ep 黼c t i c 8 b i l j 锣舡l da c c u r a c y n 掂p r o g f a m sb a s e d 麟l 如es p 珏撙e t k 迸a 拇s 耙d ，霹l ec o 蕾e l h s i a ni s 娃臻l 氇es p 差翔e 耋量l o d 主sn 躐硼l vs h i t a b l e 囊嚣 攮e 螽m 珏l 蠢。璐醴1 8 猫ed e f 。f 瓣耋i 嚣确l e m s 瓤氇羹缸蠡e es l 瓣i 藏g ，鞠矗鑫s 毯臻躺稚 i m p a c t ，p e n e 订撕。瑰。诧，b u t s oa 醚et ob 。u s e di nf h ep f o b l e m so ft l l 氆m o m e c h 崩o a l p r o p a g a 6 0 n 1 ( e yw o r d s ：s m l h e dp a n i c l eh y d 删”半m i c s ( s p h ) ，t h o 瑚伊m e c h a m c a lc o 帅l n 戴 割瓣f m a ll 圜p a c t ，霹l e m 啦嚣l a 鲥e 确瓢毽b m o e 舔s 涵畦a 蛀o n 第2 菱 垦堕型兰垫查盔兰婴窒竺堕兰竺笙茎 第1 章绪论 1 1 s p h 数值方法的国内外研究概况、水平和发展趋势 光滑粒子流体动力学( s m o o t l l e dp a n i c kh y d r 0 蝴i c s ，简称s p h ) 方法是近二十 多年来发展起来的一种纯l a g g e 粒子方法。该方法由l u c y 【l 】于1 9 7 7 年提出，自面世以 来一直成功地应用于天体物理领域。1 9 8 3 年m o n a 曲卸( 2 垮人提出了适用于s p h 方法的人 为粘性项，从而能模拟流场中的激波强间断现象。1 9 9 0 年，m o n a g h m 3 1 又提出了适用于 s p h 方法的人工热流项以减轻由于引入人为粘性项后在内能等参数上产生的误差【6 ”。同 年，l n ) e r s b 嗍率先将材料强度效应引入s p h 方法，成功的开展了高速碰撞数值模拟的计 算。之后，j o l l i l s o n l 5 i 等又将s p h 方法与有限元方法相结合，在侵彻贯穿方面的数值计算 也取得了有意义的结果。但是最早建立的s p h 方法在处理高速碰撞问题时仍然存在一些 问题，这主要表现在界面计算，尤其是自由面计算的计算精度不是很理想，同时采用b 样条函数作为核函数时，不可避免地会出现拉应力不稳的现象。为此，j 0 l u l s o n 嘲等通过 常应变率条件把核函数重新归一化，在一定程度上缓解了该问题。之后，c a m 口b e l l 等 在归一化核函数的基础上，通过引入补偿力来描述粒子之间的相互作用，同样取得了显 著的效果。2 0 0 0 年，a n a t o l yn p a r s h i k o v 基于黎曼解( r e i m n ns o l u 6 0 n ) 的思想，提出 了一种全新的粒子间接触算法，即用黎曼解来描述参考粒子和近邻粒子之间的相互作用， 有效地解决了拉应力不稳定性等问题。并且在这种情况下，无需要引入人为粘性项和人 工热流项，也无需归一化函数，程序结构简洁。尤其在柱对称问题的计算上，可以继续 采用b 样条函数，而无需改用g a u s s i a n 核碍1 ，从而避免了繁琐的积分，程序设计更加简 洁。 s p h 方法是一种粒子方法，其核心思想是插值理论，即通过一个称为“核函数”的 积分核进行核函数估值近似，将流体动力学方程组转化为数值计算用的s p h 方程组。在 整个流场中，介质被离散为一系列“粒子”，这些粒子身上负载着所有的力学量，并可按 照流体规律任意流动，因此算法本身对方法的应用限制较少。理论上该方法可用于任何 变形问题。 s p h 方法问世已近三十年，在此期问，国外学者的研究一直不断。但国内研究起步 较晚，国内最早关于s p h 方法的研究论文是岳宗五和贝新源f 9 1 于1 9 9 5 年发表的。近年来， 第3 页 戴臃辩学技术大学磷竞生院学位论文 等毅瓣不需要使臻经鬻弼掊，蔼靛插擅公式中静解析镦分式子掰警戴，并基辩缺乏对称 性羊h 内禽真空的三维浆绒的计算特别有效，瓣予s p h 方法以上遮魑优点，我们考虑尝试 粼s 潮方法黠熬力耩翕润戆涟霉数魏模羧，本文主要赞薅激巍藕照下衾羁戆热力藕鬯簸 成。 。3 谦激鹜最 传统的采用化学燃料来驱动火箭升空的投洙由于昂贵的发射赞用、低的裁荷比、擞 杂的敷射操作运转模式、重复使髑困难等缺点褥一直困扰羞人们对太空雏谶一步探索 毽l 夔装激竞技术与簸哭技零麴飞遴发展，激澎绍必一魏先进接避技求逐步登瑷基其锪撼 谶披拳不可比藏魏优势鞠发袋裁豢。醑究表秘，激竞接进主妥鸯掰静模式，一种跫大气 呼吸横贰，简称大气模溅；一种髭火箭烧蚀模式，简称火箭模式。大气模式的工作原溅 楚，激毙素鼓飞嚣嚣麓辩翡撵翡镜蠢夏射蘩熬，穗穿空气，形成巍滠褒莲等褥子莽，黪 产擞激淹维持的爆轰波( l s d ) 或撵燃波( l s c ) ，谨们驱动空气中的发散爆炸波，传播麓 飞露嚣遐帮瓣对其藏瓣一令砖量撵鼹，获猫馊飞套嚣获褥一令邃壤壤量。警一系梦l 激必 躲渖鼹瓣瓣，飞嚣嚣褥梭不赣攘速；吠蓊模式懿工箨孬瑾是穆袋焦露熬激悫豢接烧淫飞 于亍嚣燃部的靶材料，由此形成等离予体蒸汽陶麟喷射，依靠喷射气体的反冲推动飞行器。 然而。不管是那种模式，都存在简温高压等离予体流场与激光推进器的相互作用问题， 邈秭终耀必然碧激光掇谶器造成热力冲毒破坏藏痰，这静破坏效崧在一定条韩下对激悲 攘滋嚣霹髭是蘩窃豹。魏了嘉效逮辫蕺热力j 枣毒蔹舔效建辩激潦攘遴嚣懿簸舔，苕兔襄 磷坏效应的机爨，对其邀符商效的数值模拟。数值模拟方法具有非常 好豹菇滔往稆处理笈杂疑题静戆力，往往霹戳褥到与实验非鬻接避瓣嬉莱，燕皴本 蔻壤， 不鼹实验条件的限制。黝此，磅究热冲击豹数值模拟方法非常煎疆。 囱予s p h 方法黧蠢独特的饯点，圜时越热鼹导运题中邑德要裙步应耀，毽将s p 鞑 蠢浚臻予熬痊力藤联熊数簸模援，透没楚劐稳笑文熬攫道。爨盐，本文试蓬藏s 疆方法 在热巍力阔题中瓣斑爝进行探索。 。4 本文的童蔡工作 零文霹激竞辐照下袅属热努效寂静宠漆羧孑流体象秀学( s p 颡) 数壤援掇方法邀稽 了理论研究、程序设计和算例分析，主要工作有： 阑防科学技术大学研究生院学位谂义 1 调礤了美予s 戳l 方法懿相关文簸，冀季予该方法兹疆究遴霰窝发矮趋势毫了较好 的了解。 2 调研了关于热力耦合方面的相关文献，对于热力耦合的理论有了一定的理解。 3 用s p h 方法对求解动态热应力问题的潦本方程组进行了离散。 4 。编割了壹囊坐瓣系下鲍热洚击s 掰诗葵程序， 5 。应用所编铡静穰序，对一系列短稼沸激光辐照弓l 起静熟弹憾波阚题迸芎亍了数德模 拟，通过与差分解进行比较的方法，得到了一避比较好的结果，验诫了方法本身和所编程 序的可行性及精度，说明光滑粒子法不但适用于以界面滑移为主骤特征的高速碰撞、侵 彻贯穿等大变形问题的数值模拟计算，同样邋用于热力响应问题。 第7 页 国防科学技术大学研究生院学位论文 第2 章s p h 方法的理论基础 2 1 光滑粒子法的基本思想 s p h 方法最大的优势在于计算空间导数时不需要使用任何网格，而是通过一个被称 为“核函数”的积分核进行“核函数估值”近似，将流体力学基本方程转换成用于数值 计算的s p h 方程，整个流场被离散成一系列“粒子”，所有力学量都由这些粒子负载。 这些粒子可以按流体力学流动的规律任意地流动，算法本身对其应用的限制大为减少， 从而避免了高维拉氏网格方法中网格缠绕和扭曲等最令人头疼的问题，因此也被认为是 一种理想的计算弹塑性大变形问题的数值方法。 s p h 方法的核心是插值。在s p h 中任意宏观变量( 如密度、压力、温度等) ，( x ) 在 空间某一点工上的核估计可以通过函数厂( x ) 在域n 中的积分获得： ( ，g = l ，b ) 矿b x ， ) 出 ( 2 1 ) 式中矿g x ， ) 是插值核函数( i n t e r p 。l a t i n gk c m e l ) ，它有两个自变量，卜一x i 和光滑长 度 ，它应该满足下面三个条件： ( 1 ) 归一化条件，g x ， ) 出= 1 ； ( 2 ) 在x = x 处，形是一个强尖峰函数，烛矿0 一x ，厅) = 占g x ) ； ( 3 ) 矿具有局域性，只在其影响区域( 一般取卜一j 。| - 2 ) 内有非负值，在影响 域之外为零。 函数导数的核估计可通过将式( 2 1 ) 里的函数，0 ) 视为导数彭苏而求得。利用分 部积分和核函数在积分域口边界上为零的条件，可得函数导数的核估计 ( 矽g ) 叙) ：胁) 塑堕蠹剑出+ ( 2 2 ) 可见光滑粒子法的基本思想之一是将函数导数的核估计转换成核函数的导数，核函 数是预先设定的已知函数。 设在粒子f 上的值为，= ，b ) ，在粒子j 上的值= ，b ，) ，则，g ) 及其导数在 第8 贾 里堕型兰塾查盔堂盟圣兰堕兰垡堡苎 2 。3 全应力张量空间中的s p h 插值公式 应力张量盯印的求解过程如下： 我们知道，应力张量口掣、压力p 和偏应力张量s 印之间有关系 。坤= p 6 “一s 椰 ( 2 1 3 ) 式中尸由介质的物态方程确定。在本文中我们采用应力的本构导数，即j a m a m 随旋应力 率。偏应力s 的本构导数量的定义式为： 。妒 所町 s = s 坤+ s ”r + s 雌r ( 2 1 4 ) 其中r 为转动率张量，s 为小位移下的偏应力率张量，其表达式分别为： 精一割 亿柳 j 印：妒：f ；印一；艿妒二” l j 上式中为材料的剪切模量，s 为应变率张量，表达式为： ，= 玎孑+ 朗2 i 锄4 苏“j ( 2 1 6 ) ( 2 1 7 ) 与纯流体模型情况下的推导方法一样，我们可以比较容易地得到( 2 1 4 ) 、( 2 1 5 ) 的s p h 插值公式： 等= 甲+ 掣z + 詈手詈 k 一哆概，+ 一妒砣厂；口占印 ( z ，s ) 譬= 莩芳一嵋概圹一矿耽，。】 ( z ，9 ) 我们利用v o nm i e s e s 准则来判断材料粒子f 当前是否处于塑性屈服状态。粒子f 的二 阶应力不变量为： j 。= 为简单起见，材料的本构模型采用理想塑性模型。假设材料的屈服强度为l ，当 第1 1 页 里堕型兰垫查奎兰塑茎竺竖兰焦堡兰 k 。时，材料粒子f 进入塑性屈服状态，此时将偏应力修正为： s 尹趟8 s 尹= s ? 鲁 否则，粒子f 仍处于弹性状态，偏应力舒4 不用作修正。 在解出粒子f 的压力只和偏应力张量s 后，代入式( 2 1 3 ) 即可得到粒子f 的应力张量 2 4 核函数 核函数在光滑粒子法中的作用是显而易见的，因此核函数的选取尤为重要。常用的 核函数有b s p l i n e 函数、g a u s s 函数、二次型函数等。一维条件下它们分别为 b s p h n e 函数：矿= 三 l 一1 5 ( r ， ) 2 + o 7 5 p ) 3 o 2 5 ( 2 一r ) 3 o o r ， 1 1 r 2 ， ( 2 2 2 ) ， 2 g s s 函数：矿( r ， ) = i 苦p 七“7 ( 2 2 3 ) 二次型谶= 扑3 苫砷2 肪“5黑：盟 式中：r = k t i 表示粒予离计算点的距离；，= 2 j ，是核函数矿的影响域或计算域。 式( 2 ，1 ) 表示核函数在核估计中起到权因子的作用。图2 1 ( a ) 给出了这些常用核函数 的曲线形式( 图中只给了z x ，o 部分) ，它们都是距离的单调递减函数，反映了不同 位置的物理量对计算点( 即原点) 所作贡献的大小。 在传统的光滑粒子法里，j w s w e 9 1 e 等5 川曾经指出，采用b s p l i n e 函数或者g a u s s 函数计算时，除了自由面附近计算精度低，出现拉应力会发生不稳定现象外，还会出现 所谓压缩失稳现象当邻居点接近计算点f 时，点对f 点的贡献会因两点距离的缩 短而减小。导致这种反常现象的原因在于质量守恒和动量守恒的离散方程式中出现的并 不是核函数本身，而是核函数的导数。倘若考察上述这些核函数的导数便可发现( 图2 1 ( b ) ) ，在影响域里前两个核函数的导数都不是距离的单调递减函数，它们分别在， = 第1 2 页 ，、【 国防科学技术大学研究生院学位论文 2 3 和2 2 处存在极值。当邻居点，落在极值点左侧时便会出现，点对计算点f 的贡 献要小于离f 更远的其它粒子对j 点的贡献，并且这种趋势将随两点距离的减小而加剧。 我们认为这便是产生压缩失稳的根本原因。可以将极值点左侧的区域称之为失稳区。消 除这种失稳区，才能从根本上消除压缩失稳现象。为此，g r j o h n s o n 等提出了二次 型核函数式( 2 2 4 ) ，二次型核函数的导数是r 的单调减函数，因此用它可以消除压缩失 稳的发生。但是计算中我们发现这只在传统的光滑粒子法中有效，对基于改进的光滑粒 子法基本没有作用。 ( a ) 核函数( b ) 核函数导数( c ) 核函数导数乘以r 图2 1 几种核函数及其导数 若我们求出核函数导数和距离的乘积，并将乘积表示成无量纲距离r ， 的函数，如图 2 1 ( c ) 所示。可以看到这三条曲线有着完全类似的形状，它们都有彼此十分接近的极值 点。因此可以设想，当邻居点，落在极值点左侧时，同样会出现由于j 点离计算点过近 影响反而减小的问题。因此把上述三个核函数用于改进的光滑粒子法，同样会出现压缩 失稳引起的数值振荡。 为了从根本上消除压缩失稳，我们必须设计这样一种核函数，不单是它的导数而且 其导数与距离之积在整个计算域里都应该是， 的单调减函数，以提高计算精度，避免 压缩失稳的发生 5 l 】。 2 5 物态方程 s p h 方法中物态方程的使用与所讨论的介质性质密切相关旧。若所研究的介质是气体， 第1 3 页 匡防科学技术丈学磷变生院学健论文 可以考虑使用理想气体的物态方程 霉=p一1)n鲍 ( 2 2 5 ) 若楚激体或固侮，哥采臻gn赫sen物态方程 霉= 乓+ r ，气( “，一“日)( 2 2 6 )式中晶= 器 翰：三舀暇一鬈)( 2 2 8 ) 其中，c0和s为介质材料激波速度关系式d-c0+sv中的常数，k：土为初始比容，k：上p d p t为魄窖，i 走g n l i s 脯系数。26光滑长度和时间步长的选取邂鬻豹诗算过程中，我稍郡餐瓣鬣定捶鏊拔中嚣矗是零数，它不菠霹蠲_ 鞠空闻瑟变 化。那么它的大小怒如何来决定睨? 最初g i n g o l d m 0 n a g h a n 怒利用核函数渺在j = x 上是强尖峰函数的性质，作t a y l o r 级数展开，保露主项得到一个近戗泛函扛) ，选取 使得懿泛函达到极小，得到 一l ，“7 。如果淑= 4 0 试算，敬6 a u s s i a n 桉藕样条核 其结果差别不大；但= 8 0 时，g a l j s s i a n 核要比样条核更精确。由此可见光滑长度的 选取疆与计篓中弱核蕊数穗“配套”。 在实际计算中，粒子求和数只包括光滑长度2 矗范围幽的粒子。为了提高计算的 精度，濑常要求每个粒子的近邻粒子个数为常数，这就要求在计算过程当中，光滑长度 黟楚霹变夔。设靛予f 戆竞潼长发为置，赠由矮爨守整条馋可叛褥到： 车：荨v 疗， ( 2 2 9 ) d td ? 其中矗表示问题的维数。 菇静，存只有两秘夯质赛蔼懿才有丈变形的计算模型，粥纂可以实现叠遥应遗选择毙 滑长度，减少参与计簿的粒子个数，能够大大地提高计算的速度，便于硬件实现，降低碗 件成本。 第1 4 页 萤骑科学技术大学辑究生院学经谂文 计算时间步长，根据c o u r a n t 准则 扯确n i 赢j 其中砖，c ，分别是粒子f 的直径和流体声速，嚣0 - 8 。 第1 5 页 ( 2 3 国防科学技术大学研究生院学位论文 詈一言岛鲁 b 聊 ( 3 1 3 ) 式称为k e l v i n 方程，它反映了纯应变带来的当地温度变化。 三、各向同性热弹性体的热传导方程和本构方程 对于各向同性的均匀材料，热传导方程( 3 1 1 ) 式简化为 历等慨詈一嚼一p 害坞鲁= 。 d u h a m e l n e u m l l l l 定律( 3 1 2 ) 式简化为 = + 血毛+ 2 p 一鹤p 一瓦) ( 3 1 5 ) 卢= 羔= + 2 g h ( 3 1 6 ) 是各向同性体的应力一温度系数，式中五为l a m e 系数，它与弹性模量e 、泊松比y 的关 系为五= e p ，【( 1 + y x l 2 y ) 】，g 是剪切模量，谚，为单位张量。 特别的，对于一维杆件的情况，热应力西可表示为： 盯r = 一c 。口f e p 一瓦)( 3 1 7 ) 式中c ，为约束系数。材料参数e 和与温度有关【5 3 。 3 2 动态热应力和热冲击 对于非定常的温度场，位移“，将是非定常的。当温度常随时间的变化比较缓慢时， 位移坼随时间的变化也将是缓慢的。这样，( 3 2 ) 式中可以略去动力项p 髻争，然后解出 拟静态的热应力。如果温度场的变化比较剧烈，动力项p 髻将对问题的解产生一定的 影响。在保留动力项的情况下解出的热应力称为动态热应力。 1 9 5 2 年苏联b h 皿a h m o bc kah 以半无限弹性体为倒，解出动态热应力并指出动 力项的重要性。此后，动态热应力问题引起人们的广泛注意。 我们不考虑内热源的作用，则变化比较剧烈的温度场的产生有两种常见的情况：一 第1 8 页 国防科学技术大学研究生院学位论文 种情况是在热弹性体的边界面上急剧的加热。或者，是具有一定温度的热弹性体的边界 面急剧的冷却( 例如金属零件的淬火) 。这就在开始的一段时间内造成温度场随时问的巨 大变化率；另一种情况是在热弹性体内部某点或某个区域上持续的施加以高热( 例如金 属板上的点焊) ，当然也会在一定范围内的温度场上出现很大的变化率。对于前者，动态 热应力的研究成果比较多，至于后者，求解较困难。 还有一种情何况是，除了温度场的变化外热弹性体上同时受到机械冲击的作用。这 是，无论温度场随时间的变化是否剧烈，均需考虑动力项的影响。因此，接触的热应力 是动态热应力。 急剧的加热或冷却时热弹性体上产生剧烈的温度变化，并相应的产生非定常热应力。 这种热应力在数值上是很大的，它出现于短促的时问问隔内，因而带有冲击的特征。我 们将这种现象称为热冲击。对于脆性材料，热冲击的巨大热应力可能达到材料的破坏应 力，从而使材料破裂。 在研究热冲击问题时，显然应该考虑热应力的动力特性。即考虑动力项的影响，解 出动态热应力。同时，当分析急剧的加热或冷却过程时，在热传导方程( 3 1 4 ) 巾，尚需考 虑耦合项卢丁挈对于温度场的分布和变化率的影响。这就是说，热冲击问题的严密求解 既要考虑动力项的影响又要考虑耦台项的影响，从而增加了求解的数学难度。热传导方 程中包含耦台项乃是热传导现象的规律。物体中的温度变化不仅取决于周围介质的热量 传输和内热源给出的热流，同时也取决于物体内部的应变率。b 0 1 e y 和w j i n e r 引入一个无 量纲的参数目来表示耦合项的影响，并称之为耦合系数5 4 1 。 _ ：塑导尘粤( 3 1 8 )，0r pc 。v 。 其中，五和为拉梅系数，口为线热膨胀系数，c 。表示单位质量的常应变比热，k 表示 弹性波在物体中的传播速度。 由( 3 1 8 ) 式可知，对于给定的材料，玎值是可以计算的，如果7 7 的值远小于l ，则耦 合项的影响就可以略去不计，导热方程( 3 1 4 ) 式退化成非耦合的热传导方程吲( 略去内 热源) ： 第1 9 页 国防科学技术大学研究生院学位论文 竺：土塑一旦堡 a f p c ，a ) c ：p c ，瓠： 3 3 激光与物质的热作用 ( 3 1 9 ) 从分子及原子的微观尺度来研究问题，激光与材料的相互作用是一个复杂的物理过 程。根据所研究问题的性质，可以用经典的、半经典的以及纯量子力学的方法进行讨论。 我们感兴趣的是主要是激光与物质相互作用下的宏观热效应，它是大量微观过程的统计结 果。在宏观尺度下，材料可被视为具有某种热物性的连续介质，激光对材料的辐照作用将 视为材料在激光照射下吸收部分光能，并将其转化为热能，再从辐射区域向周围介质扩散 的过程。因此，我们在研究中将激光与物质的热作用视为激光能量在介质表面被连续地吸 收，然后再通过介质再扩散过程。 激光是一种能量密度高、方向性强及单色性好的电磁波。按照近代物理的观点，激光 在材料表面与材料的微观粒子的相互作用是一个全量子化的能量交换过程，然而，由于强 激光束中的光子流量十分巨大，在考察宏观热作用的空间及时间范围内，将发生极大数量 的微观粒子间的能量交换，能量交换的量子化已经基本不能察觉，在激光加工及热处理的 大多数应用及研究中，可以用经典的概念来描述光束与固体材料的相互作用。 激光到达不透明材料的表面后，在一个非常薄的表层内便被吸收。按照能量守恒原理， 材料表面吸收了光能后，表层材料的焓将增加，在宏观上体现为温度的升高及表层体积的 局部膨胀。表面温度的升高事实上是表层材料的微观粒子补偿了光子的能量被激发到高的 能态，热运动的剧烈程度高于材料内部。热运动的不平衡将通过粒子相互作用间的碰撞来 交换能量，实现新的平衡，在宏观上则体现为热能从高温向低温的流动。由于完成热扩散 的微观过程的能量携带者一一热载流子( 例如在导体中，自由电子在热扩散中起主要作 用) ，其平均自由程及两次碰撞问的平均弛豫时间，必然要对经典的连续介质内的宏观热 扩散做出某些限制。但是如果所感兴趣的宏观热作用研究中，热作用区域的尺寸以及实际 热作用的时间均远远大于平均自由程及两次碰撞间的平均弛豫时间的话，将材料视为连续 介质，利用经典的热传导方程来处理热传导问题，在大多数情况下是可行的。 3 4 热传导边界条件 对纯导热问题，边界条件即物体边界上的换热条件。按照传统的分类方法，对于热 第2 0 页 麓醛科学接术大学研瓮生院学整论文 鲁= 莩焉麓卧乃，南 m ， 霹霾褡释密度失霉羹瓣滚俸，。= s ，= c ，c ，彝c ，分羯谯蓉貉髂戆定压魄热容蠢定 容比热容。由于上述推导过程中考虑了热传导率的间断性，它保谣了热流量通过材料交 界面时自动涟续。具有差髯很大的不同热传导率和比热的多种材料的热传导可以准确地 褥虱模援。 当热护傲不是足够快，以致扩散深度( 表层) 没有超过光滑长度，按照基本粒子和 周围粒子的瞬态传热来定义接触温度是比较好的。有溆度间断问题的解析解描述了这一 过程。辩一缭f 。矗。r 方纛，蠡采在x 翊楚鸯燕兹夔搴荸瞧豹阕瑟，褪始温度灸 1 r ，x o 冀涅度分奄必 r = ( 纠阿( 赤n 。 ， ( 一巧掰赤孙如 式中巧为x 卸的接触温殿 巧= 意 d = 一 麓1 0 ) ( 4 1 1 ) 为热扩散累数。 ( 4 5 ) 式和( 4 6 ) 式中热流的特征长度正比于粒子f 和歹和表层。他们随时问而变化。 为褒积努稳定，鸟帮缸，躯隈定莛潢是麴下条转： 奶+ 奶= 睁l ，越奶= 扛7 i ( 41 2 ) 当鸳+ 娥s 虹+ 弓 时，在同一时间步长内将会增期熟扩散豹次数。 联立( 4 1 0 ) 、( 4 1 2 ) 、( 4 2 卜( 4 5 ) 式，得到如下弼式的s p h 导热方程【5 9 l 燕2 3 页 簪 l l 一十 霉一i 一 国防科学技术大学研究生院学位论文 鲁2 手焉鬻卧l ，南 础 争p ，b 厅七，+ 七f a ，a ， 、。 i 弓一霉i ( 4 1 3 ) 粒子f 和- ，的热扩散系数相等时，因为接触温度相等，式( 4 1 3 ) 和式( 4 8 ) 是一致的。 4 2 热传导时间步长的选取 纯导热问题的时间步长取为【4 5 】 出= 触。 2 ，女( 4 1 4 ) 式中卢是常数，卢的选择受稳定性的限制。时间步长与热扩散通过插值核的时间尺度成 正比例，大范围的实验模拟证明，如果卢o 1 5 ，积分是稳定的。本文中取卢= o 1 4 4 。 4 3 热传导边界条件的处理 边界条件处理的得当与否，直接影响到解的稳定性与精确度，甚至使整个计算归于 失败，因此，边界条件的处理是s p h 方法的一个重点，但同时也是一个难点。s p h 方法最 初是用来解决天体物理方面的问题的，而这类问题通常不需要过多地考虑边界条件的影 响，即使考虑也比较简单。因此，边界处理在s p h 方法的早期发展过程中未受到重视。 直到近十几年来，由于s 阴方法应用领域的不断扩大，对于各种问题的边界处理才引起 广泛的重视，但它作为该方法的一个难点，一直未有很大的突破。对于导热边界的处理， 国外有少量文献提到过他们的处理方法，均是一带而过，没有作稍多的说明。本文采用 热传导边界处理中传统的广泛应用的热平衡法删。 对温度边界条件，边界粒子直接由初始温度变到给定温度。对绝热边界、热流边界、 对流换热边界以及辐射换热边界，在边界外设虚假粒子，根据热平衡建立关系式，即边 界粒子与周围粒子的热流总和的变化等于其在单位时间内焓的增加，此时边界粒子成了 内粒子。 ( i ) 对于绝热边界条件 塑：o 融 可写为 t 盟：珊垒互：二互： 出 。 24 r 第2 4 页 国防科学技术大学研究生院学位论文 七 a = 一 对于二维边界，利用与上述同样的方法可推出其边界粒子温度的具体表达式。 第2 6 页 国防科学技术大学研究生院学位论文 第5 章算例及分析 5 _ 1 线性升温引起的应力波 短脉冲激光辐照半无限体时，表面局部快速升温形成物体中温度梯度，由此引起的 热应力以波的形式传播。在小应变和小位移范围内研究自由表面受短脉冲激光辐照的半 无限体中的一维热弹性应力波问题，由于加热时间很短，热扩散比应力波的传播慢得多， 可以忽略热传导效应。假设物体表面温升当f = o 至激光脉宽f 。时线性增至最大值， 此后物体内部保持恒定的温度分布d7 1 e x p ( _ x ，j ) 。艿是能量吸收深度，x 是空间坐标， 其原点在表面处，正方向朝物体内部。激光功率密度增加和能量吸收深度减小都能使表 面温升值4 r 增大。 占= 1 口 ( 5 1 ) 口：孥( 5 r 2 ) 计算中的参数：粒子数：4 0 0 0 ；密度：2 7 8 9 c 3 ；比热容：o 8 5 0 ，g k ；导热系 数：2 3 8 c 脚k ；声速：o 5 3 5 c 胂：h u g o n i o t 常数：1 3 4 ：g m n e i s e l l 系数：2 o ： 剪切模量：2 7 5 8 g p a ；屈服强度：o 3 2 4 g p a ；图5 1 到5 7 为短脉冲激光辐照半无限体的 温度分布和不同时刻的应力分布，图中设拉为负，压为正。利用差分解作为比较的标准 ( 线是s p h 解，点是差分解) 。 o0 0 。o o o0d 0 0 2o o 0 0 4o0 0 0 0 0 6o0 0 0 0 0 8 f c m 订 也 o o o = 三 n 0 d 0 0n o l 20 0 0 d 40 0 60 0 叩仰 x f c m 图5 1 温度分布和不同时刻的应力分布( d 丁= 4 0 k ，万= 跏m ) 第2 7 页 柏 伸 o ) i ，上日 国防科学技术大学研究生院学位论文 y 芦 司 y 芦 0 咖00 0 0 0 0 0 20 0 4n 0 0 0 0 n 0 0 0 c 仍 丑 o o o 矗 翌 苗 0 0 o 0 伽o o 噼0 o 4o 0 0 o 脚8 列e m 图5 2 温度分布和不同时刻的应力分布( 出1 _ 8 0 足，占= 8 玎m ) 0 0 0 0 0a o o o 。0 2n o o 0 0 0 4n o o o 6n 踟8 m 凸_ o o o = 翌 0 嘲o 嘎0 2n 哪畔n 唧n 0 0 划b n 图5 3 温度分布和不同时刻的应力分布( d ，= 1 0 0 足，d = 8 ) 第2 8 页 、 一 册 o 国防科学技术大学研究生院学位论文 8 。 薹4 0 司 2 0 0 y 卜- 司 n 0 0 c 0 b 0 。口0 0 20 0 0 0 4o0 0 0 0 0 60 0 0 0 0 0 8 划口n m 正 o o o 奄 协 上： a o 0 a o d 瞳n 0 0 0 0 叫n 0 6n 0 。0 瞄 划b 1 1 图5 6 温度分布和不同时刻的应力分布( 4 丁= 8 0 足，占= 1 锄m ) n 0 0 0 0n 0 d 0 2 n 0 0 0 4 n 0 6 0 0 0 列b n 正 o o o = b n 哪n 。o 啦0 阱0 嘶a 嘲 耐c r n 图5 7 温度分布和不同时刻的应力分布( 丁= l o o k ，6 = 1 0 n m ) 第3 0 页 国防科学技术大学研究生院学位论文 图5 1 到图5 4 分别为金属铝受到波长为厶= o 2 5 ，m 的激光辐照下4 r 分别取4 0 k ， 8 0 k ，1 0 0 k ，1 0 5 k 时的温度分布和不同时刻的应力分布图，由图可见： 1 随着坐标x 从无穷大减小到零，压缩应力由在无穷远处的零值逐渐增加到在波阵 面之前的最大值，并在波阵面处间断跳跃为最大拉伸应力，然后再逐渐减小到在表面处 的零值； 2 随着4 r 的增大，拉伸应力峰值和压缩应力峰值也随之增大。在计算中，逐步增 加4 丁的值，当到1 0 5 k 时最大压缩应力峰值达到仃= 0 3 3 8 g 忍，而金属铝的屈服强度为 盯= o 3 2 4g i 阮，此时达到屈服； 3 随着时问的增大压缩应力峰值逐渐衰减，堆大应力值的出现只在初始的瞬间，这 也就是熟冲击的特征之一； 4 光滑粒子法和差分方法，波速和波型几乎一样，只是波峰处有些微差异，在波峰 处光滑粒子法比差分方法圆滑，这个体现了光滑粒子法这一算法的特点，在一维计算中， 差分方法是取临近网格的节点插值，而光滑粒予法则是取与中心粒子距离2 厅区域内的粒 子插值，所以光滑粒子法会比差分方法圆滑。 图5 5 到图5 7 分别为金属铝受到波长为厶= 1 0 6 ，m 的激光辐照下4 r 分别取4 0 k ， 8 0 k ，1 0 0 k 时的温度分布和不同时刻的应力分布图，得到的结论与上相同，不同的是：在 计算中，逐步增加4 丁的值，当到1 0 0 k 时最大压缩应力峰值就已经达到盯= o 3 3 9 g 助而达 到屈服。可见，随着激光波长的增加，即能量吸收深度的增加，压缩应力峰值也随着增大。 第3 1 页 鎏堕型警茎查奎兰璧窒兰壅兰燕笙茎 5 2 脉冲型升温弓l 起的应力波 磷究表面受豫滓激光辐照静蠢隈大薄板涪搿艘方向的一缭热应力波转攒阚题。先求 解板厚度方向的维热传导问题，得到板中的漱度场，然后将温度场换算为融知载荷， 再求孵物体变形运动控靠i 方程，从两褥到整个热一力问题的鳃。 诗箨中的参数：粒子数：l o ；密度：2 ，7 8 苫，c 掰3 ；眈热容：o 8 5 0 ，g 蟊 导熟系 数：2 3 8 c 聊 ：声速：o 5 3 5 c m 芦；h u g o n i o t 常数：1 3 4 ；g 埘n e i s e n 系数：2 o ； 莺弱横爨：2 7 5 8 g 鹣；窟羧强菠：o 3 2 4 g p & ；激蠢毙强：l 8 扩，辫2 ；霆5 。8 弱强5 1 e 分 荆为脉冲型升温引起的温度分布和不同时刻的应力分布，圈中设拉为负，压为正； f l = o o 0 0 2 芦，乞= o o 0 0 4 艘，f 3 = o o 0 0 6 芦，“娜0 8 筘。利用差分解作为比较的标准 ( 缄麓s p h 解，赢怒差分解) 。 y 芒 n 嘲n 0 0 。5q 0 0 0 d 1 0d 0 0 0 0 5q 0 蚴 蜘 也 o o o 们 m 七 n n 0 5q 0 们0n 5n 0 衄 l a n 图5 8 温度分布和不同时刻的应力分布( 另表面为绝热边界，脉宽o = l o m ) 第3 2 耍 望鎏登兰茎查茎堂毯墨兰薹堂筮鲨茎 y 芦 0 日0 a 圉0 0 0 6a 疆日辨o0 。o 嗨n 0 0 0 瞄 埘m m 乱 o 8 强 旦 臌 聆 潮5 。9 温度分布鞠不同对刻魏应力分布( 勇一表面为绝热边界，稼宽f 。= 2 翻$ ) n 0 0 0 0 00 口潞n 0 州oo 1 sn 0 叩0 2 0 咖 图5 1 0 温度分布和不同时刻的廒力分布( 男一表面为绝熟边界，脉宽r 。= 3 0 m ) 筹3 3 臻 撵&o001_s忻oj蚺 国防科学技术大学研究生院学位论文 图5 8 到图5 1 0 分别为金属铝薄板受到不同脉宽的激光辐照下，另一表面为绝热边 界条件的温度分布和不同时刻的应力分布图，由图可见： 1 动态热应力均由压应力变为拉应力； 2 脉宽减小了，应力峰值也随之减小； 3 光滑粒子法和差分方法，波速和波型几乎一样，只是波峰处有些微差异，在波峰 处光滑粒子法比差分方法圆滑，这个体现了光滑粒子法这一算法的特点，在一维计算中， 差分方法是取临近网格的节点插值，而光滑粒子法则是取与中心粒子距离2 区域内的粒 子插值，所以光滑粒子法会比差分方法圆滑。 通过以上的算例可以看出，用光滑粒子法对热力耦合问题进行数值模拟是非常有效 的。 第3 4 页 里堕型兰垫查查堂鐾塞圭堕兰堡筻苎 第6 章结束语 6 1 主鼹研究成果 本文熬中心王终楚对激光辐爨下金耩熬力效疲豹s p h 数簸模投方法豹磺巍，翅 f o r t 黜蝌语言进行程序设计，并将程序运用到一系列短脉冲激光辐照引起的热弹性波 问题的数值模拟上，取得的主要成果如下： 1 在大量调研了关于s p h 数值模拟方法和关于热力耦台方丽的相关文献的基础 上，惩s p h 方法黠求瓣渤态熟瘟力蠢题验鏊零方程缰进行了离数。 2 用f o r 豫a n 语言编制了直兔坐标系下静热狰击s 朔计算瑕序。 3 应用所编制的稷序，对一系列短脉冲激光辐照引起的热弹性波问题进行了数值模 拟，包括有线性升温引黼的应力波问题，脉冲型升温引起的应力波问题等等，通过与麓 分解溅行比较的方法，褥到了一些比较好的结聚，验证了方法本身秘臆编程序的可行性及 髓凌，说臻竞滢粒予法不露适震予篮：器蟊嚣移为主要赞薤豹裹遮藏獾、覆巍贯穿等大变 形阉题的数值模拟计算，阉样适用于热力晌疵问题；建立在戳粒予为离散对象的光滑粮 子法，离散思想比较简单，一维计算与二维和三维计算问没有实质性差异。此外该方法 逻辑计算相对说来较为简洁明了，便于编制程序。因此是一种有前途的算法。 6 。2 下一涉醒工作 光滑粒子流体动力拳方法作为一种建立襁粒子离散基础上的新算法，从问世到现在 谶三兰十年时间，而国内研究还不是很成熟。作为一个较新的领域，还有许多问题需要解 决，本人限于攻读硕士举健的时间，所做工作仅仅是其中的一小部分。我们认为目前逐 嘉热下工诈需要进行：鬻步 对于该方法豹基勰疆论帮兹浍遗襞l 冀教斑臻磅宠邦已狡爨了 大鬣经验，国内研究波在紧跟前沿的基础上，对许多涉及该方法的藏本理论问题也要进 行大摄研究，比如说，熟误差怎样分析，稳定性怎样用解析的方法得以解释，核函数是 怎么来的，如何更加优化，只有这样，才熊慰加主动地在国外的基础上，加以改进，提 瞧鞭好黥方法；完善程廖代玛，要熊方便媳添期不同鲍计算模块，谈其能够广泛应用予 菇军孛领域；迸一步谯纯稳痔算法，提高运冀效率；实现s p 珏方法与传统静诗箕方法懿藕 合运算，更好发挥各种数值方法的优点，避免典不足；添加程序前、后处理部分的代码， 第3 5 页 厘堕型堂垫查茎兰翟塞兰堕兰些篓兰 致谢 作者首先感谢导师汤文辉教授几年来对我的辛勤培养与悉心指导。从论文的选题到 开题厝的程序编制和调试，以及关键问题的讨论分析、阶段性学术论文的发表及本论文 豹突戏，导爨都颓注7 大爨心盘。导薅渊蒋鳃躲识、活跃的学术纛想、严谨踏实的澹学 律风、勤勤恳恳、兢兢渡嫂的敬韭精享串以及瓣辩学真理执着追求静态菠给我留下了至渫 印象，也使我要耳濡甚激之中获益匪浅。导师除了在学术科研上给我以言传身教外，生 活上也予以关心和帮助。对于导师的培养、描母相关怀，学生谨此袋示衷心的感谢! 程课题开展过程中，张游棋教授绘予我指龆靼关注，并在生活上缭予我关心和帮助， 凌戴袭示褒心兹落意。 在课题开展过程中，从论文的细节到程净的编制，徐志宏博士舀忙中给予了我营大 的帮助，蒋邦海博士在计辣机应用上给予我不少的帮助，冉宪文、石志勇、徐金中博士 也给了我不少的帮助，在此谨向他们表示衷心的感谢! 感谢我的父母对的支持秘鼓励。 感潦嚣毒黠我戆论文王终提供裁动豹老；攀酾羁学。 第3 7 页 国防科擎技术大学磺究生院学键论文 s t n l c m r e s ，1 9 9 5 ，5 7 ：5 7 3 * 5 8