向必偉， 姜大志， 謝成龍， 胡召奇

（1.安徽大學(xué) 資源與環(huán)境工程學(xué)院， 安徽 合肥 230601； 2.加拿大西安大略大學(xué) 地球科學(xué)系， 加拿大 倫敦N6A5B7； 3.合肥工業(yè)大學(xué) 資源與環(huán)境工程學(xué)院， 安徽 合肥 230009； 4. 安徽省地質(zhì)調(diào)查院， 安徽 合肥230001）

?

構(gòu)造變形流變場(chǎng)分配多尺度數(shù)值模擬研究

向必偉1， 姜大志2， 謝成龍3， 胡召奇4

（1.安徽大學(xué) 資源與環(huán)境工程學(xué)院， 安徽 合肥 230601； 2.加拿大西安大略大學(xué) 地球科學(xué)系， 加拿大 倫敦N6A5B7； 3.合肥工業(yè)大學(xué) 資源與環(huán)境工程學(xué)院， 安徽 合肥 230009； 4. 安徽省地質(zhì)調(diào)查院， 安徽 合肥230001）

摘 要：流變場(chǎng)分配（flow field partitioning）現(xiàn)象在自然界高應(yīng)變巖石中十分常見(jiàn)。傳統(tǒng)的基于固體連續(xù)變形機(jī)制理論的變形分析中， 流變場(chǎng)分配問(wèn)題往往被忽略， 致使應(yīng)變帶中流變場(chǎng)如何分配一直都缺乏深入認(rèn)識(shí)。Eshelby闡述了嵌入均勻介質(zhì)中的橢球體內(nèi)流變場(chǎng)的數(shù)學(xué)方法， 為探討流變場(chǎng)的分配奠定了理論基礎(chǔ)。本文從固體連續(xù)變形機(jī)制入手， 重點(diǎn)介紹基于Eshelby理論的多尺度數(shù)值模擬思路和方法， 探討流變場(chǎng)分配問(wèn)題。模擬結(jié)果表明: 嵌入基質(zhì)中的橢球體內(nèi)分布的流變場(chǎng)主要取決于橢球體與基質(zhì)間的相對(duì)流變強(qiáng)度， 橢球體的相對(duì)流變強(qiáng)度越低， 其變形越接近于簡(jiǎn)單剪切， 且有限應(yīng)變積累越快。模擬還揭示， 不同流變強(qiáng)度的橢球體內(nèi)模擬拉伸線理和面理產(chǎn)狀的總體格局反映基質(zhì)流變場(chǎng)特征。由此得到以下結(jié)論: （1）在流變場(chǎng)分配現(xiàn)象顯著的區(qū)域， 局部小尺度上應(yīng)變、渦度測(cè)量分析結(jié)果無(wú)法直接揭示區(qū)域流變場(chǎng)運(yùn)動(dòng)學(xué)邊界條件， 對(duì)這些區(qū)域的構(gòu)造變形分析必須是多尺度的； （2）基于Eshelby理論的以實(shí)質(zhì)變形組構(gòu)為約束的多尺度數(shù)值模擬分析， 能更為合理地揭示高應(yīng)變巖石中流變場(chǎng)的分配。

關(guān)鍵詞：流變場(chǎng)分配； 連續(xù)變形機(jī)制； Eshelby理論； 多尺度； 數(shù)值模擬

項(xiàng)目資助： 國(guó)家自然科學(xué)基金項(xiàng)目（41472194， 41002069， 40902064）資助。

0　引 言

構(gòu)造地質(zhì)學(xué)研究的一個(gè)關(guān)鍵目的是通過(guò)實(shí)際觀察的變形構(gòu)造揭示宏觀區(qū)域變形的流變場(chǎng)及其對(duì)應(yīng)的運(yùn)動(dòng)學(xué)邊界條件。Lister和Williams （1979）引入了簡(jiǎn)單剪切和純剪切變形的概念來(lái)描述巖石變形的流變場(chǎng)， Ramberg （1975）深化了這一概念并依據(jù)流變場(chǎng)中物質(zhì)線與瞬時(shí)拉伸軸間的關(guān)系對(duì)變形形式進(jìn)行了分類。同一時(shí)期研究者先后提出， 純剪切和簡(jiǎn)單剪切是兩種端員模式流變場(chǎng)， 自然界流變場(chǎng)是由這兩個(gè)分量構(gòu)成的， 兩者的相對(duì)大小可用運(yùn)動(dòng)學(xué)渦度數(shù)（kinematic vorticity number）來(lái)衡量（Truesdell， 1953；Means et al.， 1980）。

Ramsay和Graham開(kāi)創(chuàng)了巖石韌性變形的平面應(yīng)變模式， 并將應(yīng)變帶假設(shè)為局限于剛性邊界內(nèi)的均勻連續(xù)體， 應(yīng)變帶內(nèi)流變場(chǎng)由邊界的相對(duì)運(yùn)動(dòng)產(chǎn)生（Ramsay and Graham， 1970）。在Ramsay變形分析思路的啟發(fā)下， 許多學(xué)者對(duì)自然界應(yīng)變帶的流變場(chǎng)模式開(kāi)展了深入研究， 流變場(chǎng)模式也從二維平面應(yīng)變逐漸推廣到三維單斜變形模式（Sanderson and Marchini， 1984； Fossen and Tikoff， 1993； Tikoff and Fossen， 1993）和一般三斜變形模式（Jiang and Williams，1998； Lin et al.， 1998）。Ramasy及其后續(xù)變形模型的理論基礎(chǔ)是連續(xù)統(tǒng)內(nèi)的均勻介質(zhì)變形（Jiang， 2013），其中任意一物質(zhì)點(diǎn)的運(yùn)動(dòng)學(xué)直接反應(yīng)了連續(xù)統(tǒng)的宏觀流變場(chǎng)特征?；谶@一假設(shè)， 實(shí)際構(gòu)造變形分析中研究者通常用從任意尺度構(gòu)造（如: 露頭構(gòu)造、手標(biāo)本構(gòu)造、顯微構(gòu)造）獲得的流變場(chǎng)代表宏觀區(qū)域流變場(chǎng)。

然而， 局部流變場(chǎng)與區(qū)域流變場(chǎng)往往并不一致（Kuiper et al.， 2011）， Lister and Williams （1983）稱之為流變場(chǎng)分配（flow field partitioning）。在野外觀察的基礎(chǔ)上人們對(duì)流變場(chǎng)分配取得了一些直觀的認(rèn)識(shí)，如自然界的應(yīng)變集中帶往往比圍巖“軟弱”， 易于變形； 簡(jiǎn)單剪切變形往往集中在“軟弱”的應(yīng)變局限化帶內(nèi)（Michibayashi and Murakami， 2007； Kuiper et al.，2011； Lee et al.， 2012； Bell et al.， 2013）。然而， 迄今為止對(duì)于流變場(chǎng)如何在不均勻的變形區(qū)域內(nèi)分配以及局部分配的流變場(chǎng)與宏觀區(qū)域流變場(chǎng)如何聯(lián)系仍然缺乏深入了解。正因如此， 雖然流變場(chǎng)分配的概念早已被構(gòu)造研究者認(rèn)同， 并被用來(lái)解釋一些區(qū)域內(nèi)局部的不均勻變形（Wilson and Powell， 2001；Misra et al.， 2009； Shigematsu et al.， 2009； Carreras et al.， 2013； Aydin and de Joussineau， 2014）， 然而在實(shí)際研究應(yīng)變帶的整體運(yùn)動(dòng)學(xué)邊界條件時(shí)卻極少考慮流變場(chǎng)分配因素。

Eshelby （1957， 1959）在連續(xù)變形假設(shè)基礎(chǔ)上，建立了均勻彈性基質(zhì)與嵌入其中的均勻橢球體之間流變場(chǎng)的數(shù)學(xué)聯(lián)系。Eshelby公式被拓展到探討嵌入在線性粘性牛頓體材料（Newtonian）中橢球體的變形（Bilby et al.， 1975； Ellis and Watkinson， 1987）。牛頓體的流變學(xué)強(qiáng)度一般用粘性系數(shù)表述， 而對(duì)于非線性粘性材料（power-law）， 其流變學(xué)強(qiáng)度受諸多因素影響， 如流變指數(shù)、流變速率等。為了研究非線性粘性材料基質(zhì)與嵌入的橢球體的變形， 非線性粘性材料的流變學(xué)強(qiáng)度用變形過(guò)程中的有效粘度（effect viscosity）描述， Lebensohn and Tomé （1993）、Jiang and Bentley （2012）以及Jiang （2013）對(duì)非線性粘性材料的有效粘度及其數(shù)值算法進(jìn)行了深入開(kāi)發(fā)。Eshelby公式奠定了研究宏觀區(qū)域流變場(chǎng)在局部變形體內(nèi)分配的理論基礎(chǔ)。眾多研究者利用給定的基質(zhì)流變場(chǎng)代表區(qū)域流變場(chǎng)， 基質(zhì)中嵌入的橢球體代表局部地質(zhì)體， 建立了多尺度流變數(shù)值模型（Jiang and Bentley， 2012； Jiang， 2013； Xiang and Jiang，2013）。這一模型在局部變形體與區(qū)域運(yùn)動(dòng)學(xué)邊界條件之間的鴻溝上建立了橋梁。

本文介紹基于Eshelby理論建立3D多尺度數(shù)值模型， 模擬研究流變場(chǎng)在不同尺度上的分配。我們將給定的基質(zhì)流變場(chǎng)對(duì)應(yīng)宏觀尺度上的運(yùn)動(dòng)學(xué)邊界條件， 基質(zhì)中可變形的橢球體視為局部地質(zhì)體， 橢球體內(nèi)部的有限應(yīng)變主軸模擬變形組構(gòu)。這一簡(jiǎn)單模型可以模擬三個(gè)尺度間的構(gòu)造變形關(guān)系，例如宏觀尺度、露頭尺度以及組構(gòu)尺度， 因此這一模型是多尺度的。模擬計(jì)算揭示， 橢球體相對(duì)于基質(zhì)的流變學(xué)強(qiáng)度越低， 其變形越傾向于簡(jiǎn)單剪切變形， 并且有限應(yīng)變的積累越快； 單一橢球體內(nèi)的流變場(chǎng)都不等同于宏觀流變場(chǎng)， 但全部橢球體內(nèi)發(fā)育的模擬面理和拉伸線理產(chǎn)狀的總體格局能反映基質(zhì)流變場(chǎng)的特點(diǎn)。模擬結(jié)果揭示天然巖石中流變場(chǎng)分配主要取決于巖石的相對(duì)流變強(qiáng)度。模擬結(jié)果還表明， 由于流變場(chǎng)的分配， 對(duì)自然界應(yīng)變帶的宏觀流變場(chǎng)和對(duì)應(yīng)的運(yùn)動(dòng)學(xué)邊界條件研究，利用傳統(tǒng)的以變形材料均勻?yàn)榍疤峒僭O(shè)的變形分析方法是不合適的。相對(duì)于單一尺度的連續(xù)統(tǒng)變形理論， 基于Eshelby理論的多尺度構(gòu)造變形分析思路在很大程度上避免了巖石的均一性假設(shè)， 更具普遍意義。

1　數(shù)值模擬理論背景

1.1坐標(biāo)系與參考系

數(shù)值模型的建立依賴于參考坐標(biāo)系， 在此我們首先定義下文所涉及的參考坐標(biāo)系和相關(guān)參數(shù)。在圖1a定義的右手固定坐標(biāo)系x1x2x3中， 均勻高應(yīng)變帶的走向平行于x1軸， 任意給定邊界速度v。邊界速度沿坐標(biāo)軸方向分解為v1， v2和v3（圖1a）， 相應(yīng)在應(yīng)變帶內(nèi)產(chǎn)生分別平行于x1， x2， x3的拉伸應(yīng)變率，，和平行于x1， x3的剪應(yīng)變率和（圖1b）?？臻g內(nèi)任意一個(gè)向量e的方向可用兩個(gè)角度θ和φ確定， θ為向量e在x1x2平面上的投影與x1軸正方向的夾角， φ為向量與x3軸正方向的夾角（圖1c）?？臻g內(nèi)任意一個(gè)橢球的方位由三個(gè)橢球半徑軸方位確定（三對(duì)θ和φ角）， 其中只有三個(gè)角度是獨(dú)立的。我們?nèi)∽畲髾E球半徑軸的θ和φ角和中間半徑軸的一個(gè)方位角φ （圖1d）， 當(dāng)φ = 90°時(shí)， φ為橢球中間軸在x1x2平面的投影與x1軸正方向的夾角， 但當(dāng)φ=90°時(shí)， φ為橢球中間半徑軸與x3軸正方向的夾角。

1.2均勻介質(zhì)的連續(xù)變形機(jī)制

圖1　數(shù)值模型坐標(biāo)系參數(shù)示圖Fig.1 Diagrammatic drawing for the coordinate system and parameters of the numerical model

均勻介質(zhì)的流變學(xué)在巖石流變學(xué)方面的應(yīng)用始于簡(jiǎn)單的二維變形研究。Ramsay和Graham遵循有限應(yīng)變思路闡述了剛性邊界條件下變形帶恒體積簡(jiǎn)單剪切型式的變形機(jī)制（平面應(yīng)變）（Ramsay and Graham， 1970； Ramsay， 1980）。Ramberg （1975）詳盡闡述了等體積、連續(xù)、穩(wěn)態(tài)條件下巖石有限應(yīng)變的二維一般剪切模式。Sanderson and Marchini （1984）提出了走滑擠壓變形（transpressional zones）的三維單斜分析模式。隨后大量三維分析模式陸續(xù)發(fā)表（Fossen and Tikoff， 1993； Krantz， 1995； Jones and Geoff， 1995； Teyssier et al.， 1995）。在此基礎(chǔ)上， Jiang and Williams （1998）提出了一般的三斜變形模式， 探討了三斜變形與以往研究的單斜變形的關(guān)系， 推演了三斜變形的速度梯度張量和位移梯度張量。

給定邊界速度（圖1a）均勻介質(zhì)變形的流變場(chǎng)用速度梯度張量L來(lái)表達(dá):

只要流變場(chǎng)的速度梯度張量L已知或者其增量函數(shù)已知， 就可以求解任一時(shí)刻流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度， 并可以求解經(jīng)歷任意時(shí)間段的變形后的有限應(yīng)變量大小、主應(yīng)變大小以及應(yīng)變主軸的方位。

1.3Eshelby理論

Eshelby （1957）提出了處理均勻彈性橢球嵌入體與均勻彈性基質(zhì)間應(yīng)變場(chǎng)關(guān)系的數(shù)學(xué)方法。這一方法隨后被拓展到應(yīng)用于線性牛頓體材料（Bilby et al.，1975）乃至非線性的粘性（power-law）材料（Lebensohn and Tomé， 1993； Jiang and Bentley， 2012； Jiang，2013）。牛頓體材料， 流變強(qiáng)度可由粘度表達(dá)。非線性粘性體， 其粘度由應(yīng)力指數(shù)、應(yīng)變率以及其他熱力學(xué)性質(zhì)決定， 在變形過(guò)程中呈非線性變化。設(shè)想在變形過(guò)程中的任一無(wú)限短時(shí)間內(nèi)， 非線性粘性體的實(shí)際粘度與應(yīng)變率間的非線性關(guān)系可以用一個(gè)線性方程近似擬合； 因此， 在這一時(shí)間段內(nèi)的變形也近似地滿足Eshelby理論（Lebensohn and Tomé， 1993；Jiang and Bentley， 2012）。此時(shí)， 非線性粘性體的實(shí)際粘度被稱為“有效粘度”（Lebensohn and Tomé，1993）， 在Eshelby模型中嵌入的橢球體與基質(zhì)間的有效粘度比值我們稱為有效粘度比（reff）。Jiang and Bentley （2012）系統(tǒng)地探討了非線性粘性體有效粘度的近似方法、計(jì)算精度以及數(shù)值方法。

Eshelby （1957）在數(shù)學(xué)上證明當(dāng)均勻基質(zhì)中的嵌入體為橢球體時(shí)， 變形過(guò)程中橢球體內(nèi)部的應(yīng)力和應(yīng)變都是處處一致的， 但隨著橢球的幾何形狀與方位的變化， 橢球內(nèi)部流變場(chǎng)是非穩(wěn)態(tài)的。對(duì)于各向同性非線性粘性材料， 橢球體內(nèi)流變場(chǎng)與遠(yuǎn)離嵌入體的基質(zhì)內(nèi)流變場(chǎng)滿足如下公式（Jiang， 2013）:

兩式中DE和DM分別是橢球嵌入體內(nèi)部和遠(yuǎn)離嵌入體的基質(zhì)內(nèi)的拉伸張量； WE和WM分別是橢球嵌入體內(nèi)部和遠(yuǎn)離嵌入體的基質(zhì)內(nèi)的渦度張量。JS是四階對(duì)稱單位張量； S和Π分別是對(duì)稱和反對(duì)稱的Eshelby張量。A是四階應(yīng)變率分布（Partitioning） 張量， 對(duì)各向同性材料其表達(dá)式為:

式（3）中n是非線性粘性基質(zhì)材料的應(yīng)力指數(shù)， reff是橢球嵌入體與基質(zhì)材料的有效粘度比。由于橢球嵌入體的有效粘度取決于其變形時(shí)的真實(shí)應(yīng)力（或者應(yīng)變率）， 嵌入體材料的有效粘度可通過(guò)Mancktelow和Jiang開(kāi)發(fā)的迭代程序進(jìn)行數(shù)值求解（Mancktelow，2013； Jiang， 2013）。

當(dāng)給定遠(yuǎn)離嵌入體的基質(zhì)中流變場(chǎng)速度梯度張量LM， 將LM分解為DM和WM， 我們就可以通過(guò)公式（2a）， （2b）及公式（3）求解分布于橢球嵌入體內(nèi)部流變場(chǎng)的DE和WE， 從而求解速度梯度張量LE。在以下模擬計(jì)算中我們假設(shè)基質(zhì)內(nèi)的流變場(chǎng)是均勻且穩(wěn)態(tài)的。此時(shí)， LM是恒定的， 并可由內(nèi)部運(yùn)動(dòng)學(xué)渦度刻畫(huà)。由于橢球嵌入體的形狀和方向是不斷變化的， 嵌入體內(nèi)部的流變場(chǎng)雖然是均勻的但是非穩(wěn)態(tài)的。由公式（2）確定的渦度張量WE， 包含有剛體旋轉(zhuǎn)分量。只有內(nèi)部渦度才能真正刻畫(huà)流變場(chǎng)LE的本質(zhì)， 利用Jiang （2014）開(kāi)發(fā)的應(yīng)用程序我們可以計(jì)算LE的內(nèi)部運(yùn)動(dòng)學(xué)渦度數(shù)。顯然， 橢球嵌入體內(nèi)部流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度是隨時(shí)間演化的。通過(guò)求解增量應(yīng)變， 我們可以求解橢球嵌入體經(jīng)歷任意時(shí)間段的變形后的有限應(yīng)變量大小、主應(yīng)變大小以及應(yīng)變主軸的方位。

2　數(shù)值模擬

天然構(gòu)造變形材料（巖石）組成成分往往十分多樣、結(jié)構(gòu)上具有多尺度特征。從變形巖石整體上（大尺度）看， 巖石由許多形狀各異、流變強(qiáng)度不同的次級(jí)尺度元素組成。顯然， 建立與天然巖石對(duì)等的理論模型來(lái)研究巖石變形行為， 是巖石流變學(xué)研究的理想目標(biāo)。然而， 受限于現(xiàn)有的材料力學(xué)理論以及對(duì)天然巖石力學(xué)性質(zhì)的了解， 這一目標(biāo)暫時(shí)還難以實(shí)現(xiàn)。因此， 對(duì)天然巖石構(gòu)成進(jìn)行適當(dāng)簡(jiǎn)化建立相對(duì)簡(jiǎn)單的相似模型就成為必然。目前， 材料力學(xué)領(lǐng)域在研究非均勻材料變形行為時(shí)， 往往將非均勻材料對(duì)等為一個(gè)力學(xué)性質(zhì)十分相近的均勻介質(zhì)， 這個(gè)介質(zhì)被稱為“均勻?qū)Φ冉橘|(zhì)（HEM）”（Mura， 1987； Lebensohn et al.， 1998）。借鑒這一思路， 我們?cè)谘兄茙r石整體（大尺度）變形時(shí)也將其視為一個(gè)均勻介質(zhì)， 在研究巖石中某一元素（次級(jí)尺度）時(shí)將其視為嵌入均勻介質(zhì)中的一個(gè)非均勻體。這樣， 我們就可以依據(jù)Eshelby理論建立簡(jiǎn)單的多尺度數(shù)值模型， 來(lái)模擬研究巖石的變形行為。

2.1數(shù)值模型

我們構(gòu)建一個(gè)由單一橢球體嵌入基質(zhì)材料的數(shù)值模型來(lái)模擬宏觀區(qū)域構(gòu)造與局部地質(zhì)體。數(shù)值模型中， 橢球體與基質(zhì)都是均勻的各向同性的非線性粘性（power-law）材料， 兩者的有效粘度不同。以同一應(yīng)變率狀態(tài)為參照， 我們用橢球體材料與基質(zhì)材料的流變應(yīng)力指數(shù)ne和nm， 以及橢球體與基質(zhì)材料的有效粘度比（reff）定義橢球體與基質(zhì)材料的相對(duì)流變學(xué)性質(zhì)（Jiang and Bentley， 2012）。橢球體的方向由三個(gè)球面角θ， φ和φ確定 （圖1d）， 橢球幾何形狀由橢球半徑比值（a1: a2: a3）確定（a1， a2， a3分別是最大半徑、中間半徑和最小半徑）。

在圖1a坐標(biāo)系內(nèi)給定恒定的基質(zhì)流變場(chǎng)LM。除橢球嵌入體外， 基質(zhì)的變形可視為均勻介質(zhì)的穩(wěn)態(tài)變形， LM就等同于等式1所定義的L。我們假設(shè)變形是恒體積的， 即。變形過(guò)程中跟蹤計(jì)算橢球嵌入體內(nèi)流變場(chǎng)速度梯度張量LE， 從而跟蹤流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)、有限應(yīng)變主軸的方位以及有限應(yīng)變量大小。顯然， 所有相關(guān)的參數(shù)都會(huì)影響橢球嵌入體內(nèi)的流變場(chǎng)， 包括:

通過(guò)系統(tǒng)地調(diào)整上述參數(shù)， 運(yùn)行大量不同橢球體初始形態(tài)和方位模擬計(jì)算， 考察基質(zhì)流變場(chǎng)在橢球體內(nèi)的分配特征以及影響流變場(chǎng)分配的因素。模擬過(guò)程中， 首先系統(tǒng)地考察橢球體在平面應(yīng)變型式下的基質(zhì)流變場(chǎng)中發(fā)生平面變形時(shí)流變場(chǎng)的分配特征（圖2）， 然后逐步考察一般三斜的基質(zhì)流變場(chǎng)在橢球體內(nèi)的分配。本次研究中， 我們主要用橢球體中流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)和有限應(yīng)變量大小來(lái)刻畫(huà)流變場(chǎng)的分配特征。

圖2　平面應(yīng)變橢球的空間方位Fig.2 Orientation of the plane straining ellipsoid

2.2橢球體內(nèi)的流變場(chǎng)特征

圖3和圖4分別揭示了平面應(yīng)變的橢球嵌入體內(nèi)流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)和有限應(yīng)變量的演化歷史。在相對(duì)粘度低的橢球嵌入體內(nèi)（reff＜ 1）， 流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)經(jīng)歷了初始階段的起伏變化后很快趨于穩(wěn)定（NM≈ 2）。reff越低， wk越易于達(dá)到穩(wěn)定， 且最終穩(wěn)定的wk值越接近于1（圖3）。本次模擬計(jì)算設(shè)定的橢球體與基質(zhì)流變強(qiáng)度差異處于自然界巖石強(qiáng)度差異范圍內(nèi)， 低流變學(xué)強(qiáng)度的橢球體的變形接近于簡(jiǎn)單剪切。例如， 當(dāng)基質(zhì)流變場(chǎng)的wk=0.6， 在有效粘度是基質(zhì)粘度的五分之一的橢球嵌入體內(nèi)， 流變場(chǎng)的wk約為0.95。橢球體與基質(zhì)積累有限應(yīng)變的速率也有明顯差異（圖4）。盡管基質(zhì)和橢球嵌入體內(nèi)有限應(yīng)變的積累都是非線性的， 但兩者的比值近似線性。有效粘度越低， 橢球嵌入體積累有限應(yīng)變的速度越快。當(dāng)有效粘度比并不十分顯著時(shí)， 低流變學(xué)強(qiáng)度的橢球嵌入體積累的有限應(yīng)變就可達(dá)到基質(zhì)材料的數(shù)倍， 如: 當(dāng)reff≈0.2時(shí)， NE/ NM≈3.5（圖4f）。

平面應(yīng)變流變場(chǎng)的分配規(guī)律也適合于一般三維流變場(chǎng)。如圖5所示， 不同型式的基質(zhì)流變場(chǎng)在一任意初始方位的三維橢球體內(nèi)的分配特點(diǎn)基本相同。橢球體內(nèi)流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)經(jīng)歷短暫的波動(dòng)后很快趨于穩(wěn)定， 且穩(wěn)定的渦度數(shù)完全取決于橢球體與基質(zhì)間的相對(duì)有效粘度。當(dāng)給定的基質(zhì)流變場(chǎng)的渦度數(shù)一致（圖5所示， wk=0.6）， 無(wú)論基質(zhì)流變場(chǎng)是平面應(yīng)變型式（圖5a）或一般單斜型式（圖5c）還是一般三斜型式（圖5e）， 有效粘度比近似相等的橢球體內(nèi)流變場(chǎng)的最終穩(wěn)定渦度數(shù)近似相等。在不同類型基質(zhì)流變場(chǎng)中， 基質(zhì)與橢球體內(nèi)的有限應(yīng)變量仍然呈現(xiàn)了近似的線性關(guān)系， 并且有效粘度越低的橢球體內(nèi)有限應(yīng)變積累速率越快（圖5b， d 和f）。

圖3　平面應(yīng)變橢球體內(nèi)流變場(chǎng)運(yùn)動(dòng)學(xué)渦度演化圖Fig.3 Evolution of the kinematical vorticity number of the flow field in plane straining ellipsoids

上述模擬計(jì)算表明， 流變場(chǎng)在不同流變學(xué)強(qiáng)度元素組成的介質(zhì)中呈現(xiàn)差異分布。介質(zhì)的相對(duì)流變學(xué)強(qiáng)度（粘度）越低， 其變形越接近于簡(jiǎn)單剪切變形，且有限應(yīng)變積累得越快。自然界中普遍存在的應(yīng)變非均勻分配（strain partition）的地質(zhì)現(xiàn)象（Lister and Williams， 1983； Ramsay， 1980； Hudleston， 1999；Hobbs et al.， 2011）與我們的理論結(jié)果是完全吻合的。

2.3組構(gòu)與區(qū)域流變場(chǎng)

圖4　平面應(yīng)變橢球有限應(yīng)變演化圖Fig.4 Evolution of the finite strain in plane straining ellipsoids

我們假設(shè)自然界高應(yīng)變帶由方向不同、形狀各異的地質(zhì)體構(gòu)成。在圖1a所示的參考坐標(biāo)系下， 將 300個(gè)在三維空間內(nèi)隨機(jī)定向的橢球體嵌入到同一均勻基質(zhì)中， 且假設(shè)任意兩個(gè)橢球體之間的距離足夠遠(yuǎn)， 變形過(guò)程中不會(huì)產(chǎn)生相互影響。橢球的長(zhǎng)半徑與中間半徑的比值在1到10之間隨機(jī)變化， 短半徑給定為1。這些橢球能夠覆蓋很大范圍的形態(tài)變化， 從近等軸的球狀（1∶1∶1）到層狀（10∶10∶1），再到桿狀（10∶1∶1）。我們同樣假設(shè)橢球嵌入體與基質(zhì)都是非線性粘性材料， 給定材料的應(yīng)力指數(shù)均為3。橢球嵌入體與基質(zhì)間的相對(duì)粘度比在0.1到0.5范圍內(nèi)隨機(jī)給定。以一般單斜形式的基質(zhì)流變場(chǎng)為例模擬橢球嵌入體內(nèi)拉伸線理與面理產(chǎn)狀的總體演化。此時(shí)， LM可寫(xiě)為:

等式（4）的LM刻畫(huà)了一個(gè)垂直的右行走滑應(yīng)變帶。其中，是平行于x2軸的縮短應(yīng)變率， μ從0到1之間取值，˙可由LM的運(yùn)動(dòng)學(xué)渦度數(shù)wk和μ確定。

圖5　三維橢球體內(nèi)流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)和有限應(yīng)變演化圖Fig.5 Evolution of the kinematical vorticity number and finite strain in a 3D ellipsoid embedded in different imposed flow field

基質(zhì)流變場(chǎng)特征決定了組構(gòu)產(chǎn)狀的整體格局（圖6）。橢球嵌入體內(nèi)有限應(yīng)變最小主軸S3的方位幾乎不受wk和μ的影響， 總是近平行于x2軸； 而最大主應(yīng)變軸S1的方位則由wk和μ控制。S1的方位可以從以x3軸為中心的點(diǎn)集密經(jīng)大圓環(huán)帶過(guò)渡到以x1軸為中心的點(diǎn)集密。其代表的拉伸線理產(chǎn)狀投影從垂直的點(diǎn)集密經(jīng)過(guò)大圓環(huán)帶過(guò)渡到水平點(diǎn)集密（圖6b）。

圖6　理論模擬變形組構(gòu)演變圖Fig.6 Evolution pattern of deformation fabrics from numerical modeling

3　討 論

3.1Eshelby理論在韌性變形研究中的優(yōu)勢(shì)

與傳統(tǒng)的韌性變形分析理論相比較， 利用Eshelby理論開(kāi)展多尺度韌性變形分析的最主要優(yōu)勢(shì)在于其假設(shè)了天然巖石構(gòu)成的非均勻性， 能夠合理地解釋自然界流變場(chǎng)的不均勻分配現(xiàn)象。自從Ramsay利用連續(xù)變形理論分析巖石構(gòu)造變形以來(lái)，這一理論在巖石韌性變形分析領(lǐng)域一直占有統(tǒng)治地位。巖石材料均勻是Ramsay理論的基本假設(shè)之一。在這一前提下， 自然界構(gòu)造變形區(qū)域流變場(chǎng)的不均勻分配現(xiàn)象始終難以得到合理的解釋。例如在均勻性假設(shè)條件下， 人們通常利用材料彈性失穩(wěn)模式來(lái)解釋?xiě)?yīng)變局限化帶的形成機(jī)制。而按照這一觀點(diǎn)，高應(yīng)變帶的尾端相容性問(wèn)題（Hudleston， 1999）始終難以解決； 盡管Ramsay （1980）提出了兩種解釋模式，但到目前為止未見(jiàn)實(shí)際資料的報(bào)道。Eshelby理論的初始邊界條件是均勻介質(zhì)中嵌入可變形的均勻橢球體， 基質(zhì)和橢球體具有不同的流變強(qiáng)度。這就強(qiáng)調(diào)了介質(zhì)與嵌入其中的橢球體間的不均勻特征。隨著均勻?qū)Φ冉橘|(zhì)（HEM）思路的提出以及可行的數(shù)值求解方法的開(kāi)發(fā)， 基于Eshelby理論的模型可以視為由流變強(qiáng)度不同的橢球體構(gòu)成（Jiang， 2014）； 研究其中一個(gè)橢球體內(nèi)的流變場(chǎng)時(shí)， 其他橢球體整體視為基質(zhì)， 基質(zhì)的流變性質(zhì)用實(shí)時(shí)求解的HEM的流變性質(zhì)近似（Jiang and Bentley， 2012）。顯然， 通過(guò)改進(jìn)的Eshelby理論模型進(jìn)一步強(qiáng)調(diào)了變形材料的不均勻性， 與天然巖石的構(gòu)成也就更加相似。上述我們利用簡(jiǎn)單的Eshelby理論模型的模擬計(jì)算， 揭示了天然巖石變形流變場(chǎng)的不均勻分配主要受巖石相對(duì)流變強(qiáng)度制約， 這一規(guī)律與已有的觀察結(jié)論一致。

與其他數(shù)值模擬方法相比， 基于Eshelby理論的數(shù)值模擬的優(yōu)勢(shì)在于大應(yīng)變和多對(duì)象同時(shí)研究。目前用于變形模擬研究的常用數(shù)值模擬方法主要有有限元法和有限差分法， 兩種方法都有成熟的商業(yè)軟件， 也都應(yīng)用于模擬巖石的韌性變形（Mancktelow，2011）。兩種方法的一個(gè)共同特點(diǎn)是將模型網(wǎng)格化，通過(guò)求解網(wǎng)格節(jié)點(diǎn)的運(yùn)動(dòng)學(xué)和力學(xué)參數(shù)來(lái)描述變形過(guò)程。通過(guò)對(duì)模型中不同區(qū)域內(nèi)的網(wǎng)格設(shè)定不同的力學(xué)參數(shù)， 也可以建立由不均勻材料構(gòu)成的數(shù)值模型。模擬過(guò)程中， 當(dāng)任意相鄰的兩個(gè)網(wǎng)格交差產(chǎn)生幾何破壞時(shí)， 模型材料變形的應(yīng)力、應(yīng)變連續(xù)性被破壞， 運(yùn)算就無(wú)法繼續(xù)進(jìn)行。已有的模擬研究經(jīng)驗(yàn)表明， 網(wǎng)格密度越大， 材料的力學(xué)性質(zhì)刻畫(huà)得就越精確， 計(jì)算精度也就越高， 但此時(shí)運(yùn)算速率也越慢且網(wǎng)格越容易發(fā)生幾何破壞， 模型材料經(jīng)歷的有限應(yīng)變量就越小。因此， 無(wú)法實(shí)現(xiàn)模型與天然巖石應(yīng)變相當(dāng)?shù)拇髴?yīng)變是基于有限元理論和有限差分理論的模擬方法的主要困難之一， 特別是刻畫(huà)具有較大流變學(xué)強(qiáng)度差異的非均勻材料模型尤為如此。理論上， 基于有限元理論和有限差分理論的數(shù)值模型可以無(wú)限包含非均勻元素， 但實(shí)際模擬實(shí)驗(yàn)過(guò)程中，這類模型能夠發(fā)生的有限應(yīng)變量有限而且運(yùn)算量也無(wú)法由一般硬件設(shè)備承擔(dān)。因此， 已有的基于這兩種數(shù)學(xué)理論在巖石變形方面的研究模型材料都較為單一， 往往僅包含一個(gè)或有限的幾個(gè)非均勻元素（Mancktelow， 2011）。已有研究表明（Jiang and Bentley，2012； Jiang， 2013； Xiang and Jiang， 2013）， 基于Eshelby理論的數(shù)值模型可以開(kāi)展大應(yīng)變模擬計(jì)算，并且可以創(chuàng)建由大量流變學(xué)性質(zhì)不同的元素構(gòu)成的數(shù)值模型。

3.2Eshelby數(shù)值模型在巖石變形中的局限性

同傳統(tǒng)的巖石韌性變形分析理論一樣， Eshelby理論也以連續(xù)變形假設(shè)為前提， 這一假設(shè)對(duì)于天然巖石變形可能過(guò)于嚴(yán)格。天然巖石在構(gòu)造變形過(guò)程中， 不同塊體（或元素）之間不僅發(fā)生連續(xù)變形， 還可能會(huì)發(fā)生相對(duì)滑移的非連續(xù)變形。顯然， 如果巖石變形過(guò)程中非連續(xù)變形十分顯著， 用Eshelby理論建立的連續(xù)變形模式來(lái)研究其變形行為就不合適了。然而， 在考慮引入新的理論建立包含非連續(xù)變形基質(zhì)的模型之前， 我們必須對(duì)巖石變形行為開(kāi)展系統(tǒng)的研究， 揭示巖石發(fā)生非連續(xù)變形的標(biāo)志， 非連續(xù)變形對(duì)巖石整體變形的貢獻(xiàn)大小等等。如果包含非連續(xù)變形的理論模型是研究天然巖石變形的最終理想模型， 那么基于Eshelby理論的數(shù)值模型只是人們邁向最終目標(biāo)的一個(gè)中間過(guò)程。

Eshelby理論對(duì)嵌入基質(zhì)中的次級(jí)元素的假設(shè)條件可能過(guò)于嚴(yán)格。一方面， 地質(zhì)體的初始幾何形態(tài)可能與標(biāo)準(zhǔn)橢球體有一定出入， 這樣在實(shí)際變形過(guò)程中其內(nèi)部的流變場(chǎng)并非嚴(yán)格的處處均一。另一方面， 利用均勻橢球體模擬地質(zhì)體也是簡(jiǎn)化地近似；天然地質(zhì)體的流變性質(zhì)往往并不完全均勻， 變形過(guò)程中， 流變學(xué)性質(zhì)的不均勻性也會(huì)導(dǎo)致實(shí)際流變場(chǎng)比理想的模擬流變場(chǎng)復(fù)雜得多。因此， 實(shí)際天然變形巖石中記錄流變場(chǎng)特征的組構(gòu)可能比數(shù)值模擬預(yù)測(cè)結(jié)果更為復(fù)雜， 數(shù)值模擬結(jié)果更多地揭示了巖石的總體變形行為（Xiang and Jiang， 2013）。此外， 天然地質(zhì)體還可能具有各向異性的流變學(xué)性質(zhì)。盡管對(duì)一些相對(duì)簡(jiǎn)單的各向異性體的變形已開(kāi)展了一些初步的模擬研究工作， 如云母殘斑晶的變形和定向（Chen et al.， 2014）， 但是多數(shù)情況下地質(zhì)體的各向異性特征并不被人們所掌握， 對(duì)構(gòu)建具有各向異性特征的數(shù)值模型缺乏基本的約束。因此， 地質(zhì)體的各向異性特征也會(huì)影響地質(zhì)體內(nèi)部的流變場(chǎng)以及變形組構(gòu)發(fā)育。

3.3韌性變形分析思路

本次模擬揭示， 構(gòu)造變形區(qū)域內(nèi)分配在不同流變學(xué)強(qiáng)度的局部巖石內(nèi)的流變場(chǎng)往往不同。因此， 傳統(tǒng)的以均勻性假設(shè)為前提的變形分析方法就難以勝任對(duì)流變場(chǎng)不均勻分配區(qū)域的宏觀流變場(chǎng)特征的研究了。首先， 以往的韌性構(gòu)造變形分析中， 拉伸線理方向或優(yōu)選方位直接代表宏觀流變場(chǎng)的剪切方向。本次模擬研究表明， 變形區(qū)域內(nèi)拉伸線理產(chǎn)狀的總體分布格局是由宏觀流變場(chǎng)形式?jīng)Q定的， 隨著流變場(chǎng)特征的變化， 拉伸線理優(yōu)選方位可以從與剪切方向一致漸變?yōu)榕c剪切方向垂直。模擬實(shí)驗(yàn)表明， 可以利用變形區(qū)域內(nèi)拉伸線理產(chǎn)狀的總體分布格局來(lái)揭示宏觀流變場(chǎng)形式。其次， 傳統(tǒng)研究方法中通過(guò)運(yùn)動(dòng)學(xué)渦度測(cè)量來(lái)直接指示區(qū)域上的變形程度和運(yùn)動(dòng)學(xué)邊界條件， 在流變場(chǎng)不均勻分配的區(qū)域就不合適了。本次模擬研究表明， 局部巖石內(nèi)流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)不僅由宏觀流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度決定， 而且還受局部巖石的相對(duì)流變學(xué)強(qiáng)度影響。對(duì)于如何利用從不同局部變形巖石中測(cè)量的運(yùn)動(dòng)學(xué)渦度數(shù)來(lái)揭示宏觀區(qū)域流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)， 從而揭示宏觀運(yùn)動(dòng)學(xué)邊界條件， 還需要進(jìn)一步研究， 我們將在后續(xù)研究中進(jìn)一步探討。

4　結(jié) 論

本次工作系統(tǒng)地綜述了基于Eshelby理論的多尺度數(shù)值模擬研究思路和方法。通過(guò)多尺度數(shù)值模擬， 我們對(duì)天然巖石的構(gòu)造變形研究取得了兩點(diǎn)初步認(rèn)識(shí): （1）構(gòu)造變形區(qū)域， 流變場(chǎng)的分配主要受巖石的相對(duì)流變學(xué)強(qiáng)度影響， 流變強(qiáng)度越低的巖石內(nèi)流變場(chǎng)的剪應(yīng)變?cè)郊星矣邢迲?yīng)變累積速度越快；（2）在流變場(chǎng)不均勻分配的區(qū)域開(kāi)展構(gòu)造變形分析，傳統(tǒng)的以變形巖石均勻性假設(shè)為前提的韌性變形分析方法難以勝任； 變形區(qū)域上宏觀流變場(chǎng)可以通過(guò)面理和拉伸線理組構(gòu)產(chǎn)狀的總體分布格局并結(jié)合可能的變形形式來(lái)推斷； 構(gòu)造變形區(qū)域上流變場(chǎng)的運(yùn)動(dòng)學(xué)渦度數(shù)無(wú)法從測(cè)量局部變形巖石中的運(yùn)動(dòng)學(xué)渦度數(shù)直接獲取。

致謝： 兩位審稿專家對(duì)本文的最終成稿提出了科學(xué)嚴(yán)謹(jǐn)?shù)男薷囊庖?jiàn)， 讓筆者獲益匪淺， 在此表示由衷的感謝。

參考文獻(xiàn)（References）:

Aydin A and de Joussineau G. 2014. The relationship between normal and strike-slip faults in Valley of Fire State Park， Nevada， and its implications for stress rotation and partitioning of deformation in the east-central Basin and Range. Journal of Structural Geology， 63: 12-26.

Bell T H， Rieuwers M T， Cihan M， Evans T P， Ham A P and Welch P W. 2013. Inter-relationships between deformation partitioning， metamorphism and tectonism. Tectonophysics， 587: 119-132.

Bilby B A， Eshelby J D and Kundu A K. 1975. The change of shape of a viscous ellipsoidal region embedded in a slowly deforming matrix having a different viscosity. Tectonophysics， 28（4）: 265-274.

Carreras J， Cosgrove J W and Druguet E. 2013. Strain partitioning in banded and/or anisotropic rocks: Implications for inferring tectonic regimes. Journal of Structural Geology， 50: 7-21.

Chen Y， Jiang D， Zhu G and Xiang B. 2014. The formation of micafish: A modeling investigation based on micromechanics. Journal of Structural Geology，68（Part B）: 300-315.

Elliott D. 1972. Deformation paths in structural geology. Geological Society of America Bulletin， 83（9）: 2621-2638.

Ellis M and Watkinson A J. 1987. Orogen-parallel extension and oblique tectonics: The relation between stretching lineations and relative plate motions. Geology， 15（11）: 1022-1026.

Eshelby J D. 1957. The determination of the elastic field of an ellipsoidal inclusion， and related problems. Proceedings of the Royal Society of London. Series A，Mathematical and Physical Sciences， 241（1226）: 376-396.

Eshelby J D. 1959. The elastic field outside an ellipsoidal inclusion. Proceedings of the Royal Society of London. Series A， Mathematical and Physical Sciences，252（1271）: 561-569.

Fossen H and Tikoff B. 1993. The deformation matrix for simultaneous simple shearing， pure shearing and volume change， and its application to transpressiontranstension tectonics. Journal of Structural Geology，15（3-5）: 413-422.

Goscombe B D and Gray D R. 2008. Structure and strain variation at mid-crustal levels in a transpressional orogen: A review of Kaoko Belt structure and the character of West Gondwana amalgamation and dispersal. Gondwana Research， 13（1）: 45-85.

Hobbs B E， Ord A and Regenauer-Lieb K. 2011. The thermodynamics of deformed metamorphic rocks: A review. Journal of Structural Geology， 33（5）: 758-818.

Hudleston P. 1999. Strain compatibility and shear zones: Is there a problem？ Journal of Structural Geology，21（8-9）: 923-932.

Jiang D. 2007. Sustainable transpression: An examination of strain and kinematics in deforming zones with migrating boundaries. Journal of Structural Geology，29（12）: 1984-2005.

Jiang D. 2013. The motion of deformable ellipsoids in power-law viscous materials: Formulation and numerical implementation of a micromechanical approach applicable to flow partitioning and heterogeneous deformation in Earth's lithosphere. Journal of Structural Geology， 50: 22-34.

Jiang D. 2014. Structural geology meets micromechanics: A self-consistent model for the multiscale deformationand fabric development in Earth's ductile lithosphere. Journal of Structural Geology， 68（Part B）: 247-272.

Jiang D and Bentley C. 2012. A micromechanical approach for simulating multi-scale fabrics in large-scale high-strain zones: Theory and application. Journal of Geophysical Research， 117（12）: 12201-12216.

Jiang D and Williams P F. 1998. High-strain zones: A unified model. Journal of Structural Geology， 20（8）: 1105-1120.

Jones R and Geoff T P W. 1995. Strain partitioning in transpression zones. Journal of Structural Geology，17（6）: 793-802.

Krantz R W. 1995. The transpressional strain model applied to strike-slip， oblique-convergent and obliquedivergent deformation. Journal of Structural Geology，17（8）: 1125-1137.

Kuiper Y D， Lin S and Jiang D. 2011. Deformation partitioning in transpressional shear zones with an along-strike stretch component: An example from the Superior Boundary Zone， Manitoba， Canada. Journal of Structural Geology， 33（3）: 192-202.

Lebensohn R A and Tomé C N. 1993. A self-consistent anisotropic approach for the simulation of plastic deformation and texture development of polycrystals: Application to zirconium alloys. Acta Metallurgica et Materialia， 41（9）: 2611-2624.

Lebensohn R A， Turner P A， Signorelli J W， Canova G R and Tomé C N. 1998. Calculation of intergranular stresses based on a large-strain viscoplastic self-consistent polycrystal model. Modelling and Simulation in Materials Science and Engineering， 6（4）: 447-465.

Lee P E， Jessup M J， Shaw C A， Hicks Iii G L and Allen J L. 2012. Strain partitioning in the mid-crust of a transpressional shear zone system: Insights from the Homestake and Slide Lake shear zones， central Colorado. Journal of Structural Geology， 39: 237-252.

Lin S， Jiang D and Williams P F. 1998. Transpression （or transtension） zones of triclinic symmetry: Natural example and theoretical modelling. Geological Society，London， Special Publications， 135（1）: 41-57.

Lister G S and Williams P F. 1979. Fabric development in shear zones: Theoretical controls and observed phenomena. Journal of Structural Geology， 1（4）: 283-297.

Lister G S and Williams P F. 1983. The partitioning of deformation in flowing rock masses. Tectonophysics， 92（1-3）: 1-33.

Mancktelow N S. 2011. Deformation of an elliptical inclusion in two-dimensional incompressible power-law viscous flow. Journal of Structural Geology， 33（9）: 1378-1393.

Mancktelow N S. 2013. Behaviour of an isolated rimmed elliptical inclusion in 2d slow incompressible viscous flow. Journal of Structural Geology， 46: 235-254.

Means W D， Hobbs B E， Lister G S and Williams P F. 1980. Vorticity and non-coaxiality in progressive deformations. Journal of Structural Geology， 2（3）: 371-378.

Michibayashi K and Murakami M. 2007. Development of a shear band cleavage as a result of strain partitioning. Journal of Structural Geology， 29（6）: 1070-1082.

Misra S， Burlini L and Burg J P. 2009. Strain localization and melt segregation in deforming metapelites. Physics of The Earth and Planetary Interiors， 177（3-4）: 173-179.

Mura T. 1987. Micromechanics of Defects in Solids. Martinus Nijhoff Publishers: 388.

Ramberg H. 1975. Particle paths， displacement and progressive strain applicable to rocks. Tectonophysics，28（1-2）: 1-37.

Ramsay J G. 1980. Shear zone geometry: A review. Journal of Structural Geology， 2（1-2）: 83-99.

Ramsay J G and Graham R H. 1970. Strain variation in shear belts. Canadian Journal of Earth Sciences， 7（3）: 786-813.

Sanderson D J and Marchini W R D. 1984. Transpression. Journal of Structural Geology， 6（5）: 449-458.

Shigematsu N， Fujimoto K， Ohtani T， Shibazaki B， Tomita T，Tanaka H and Miyashita Y. 2009. Localisation of plastic flow in the mid-crust along a crustal-scale fault: Insight from the Hatagawa Fault Zone， NE Japan. Journal of Structural Geology， 31（6）: 601-614.

Teyssier C， Tikoff B and Markley M. 1995. Oblique plate motion and continental tectonics. Geology， 23（5）: 447-450.

Tikoff B and Fossen H. 1993. Simultaneous pure and simple shear: The unifying deformation matrix. Tectonophysics，217（3-4）: 267-283.

Truesdell C A. 1953. Two measures of vorticity. Journal of Rational Mechanics and Analysis， 2: 173-217.

Truesdell C A and Toupin R A. 1960. The classic field theories // Flügge S. Encyclopedia of Physics， Vol. II: Principles of Classical Mechanics and Field Theory. Springer-Verlag， Berlin: 226-793.

Wilson C J L and Powell R. 2001. Strain localisation and high-grade metamorphism at Broken Hill， Australia: A view from the Southern Cross area. Tectonophysics，335（1-2）: 193-210.

Xiang B and Jiang D. 2013. Small-scale ductile shear zones as transposed rheologically weak domains: A numerical modeling investigation and practical application. Journal of Structural Geology， 54: 184-198.

Multi-scale Numerically Modeling Flow Field Partitioning During Structural Deformation

XIANG Biwei1， JIANG Dazhi2， XIE Chenglong3and HU Zhaoqi4
（1. School of Resources and Environmental Engineering， Anhui University， Hefei 230601， Anhui， China； 2. Department of Earth Sciences， University of Western Ontario， London N6A5B7， Canada； 3. School of Resources and Environmental Engineering， Hefei University of Technology， Hefei 230009， Anhui， China； 4. Geological Survey of Anhui Province， Hefei 230001， Anhui， China）

Abstract：Flow field partitioning is common in inhomogeneous rocks in high strain area. However， the classic structural analysis based on the continuum mechanics usually ignores the ways of flow field partition. Eshelby's formulation deals with the flow field in an ellipsoid embedded in a homogeneous matrix， which established the foundation for the study of flow field partitioning. In this paper， we begin with the continuum mechanism and focus on the multi-scales structural analysis concept and method by numerical modeling based on the Eshelby's formulation. Numerical modeling results show that the flow field in embedding ellipsoids depends on the relative rheology between the ellipsoid and matrix. Under an imposed matrix flow field， the weaker the ellipsoid is， the more noncoaxial the partitioned flow field is and the more dramatic finite strain increase. Modeling results also show that the imposed flow field can be characterized by the whole pattern composed of fabrics in all embeddings with different rheology. Thereafter， we got two conclusions: （1） in a flow field partitioning area， the traditional structural analysis by the finite and vorticity number measurements can not be directly used to explain the regional kinematical boundary conditions； （2） the traditional single-scale method is not enough to explain flow field partitioning and a multi-scale modeling limited by natural fabrics should be a feasible way to understand the kinematic boundary conditions of a regional deformation.

Keywords:flow field partitioning； continuum mechanism； Eshelby's theory； multi-scale； numerical modeling

中圖分類號(hào)：P542

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1001-1552（2016）01-0001-013

收稿日期：2014-09-27； 改回日期： 2015-05-08

第一作者簡(jiǎn)介：向必偉（1976-）， 男， 副教授， 從事造山帶構(gòu)造變形分析等方面研究。Email: xbw1977@163.com