郭一諺 尹 飛 姜玉海 楊 林 李國飛 童 鶴

（1．中國兵器集團(tuán)第52 研究所煙臺事業(yè)部，山東 煙臺 264003；2．駐南京地區(qū)軍事代表局駐煙臺地區(qū)軍事代表室，山東煙臺 264003；3．陸軍裝備部駐包頭地區(qū)第一軍代室，內(nèi)蒙古 包頭 014033）

近場動力學(xué)是由美國Sandia國家實驗室的杰出研究員Stewart Silling 博士于2000 年提出的[1]。其理論是基于非局部作用思想的一整套力學(xué)理論體系，通過求解空間積分方程描述物質(zhì)力學(xué)行為，可以解決傳統(tǒng)基于連續(xù)性假設(shè)的網(wǎng)格劃分模型所面臨不連續(xù)問題的奇異性和復(fù)雜性，是當(dāng)前國際計算力學(xué)及相關(guān)領(lǐng)域的研究熱點(diǎn)。

1 “態(tài)”基近場動力學(xué)理論的背景

近場動力學(xué)理論以求微分—積分的解來解決傳統(tǒng)構(gòu)建模型所面臨的位移不連續(xù)（例如有裂紋的材料受力后斷裂）、應(yīng)變場不連續(xù)以及在不連續(xù)處存在導(dǎo)數(shù)沒有定義的問題。最初提出的中心力本構(gòu)模型的“鍵”基PD 模型有限制條件，即材料的等效泊松比只能為定值。為了解決這個限制問題，2007 年Silling 博士及其合作者提出了近場動力學(xué)的新本構(gòu)模型：“態(tài)”基PD 模型（State-based Peridynamic models）。用“態(tài)”替換經(jīng)典彈性理論中的“張量”，提出了一種同時兼顧連續(xù)（光滑）和不連續(xù)（不光滑）的空間映射?！皯B(tài)”的基本數(shù)學(xué)定義如下：m階的態(tài)為函數(shù)A〈·〉：H→Lm，因此，向量ξ∈H在態(tài)A的運(yùn)算下的映射就是m階的張量，記作A〈ξ〉。理論上講，“態(tài)”的數(shù)學(xué)定義可以到任意階張量[1-4]。

2 近場動力學(xué)的三類PD 模型

近場動力學(xué)理論是研究質(zhì)子與一定范圍內(nèi)的其他質(zhì)子之間互相作用的物理現(xiàn)象。按照向量值函數(shù)的不同定義可以分為“鍵”基近場動力學(xué)和“態(tài)”基近場動力學(xué)?！版I”基描述的是2 個質(zhì)子點(diǎn)變形前、后大小相等、方向相對的相互作用力矢量，這2 個質(zhì)子點(diǎn)間位移差允許位移場不連續(xù)，其相互作用取決于“鍵”的變形，因此泊松比被限制在一個獨(dú)立的材料常數(shù)?！皯B(tài)”基指近場動力學(xué)的力態(tài)，它與其向量函數(shù)的變形位移及周圍一定區(qū)域的物質(zhì)點(diǎn)都有統(tǒng)計學(xué)意義。

基于“鍵”和“態(tài)”的概念，Silling 博士等在2007 年定義了2 個質(zhì)子點(diǎn)之間的3 種不同本構(gòu)關(guān)系，即常規(guī)“態(tài)”基PD 模型（Odinaray State-based Peridynamic Model）、“鍵”基PD 模型（Bond-based Peridynamic Model）和非常規(guī)“態(tài)”基PD 模型（Nonordinary State-based Peridynamic Model）[2-5]（3 種基本的本構(gòu)關(guān)系模型如圖1 所示）。

圖1 3 種基本的本構(gòu)關(guān)系模型

常規(guī)“態(tài)”基的PD 模型是一個標(biāo)量態(tài)與一個向量態(tài)的乘積，力矢量的方向與“鍵”的方向平行。

“鍵”基的PD 模型只包括彈性常數(shù)，因此僅針對物質(zhì)的整體變形，無法區(qū)分剪切變形和體積膨脹，可以看作是線性變化的常規(guī)“態(tài)”基模型。

非常規(guī)“態(tài)”基的PD 模型是變形前物質(zhì)構(gòu)造的應(yīng)力與形狀張量的乘積，可用力態(tài)2 個端點(diǎn)各自的應(yīng)力張量表達(dá)（其力態(tài)方向和位移的形變方向可以不一致），且沒有材料泊松比的限制問題，更適合計算材料斷裂點(diǎn)及自然斷裂紋。

基于這3 種基本的本構(gòu)模型，Silling 博士針對裂紋、損傷過程的計算提出了材料質(zhì)子間力矢量的斷裂常用準(zhǔn)則，計算裂紋面上單位面積的等效斷裂功G0如公式（1）所示[2-6]。

式中：c為材料系數(shù)；s0為鍵伸長率的閥值（通常鍵的伸長率小于s0）；δ為物質(zhì)中1 個點(diǎn)的周圍領(lǐng)域的作用。

隨著質(zhì)子間相互作用的消失，質(zhì)子之間的力矢量斷開，連續(xù)質(zhì)子間的相互作用斷開就形成了裂紋（如圖2 所示），因此公式（1）可以作為模擬運(yùn)行過程中斷或截止的判斷準(zhǔn)則。

圖2 裂隙擴(kuò)展

3 近場動力學(xué)模擬計算彈靶作用的原理

荷蘭TenCate 公司提到，傳統(tǒng)的鋼質(zhì)組件可以抵抗質(zhì)量小于10 g、速度約為800 m/s 的步槍子彈，質(zhì)量相當(dāng)或更小的多層材料復(fù)合裝甲則可以對質(zhì)量約為65 g、速度超過1 000 m/s 的穿甲彈或重機(jī)槍子彈進(jìn)行防護(hù)。為了達(dá)到裝甲輕量化和最佳防護(hù)效果，復(fù)合裝甲材料的選擇完全依賴實彈試驗，費(fèi)用、時間成本消耗過大，還會受地域環(huán)境的影響。因此，現(xiàn)代多層材料復(fù)合設(shè)計所采用的方法是逐漸形成基于實驗研究的數(shù)據(jù)、歷史經(jīng)驗數(shù)據(jù)及其試驗結(jié)果，再結(jié)合模擬計算。常用的模擬本構(gòu)模型理論基礎(chǔ)是假設(shè)物質(zhì)的質(zhì)子位移變化是連續(xù)不斷的，如果出現(xiàn)不連續(xù)的情況，就在原網(wǎng)格的基礎(chǔ)上細(xì)劃局部網(wǎng)格，如圖3 所示。

圖3 有限元局部細(xì)化網(wǎng)格后的彈靶模型

復(fù)合材料在受到彈丸沖擊的過程中，從接觸到結(jié)束的時間很短，彈丸侵徹靶板的過程如圖4 所示。

圖4 彈丸侵徹靶板的過程

在整個沖擊過程中，先出現(xiàn)高速區(qū)，當(dāng)拉伸波產(chǎn)生的應(yīng)力超過材料的抗拉強(qiáng)度極限（Sc）時，材料易產(chǎn)生裂紋。隨著過程的進(jìn)行，彈丸的速度迅速下降，裂紋的擴(kuò)展速率也開始下降。直到應(yīng)力強(qiáng)度趨于0 MPa，裂紋的擴(kuò)展速率Sc也趨于0。傳統(tǒng)的有限元模擬采用局部網(wǎng)格細(xì)劃，無法確切模擬高速沖擊過程中的奇異點(diǎn)（此處不存在導(dǎo)數(shù)的解）。近場動力學(xué)理論采用的是微分-積分求解公式，能夠克服在斷裂處不存在導(dǎo)數(shù)的問題，其模型不再將材料中每個物質(zhì)粒的移動看作獨(dú)立的粒子移動，而是相對領(lǐng)域的粒子間形成的向量的相對位移量及向量兩端力態(tài)的乘積，在近場動力學(xué)理論的公式中使用的是位移而不是位移導(dǎo)數(shù)，采用可應(yīng)用于不連續(xù)體的空間積分方程。

最初，近場動力學(xué)理論的“鍵”基 PD 模型和“態(tài)”基PD 模型的理論都基于低速（v≤200 m/s）沖擊試驗，彈靶相互作用的沖擊破壞的過程快且復(fù)雜，模擬裂紋萌生和擴(kuò)展機(jī)理的難度很大。為了解決速度高于300 m/s 的沖擊問題（尤其是構(gòu)建沖擊模型的問題），2010 年Silling 博士引入“雙重狀態(tài)”的概念[2-6]，將“態(tài)”近場動力學(xué)理論擴(kuò)展到有限的區(qū)域內(nèi)：成對質(zhì)子相互作用和質(zhì)子之間的間接相互作用對其他質(zhì)子的影響，因此在分析復(fù)雜的多裂紋擴(kuò)展的過程中，可以模擬和分析具有復(fù)雜破壞機(jī)理和多裂紋共存的不連續(xù)問題。在實際應(yīng)用中，通過設(shè)置不同的材料屬性參數(shù)能夠更好地以積分形式重構(gòu)高速沖擊下，由各類同向性或異向性的材料粘合的連續(xù)體復(fù)合材料的運(yùn)動方程，可以獲得沖擊荷載下質(zhì)子跨界面運(yùn)動及跨界面區(qū)域中質(zhì)子的運(yùn)動軌跡的整個過程（連續(xù)的過程）。

在PD 模型本構(gòu)關(guān)系中，可以將材料看作是由無限物質(zhì)粒構(gòu)成的。材料發(fā)生形變有很多原因，此處的PD 模型中，當(dāng)單看材料應(yīng)力形變時，物質(zhì)粒變形前到變形后的位移與緊接的下次位移是斷裂開的，在相對領(lǐng)域內(nèi)，材料的形變可以看作是物質(zhì)粒間的非線性連續(xù)變化。基于非局部連續(xù)介質(zhì)力學(xué)的近場動力學(xué)理論， 單質(zhì)材料（例如鋼板）可以理想化為二維PD 模型，單層厚度不同“三明治”復(fù)合材料的每層都可以設(shè)定獨(dú)自的材料參數(shù)，當(dāng)理想化為二維結(jié)構(gòu)時，從厚度方向上就可以離散成單層的粒子，二維粒子模型（如圖5（a）所示），粒子間的間距為Δx。如果理想化為三維模型（如圖5（b）所示），那么其PD 的運(yùn)動模型可以將復(fù)合材料看作不同材料層內(nèi)單個領(lǐng)域（Hx）內(nèi)質(zhì)子的數(shù)量、位移的變化，在總體上可以看作多個領(lǐng)域Hx組成的整體，模擬的效果更接近試驗的預(yù)期結(jié)果。

圖5 PD 模型及非局部作用Hx 領(lǐng)域和子域

假設(shè)取整體復(fù)合材料中域H（xHx內(nèi)包括多個間距為Δx的質(zhì)子和作用域半徑為δ構(gòu)成的子域，如圖5（c）所示），隨著時間的變化，單質(zhì)子運(yùn)動產(chǎn)生的位移建立PD 本構(gòu)關(guān)系還需要考慮近場動力學(xué)理論中的運(yùn)動（平衡）方程，如公式（2）所示[3-9]。

式中：ρ為密度；u為位移；b為體力密度；T為位置和時間的函數(shù)；x為位置；t為時間；X和X'為相互作用的2 個質(zhì)子。

在公式（2）中，T[X，t]〈X'-X〉為在以質(zhì)子X為中心半徑δ的子域內(nèi)質(zhì)子X′對質(zhì)子X的作用；T[X，t]〈X'-X〉為在以質(zhì)子X′為中心半徑δ子域內(nèi)質(zhì)子X對質(zhì)子X′的作用。因此，PD 運(yùn)動方程的質(zhì)子拉格朗日方程如公式（3）所示。

式中：T[X（k），t]〈X（j）-X（k）〉為t時間X（k）的力密度矢量T（k）（j）；T[X（k），t]〈X（j）-X（k）〉為X（j）的力密度矢量T（j）（k）；v（j）為每個質(zhì)子的體積修正因子；b（k）為體積力密度矢量。

常見的“三明治”形式的多層復(fù)合材料由表面包覆層、阻燃膠、防彈層以及背板組成，如圖6 所示。

圖6 多層復(fù)合材料

PD 模型Hx領(lǐng)域的粒子化如圖7（a）所示。當(dāng)Hx領(lǐng)域受到高速沖擊時，粒子化的多層復(fù)合結(jié)構(gòu)的材料在PD模型中粒子可以PD 的質(zhì)子，其質(zhì)子運(yùn)動軌跡如圖7（b）所示。

圖7 材料粒子化

當(dāng)進(jìn)行PD 模型運(yùn)算時，復(fù)合結(jié)構(gòu)的材料的總厚度h如公式（4）所示。

式中：N為總層數(shù)；hn為第n層的厚度。

當(dāng)受到高速沖擊時，PD 模型的計算過程是將多層復(fù)合材料承載的外載荷以體力密度（單位體積的物質(zhì)所受的力）的形式施加在臨界材料邊界上的非零體積真實材料層上，對相鄰邊界層R1表面S1上的分布壓力p（x，t）或集中力P（t）的外載荷來說，它們的體力密度矢量b（x，t）如公式（5）所示。

最終，采用粒子離散方法建立近場動力學(xué)PD 模型，當(dāng)模擬受到?jīng)_擊時，受力粒子的運(yùn)動軌跡如圖7（b）所示，圖中是由3 種不同物質(zhì)材料復(fù)合而成“三明治”結(jié)構(gòu)。

在近場動力學(xué)中，材料損傷是通過質(zhì)點(diǎn)之間的相互作用（微勢能）的截斷來描述的。當(dāng)2 個質(zhì)點(diǎn)x（k）和x（j）之間的伸長率s（k）（j）超過它受力的臨界值sc時，就會產(chǎn)生損傷。當(dāng)出現(xiàn)損傷時，運(yùn)動方程中2 個質(zhì)點(diǎn)之間的作用力將不可逆的永久消失。質(zhì)子對受沖擊的裂隙尖的形成、裂縫的擴(kuò)展以及邊界效應(yīng)、材料崩落的形式來說，非連續(xù)的物質(zhì)的質(zhì)子位移向量和端點(diǎn)矢量力都是動態(tài)的、連續(xù)變化的，避免了基于連續(xù)性假設(shè)建模和求解空間微分方程的傳統(tǒng)宏觀方法在面臨不連續(xù)問題時的奇異性。因此，當(dāng)通過計算收斂速度估算近場動力學(xué)微分時，可以用公式（1）作為模擬計算停止的判斷準(zhǔn)則。

模擬計算中近場動力學(xué)理論的運(yùn)動（平衡）方程可以通過Taylor 級數(shù)展開近場動力學(xué)函數(shù)gp N（ξ），如公式（6）所示。

式中：ap q為a 展開的級數(shù)為p、項數(shù)q（a為多項式展開中的系數(shù)）；ω（ξ）為權(quán)函數(shù)；ξq為ξ展開q項。

權(quán)函數(shù)ω（ξ）指影響半徑內(nèi)控制相互作用的程度。根據(jù)Δx，權(quán)函數(shù)被指定為ω（ξ）=Ae-（ξ/δ）2或ω（ξ）=1，其中ξ=|x（j）-x（k）|，影響半徑通常被指定為δ=mΔx（δ為影響半徑）。

影響動力學(xué)微分與局部精確微分的偏差的因素除了權(quán)函數(shù)ω（ξ）外，還有其相互作用的域Hx、近場動力學(xué)函數(shù)gNp（ξ）的項數(shù)以及Taylor 級數(shù)展開式的剩余項R（N，x）。

在實彈試驗的過程中，隨著海拔高度的變化，大氣壓力、溫度對彈道參數(shù)影響較大，對影響試驗結(jié)果與依靠射表計算的預(yù)期結(jié)果之間的偏差影響的因素如下： 1） 空氣阻力Rx（如公式（7）所示）。2） 阻力系數(shù)Cx0（M）（如公式（8）所示）。

式中：ρ為氣體（或流體）密度；v為氣流速度；S為彈丸橫截面積；M為馬赫數(shù)；Cxj為模組系數(shù)；Cxb為渦阻系數(shù)（或地主系數(shù)）；Cxw為波阻系數(shù)。

與槍彈相比，當(dāng)炮彈飛行時，其空氣在彈丸表面形成附面層，由外到內(nèi)分別為邊界層、層流附面層和紊流附面層，計算過程還應(yīng)考慮影響炮彈彈速、著彈點(diǎn)的因素——雷諾數(shù)（附面層由層流向紊流轉(zhuǎn)變的因子），如公式（9）所示。

式中：l為彈長，有時也可以用彈丸的直徑表示；η為氣體的黏性系數(shù)。

空氣的黏性系數(shù)與高度或溫度的經(jīng)驗公式如公式（10）所示。

式中：τ為按溫度隨高度分布的標(biāo)準(zhǔn)定律，τ=τ（y）（y為高度）。

為了簡化計算過程，海拔高度的變化對彈道參數(shù)的影響因素作為約束條件可以設(shè)置為常數(shù)，也可以設(shè)置為輸入項。當(dāng)應(yīng)用微分算子運(yùn)算規(guī)則求解微分-積分時，可以放寬連續(xù)性假設(shè)條件，近場動力學(xué)的微分算子對函數(shù)進(jìn)行泰勒展開，其等價公式如公式（11）所示。

式中：f（x）和f（x+ξ）為求微分的目標(biāo)函數(shù)（ξ是x的增量）；p為階數(shù)；g為近場動力學(xué)函數(shù)；f為求微分的目標(biāo)函數(shù)；gNp為函數(shù)g展開的級數(shù)為p、項數(shù)為N。

公式（11）可以把具有p階連續(xù)性要求的目標(biāo)函數(shù)f（x）的微分運(yùn)算等效為將目標(biāo)函數(shù)f（x+ξ）乘以一個特定的近場動力學(xué)函數(shù)g的積分運(yùn)算來求解，利用這種解析方法來求解PD 運(yùn)動方程的質(zhì)子拉格朗日方程（公式（3））。

隨著中印邊界問題的凸顯，對武器裝備在不同海拔的實際防護(hù)能力的研究日益迫切。在不同海拔高度的實彈試驗過程中，會受到許多條件的限制（隨著海拔高度的變化，影響彈丸速度、著靶和威力等性能的條件（例如大氣壓、溫度等）也會發(fā)生變化）。依靠現(xiàn)代計算機(jī)運(yùn)算的強(qiáng)大功能，合理選擇、利用數(shù)值計算的理論方法，根據(jù)給出的測度值和需求條件計算，高速沖擊對目標(biāo)體的運(yùn)動過程及結(jié)果，與選定的特定海拔試驗得出的數(shù)據(jù)進(jìn)行對比，可以得到更準(zhǔn)確的研究資料。

利用近場動力理論對作用域的積分求解可以避免有限元方法對導(dǎo)數(shù)求解的奇異性，可以更準(zhǔn)確地模擬高速沖擊對目標(biāo)體的運(yùn)動過程。計算近場動力學(xué)代碼包括時間步長、最小時間步長、遞增時間、迭代次數(shù)、初始條件邊界條件的設(shè)定、域內(nèi)面力以及界面修正補(bǔ)償（多種復(fù)合材料粒子化材料邊界的非完整粒子的修正）的非線性的求解。

二維PD 模型計算量與同材料的有限元的計算量近似，在裂紋自然形成的模擬計算中，二維PD 模型比有限元模型精確，其原因是2 種模型的計算方法不同。

三維PD 模型計算量非常大，但是非常態(tài)基的近場動力學(xué)分模塊化計算的方式與現(xiàn)代計算機(jī)的并行結(jié)構(gòu)更匹配，PD 模型的可分解成多模塊的結(jié)構(gòu)適用于多點(diǎn)觸發(fā)的并行計算模式（多路處理器并行處理的計算機(jī)），與二維PD模型相比，其更能體現(xiàn)多路處理器的計算機(jī)優(yōu)勢，既能節(jié)約計算時間，又能很好地模擬材料的動態(tài)、連續(xù)變化。

4 結(jié)語

彈丸侵徹形成了對多層復(fù)合材料的高速沖擊，時間短且速度快。彈丸和多層復(fù)合材料之間相互作用的動作復(fù)雜、過程暫短，很難捕捉材料從變形、裂隙出現(xiàn)和裂紋擴(kuò)展的活動過程?；诜蔷植孔饔盟枷虢⒌慕鼒鰟恿W(xué)理論是國際上興起的、適用于模擬材料的損傷和斷裂過程的理論。該理論通過求解基于時間和位移的空間積分的運(yùn)動方程，重現(xiàn)高速沖擊下多層材料物質(zhì)粒的活動，避免了基于連續(xù)性假設(shè)建模和求解空間微分方程的傳統(tǒng)宏觀方法在面臨不連續(xù)問題時的奇異性。Silling 博士于2000 年首次提出近場動力學(xué)理論，突破了經(jīng)典分子動力學(xué)方法在計算尺度上的局限性，對宏觀、微觀不連續(xù)力學(xué)問題進(jìn)行分析，包括均勻與非均勻材料的結(jié)構(gòu)變形、損傷、斷裂、沖擊、穿透和失穩(wěn)問題。至今，近場動力學(xué)理論的研究范圍已經(jīng)擴(kuò)展至結(jié)晶相變動力學(xué)問題以及納米材料和結(jié)構(gòu)的破壞問題。研究的材料和結(jié)構(gòu)已涉及金屬、混凝土、多種復(fù)合材料、層合板結(jié)構(gòu)、玻璃、顆粒材料、木材以及納米纖維結(jié)構(gòu)等。

在實際工程中，近場動力學(xué)理論為破壞機(jī)理的研究提供了強(qiáng)有力的解決理論基礎(chǔ)。作為一種非局部理論，近場動力學(xué)避免了傳統(tǒng)方法在解決不連續(xù)問題過程中存在的奇異性，成為研究脆性材料破壞行為的一種新興研究理論，在研究彈丸高速沖擊下，近場動力學(xué)在復(fù)合材料損傷方面具有以下3 個特點(diǎn)： 1） 與經(jīng)典連續(xù)介質(zhì)力學(xué)不同，近場動力學(xué)是采用積分方程而不是位移分量的導(dǎo)數(shù)來表示的。近場動力學(xué)允許材料內(nèi)部自發(fā)的裂紋萌生與擴(kuò)展路徑自由，不需要引入額外的裂紋擴(kuò)展準(zhǔn)則。在近場動力學(xué)理論中，內(nèi)力是由連續(xù)體中任意2 個物質(zhì)點(diǎn)之間的本構(gòu)關(guān)系來表示的，損傷是本構(gòu)模型的一部分。2） 近場動力學(xué)理論研究基于有限的區(qū)域內(nèi)質(zhì)子間的相對位移和相互作用，結(jié)合物質(zhì)離散粒子理論，以積分形式重構(gòu)高速沖擊下連續(xù)體的運(yùn)動方程，允許位移場不連續(xù)，在處理不連續(xù)問題的過程中有數(shù)值優(yōu)勢，可以更好地解決不同物質(zhì)黏合后在沖擊荷載下質(zhì)子跨界面運(yùn)動及跨界面區(qū)域中質(zhì)子的運(yùn)動軌跡的連續(xù)性問題。3） 近場動力學(xué)的PD 計算模型的無網(wǎng)格架構(gòu)非常適合進(jìn)行并行計算，根據(jù)使用的處理器數(shù)量，并行計算可顯著提高計算效率。