唐偉琴，李大永，彭穎紅

(上海交通大學(xué) 機(jī)械與動(dòng)力工程學(xué)院，上海 200240)

塑性各向異性是造成金屬板材沖壓制耳現(xiàn)象的主要根源[1]。晶體塑性理論以晶體學(xué)織構(gòu)及其在成形過程中的不斷演化作為分析塑性應(yīng)變問題的基礎(chǔ)，比經(jīng)典塑性理論更忠實(shí)于各向異性性能的物理本質(zhì)，能夠更好地預(yù)測(cè)金屬板材沖壓制耳行為。隨著計(jì)算機(jī)技術(shù)的發(fā)展，晶體塑性理論已被廣泛應(yīng)用于模擬金屬板材的沖壓成形過程。董湘懷課題組[2-4]開發(fā)了相關(guān)晶體塑性模型，研究了不同初始織構(gòu)對(duì)面心立方結(jié)構(gòu)(FCC)板材深沖成形制耳的影響。NAKAMACHI等[5-6]利用彈-黏塑性有限單元方法成功模擬了體心立方結(jié)構(gòu)(BCC)板材深沖時(shí)的4個(gè)和6個(gè)制耳現(xiàn)象。KIM等[7]、ENGLER 等[8]、耿小亮等[9]和李大永等[10]分別將晶體塑性本構(gòu)模型與有限元相結(jié)合，研究FCC、BCC金屬板材沖壓成形過程中的織構(gòu)演化及其對(duì)制耳的影響。

然而，目前晶體塑性理論在沖壓成形制耳方面的應(yīng)用大都集中在立方結(jié)構(gòu)(FCC和BCC)金屬，其主要塑性變形機(jī)制為滑移，且滑移系和滑移方向確定。鎂合金為密排六方結(jié)構(gòu)(HCP)，主要滑移系有基面滑移、柱面滑移、錐面滑移和錐面滑移，中常溫下孿生也是協(xié)調(diào)c軸方向應(yīng)變的主要機(jī)制。孿生為單向切變機(jī)制，其激活不僅與加載方向有關(guān)，而且還與變形溫度、合金成分和晶格結(jié)構(gòu)等密切相關(guān)，使得滑移與孿生的相互作用十分復(fù)雜[11]。AL-SAMMAN等[12]和HUANG等[13]分別對(duì)不同織構(gòu)的鎂合金進(jìn)行了平面應(yīng)變壓縮實(shí)驗(yàn)和拉深實(shí)驗(yàn)，結(jié)果表明，不同的初始織構(gòu)將導(dǎo)致塑性變形過程中啟動(dòng)的變形機(jī)制不同。近年來，基于晶體塑性模型模擬鎂合金的變形行為及其在變形過程中的織構(gòu)特性得到越來越多的關(guān)注。對(duì)于鎂合金多晶體模擬，主要采用基于Taylor模型的模擬[14-16]和基于自洽模型(Self consistent, SC)的模擬[17-20]兩種方法。Taylor模型比較簡(jiǎn)單，數(shù)值上容易實(shí)現(xiàn)，但是該模型沒有考慮晶粒間的相互作用，不能描述變形過程中由于各向異性引起的變形不均勻現(xiàn)象[21]。自洽模型在處理孿生的問題上具有物理意義上的優(yōu)勢(shì)，至今已發(fā)展成為彈塑性自洽模型(Elastic-plastic self-consistent, EPSC)、黏塑性自洽模型(Viscoplastic self-consistent, VPSC)和彈黏塑性自洽模型(Elastic-viscoplastic self-consistent, EVPSC)，適用于小變形、大變形、蠕變行為等多種變形條件，但是當(dāng)本構(gòu)模型比較復(fù)雜時(shí)，需要反復(fù)迭代求解應(yīng)力-應(yīng)變關(guān)系，計(jì)算過程繁雜[22]。另一種方法是將晶體塑性模型與商用有限元軟件(如ABAQUS)結(jié)合，即采用晶體塑性有限元(Crystal plasticity finite element method,CPFEM)對(duì)鎂合金的塑性變形過程進(jìn)行模擬[23-25]，該方法為解決復(fù)雜邊界條件下的織構(gòu)特性及對(duì)宏觀力學(xué)行為的影響提供了有效途徑。鎂合金變形機(jī)制的多樣性和不確定性，以及沖壓過程本身的復(fù)雜性，使得從晶體塑性理論角度預(yù)測(cè)鎂合金的沖壓成形比立方結(jié)構(gòu)金屬復(fù)雜得多。目前，僅WALDE等[23]建立了黏塑性自洽多晶體本構(gòu)模型(VPSC)，并嵌入 ABAQUS有限元對(duì)鎂合金沖壓過程進(jìn)行了模擬，但是該模型假設(shè)初始織構(gòu)為簡(jiǎn)單的理想織構(gòu)，如(0°, 0°, 0°)或(0°, 90°,0°)，沒有對(duì)實(shí)際軋制鎂板的沖壓性能進(jìn)行研究。

本文作者首先將滑移主導(dǎo)的晶體塑性理論擴(kuò)展到耦合滑移、孿生的晶體塑性理論，建立適合密排六方金屬的晶體塑性本構(gòu)模型；在此基礎(chǔ)上，構(gòu)建鎂合金板材大變形晶體塑性有限元模型，預(yù)測(cè)AZ31軋制板材的杯形件沖壓制耳行為，并對(duì)該預(yù)測(cè)方法進(jìn)行試驗(yàn)驗(yàn)證。

1 密排六方單晶體塑性模型

1.1 滑移與孿生的耦合方法

建立 HCP晶體塑性本構(gòu)模型主要是要解決滑移與孿生的耦合方法，而滑移和孿生兩種變形方式耦合的關(guān)鍵問題是確定孿生對(duì)塑性應(yīng)變的貢獻(xiàn)以及孿生引起的晶體轉(zhuǎn)動(dòng)。

1.1.1 孿生對(duì)塑性應(yīng)變的貢獻(xiàn)

由于孿晶的極性，所以孿生只能沿著孿晶方向進(jìn)行。本模擬中假定孿生為特定的滑移系，將塑性應(yīng)變?cè)隽?Δ εp分解為孿生和滑移變形兩部分的加權(quán)平均，單個(gè)晶粒的塑性應(yīng)變?cè)隽勘磉_(dá)式如下[26]：

式中：∑ gβ為孿晶體積分?jǐn)?shù)，Ns和 Ntw分別為開動(dòng)滑移系和孿晶系的數(shù)目，Δ γα和分別為滑移系和孿晶系的剪切應(yīng)變?cè)隽?，采用TAKAHASHI等[27]提出的應(yīng)變速率無關(guān)的“連續(xù)迭代法”求解，P為Schmid張量，它根據(jù)滑移面法向n和滑移方向m定義如下：

1.1.2 孿生引起的晶體轉(zhuǎn)動(dòng)

孿生導(dǎo)致晶體轉(zhuǎn)動(dòng)的模式與滑移不同，本文作者采用主導(dǎo)孿晶旋轉(zhuǎn)(PTR)模型[28]，即根據(jù)占主導(dǎo)地位的孿晶系來計(jì)算孿生對(duì)于晶體轉(zhuǎn)動(dòng)的貢獻(xiàn)。計(jì)算步驟歸納如下。

在每一個(gè)變形步，第n個(gè)晶粒第β個(gè)孿晶系的孿晶體積分?jǐn)?shù)增量為

對(duì)所有晶粒、所有孿晶系進(jìn)行求和得到總的孿晶分?jǐn)?shù)(FR)：

式中：n

f為晶粒n的體積分?jǐn)?shù)。對(duì)所有已發(fā)生孿晶旋

轉(zhuǎn)的晶粒進(jìn)行求和得到有效孿晶分?jǐn)?shù)(FE)：

孿生引起旋轉(zhuǎn)的臨界值定義為

式中：A1、A2是材料常數(shù)。在每一個(gè)變形步里，若最活躍的孿晶系的累積孿晶分?jǐn)?shù)β,ng 超過臨界值TF，則整個(gè)晶粒就根據(jù)該孿晶系的取向進(jìn)行旋轉(zhuǎn)，同時(shí)將已發(fā)生旋轉(zhuǎn)的孿晶分?jǐn)?shù)加入有效孿晶分?jǐn)?shù)。

根據(jù)下式更新發(fā)生孿晶旋轉(zhuǎn)的取向：

式中：R和R′分別為發(fā)生孿晶旋轉(zhuǎn)前后的晶體旋轉(zhuǎn)矩陣；Rtw為孿晶與母體的旋轉(zhuǎn)關(guān)系；n為孿晶切變方向。

1.2 單晶體本構(gòu)方程

首先考慮單晶體的變形，現(xiàn)時(shí)構(gòu)型即t+Δt時(shí)刻的Cauchy應(yīng)力為

應(yīng)力增量σΔ由材料本構(gòu)關(guān)系確定：

式中： Ce為4階彈性張量；總應(yīng)變?cè)隽喀う趴捎珊暧^變形轉(zhuǎn)換得到，塑性應(yīng)變?cè)隽?Δ εp由式(1)計(jì)算得到。

1.3 開動(dòng)滑移/孿生系的確定

判斷滑移/孿生系是否開動(dòng)，即判定某一晶粒是發(fā)生滑移變形還是孿生變形，是計(jì)算晶體塑性變形的前提。根據(jù)拉伸實(shí)驗(yàn)，晶體滑移/孿生系能否開動(dòng)取決于拉應(yīng)力在滑移/孿生系上的切應(yīng)力分量的大小。第α滑移/孿生系上的剪切應(yīng)力為

根據(jù) Schmid臨界切應(yīng)力定律：切應(yīng)力分量達(dá)到臨界值是滑移/孿生系處于塑性屈服狀態(tài)的必要條件。但是臨界滑移/孿晶系并非一定是活動(dòng)滑移/孿晶系，塑性流動(dòng)的產(chǎn)生不僅要求滑移/孿生系處于臨界狀態(tài)，而且要求一致性條件得到滿足。因此，滑移/孿生系是否開動(dòng)可由以下條件進(jìn)行判斷：

1.4 硬化模型

式中：N為晶體滑移系的數(shù)目；αβh 為硬化模量，表

示滑移系β中的滑移剪切應(yīng)變對(duì)滑移系α所產(chǎn)生的硬化，可表示為

式中：αβ

q 為描述自硬化和潛在硬化行為的矩陣，可表示為

式中：q為通過潛在硬化實(shí)驗(yàn)確定的潛在硬化與自硬化的比值。δαβ為Kronecher符號(hào)。hα可以通過求解：

1.5 晶體取向的確定

晶體在塑性變形過程中的取向轉(zhuǎn)動(dòng)導(dǎo)致織構(gòu)演化，它由兩部分組成：晶格旋轉(zhuǎn)和孿生引起的晶體轉(zhuǎn)動(dòng)，后者根據(jù)式(8)和(9)進(jìn)行求解，而晶格旋轉(zhuǎn)引起的取向變化為

式中：R為發(fā)生晶格旋轉(zhuǎn)前的晶體旋轉(zhuǎn)矩陣，Δω*為晶格旋轉(zhuǎn)增量，可由下式計(jì)算：

式中：Δω為總體旋轉(zhuǎn)增量，可由宏觀變形轉(zhuǎn)換得到，Δωp為塑性旋轉(zhuǎn)增量，可由下式求解：

式中：W為Schmid張量，它根據(jù)滑移面法向n和滑移方向m定義如下：

2 鎂合金多晶體塑性模型

HCP多晶體模型主要有泰勒模型、自洽模型和有限單元模型(Finite element model, FEM)。文獻(xiàn)[29]通過對(duì)比模型單元數(shù)及每個(gè)單元代表晶粒數(shù)的不同組合分析得出如下結(jié)論：只要總的晶粒個(gè)數(shù)(單元個(gè)數(shù)和每個(gè)單元代表的晶粒個(gè)數(shù)的乘積)達(dá)到一定數(shù)量，就可以獲得比較穩(wěn)定的統(tǒng)計(jì)結(jié)果；該文獻(xiàn)也指出耦合有限元的彈塑性自洽模型比較耗時(shí)。泰勒模型雖然簡(jiǎn)單，但將孿生作為一個(gè)連續(xù)體，與實(shí)際情況不相符合，對(duì)孿生進(jìn)行模擬時(shí)需采取一些特殊處理手段，如考慮晶界滑移等晶間相互作用，以提高模擬精度。因此，本模擬中采用有限元方法耦合晶體塑性模擬鎂合金板材大變形，其中每個(gè)單元代表一個(gè)晶粒，因?yàn)橛糜邢迒卧椒M多晶體能夠同時(shí)滿足相容條件和平衡條件，而無需泰勒或自洽多晶體假設(shè)，降低了對(duì)計(jì)算機(jī)的性能要求，模擬精度也得到了保證。

鎂作為 HCP金屬，塑性變形時(shí)啟動(dòng)的滑移/孿生系組合及其各自的臨界剪切應(yīng)力(CRSS)存在不確定性及復(fù)雜性，并且隨著合金成分及溫度的變化而不同[11,15]，因此，對(duì)于給定的鎂合金材料及變形溫度，變形機(jī)制及模型參數(shù)的確定是建立晶體塑性有限元模型的關(guān)鍵。

本文作者基于以上 HCP單晶體塑性模型編制了用戶材料子程序，嵌入到 ABAQUS有限元中，根據(jù)各變形機(jī)制的特性及鎂板變形過程中的織構(gòu)演化特征確定微觀變形機(jī)制，通過擬合拉伸或壓縮變形的應(yīng)力-應(yīng)變曲線及各向異性參數(shù)r值獲得模型參數(shù)，從而建立了AZ31板材大變形晶體塑性有限元模型。

AZ31軋制板材的織構(gòu)如圖 1所示，初始織構(gòu)為大部分晶粒c軸平行于板材法向的基面織構(gòu)。在進(jìn)行有限元模擬時(shí)將測(cè)得的反映該織構(gòu)的一組歐拉角(Φ1,Φ,Φ2)作為鎂板多晶體初始取向賦給各個(gè)晶粒。

2.1 變形機(jī)制的確定

對(duì)于純鎂，可啟動(dòng)的變形機(jī)制有基面滑移，柱面滑移，錐面滑移，錐面滑移和拉伸孿生[28]，如圖2所示。

圖1 AZ31鎂合金軋制板材的初始織構(gòu)Fig. 1 Initial texture of hot rolled magnesium alloy sheet: {0001} and {1 010} pole figures

圖2 鎂合金的主要變形機(jī)制Fig. 2 Main slip/twinning systems in magnesium alloy: (a) Basal slip; (b) Prismatic slip; (c) Pyramidal slip;(d) Pyramidal slip; (e) Tensile twinning

由拉伸孿生啟動(dòng)條件及AZ31板材的初始織構(gòu)可知：當(dāng)沿著鎂板平面內(nèi)任何方向進(jìn)行拉伸時(shí)，拉伸孿晶處于硬取向而幾乎不能啟動(dòng)，因此，晶粒c軸方向的變形只能由錐面滑移來協(xié)調(diào)，即錐面滑移在拉伸變形過程中必不可少；基面滑移的CRSS比其他滑移系的小很多，在變形過程中通常都會(huì)啟動(dòng)；根據(jù)文獻(xiàn)[19, 30]可知，鎂合金板材單向拉伸后轉(zhuǎn)變?yōu)閧0 001}〈10 1 0〉單一織構(gòu)，而柱面滑移是導(dǎo)致鎂板該織構(gòu)形成的因素。

鑒于此，本文作者將AZ31板材150 ℃變形時(shí)的變形機(jī)制取為基面滑移、柱面滑移、錐面滑移和拉伸孿生。此處沒有考慮錐面滑移，因?yàn)殄F面滑移引起的變形和織構(gòu)演化可由基面滑移和柱面滑移組合作用代替。

2.2 模型參數(shù)的擬合

鎂合金作為各向異性材料，硬化參數(shù)的擬合須同時(shí)滿足各個(gè)方向的應(yīng)力-應(yīng)變曲線和表征塑性各向異性特征的塑性應(yīng)變比r值。本模擬中分別沿著軋向(RD)、45°方向、橫向(TD) 3個(gè)方向進(jìn)行單軸拉伸或壓縮，根據(jù)拉伸方向8%的伸長(zhǎng)量計(jì)算塑性應(yīng)變比r值，計(jì)算公式如下：

式中：lε、wε和tε分別為長(zhǎng)度、寬度和厚度方向的塑性應(yīng)變。

由于AZ31板材沿著面內(nèi)單向拉伸時(shí)拉伸孿生不啟動(dòng)，所以一方面可排除孿生的作用而先確定滑移系的參數(shù)，從而使得問題簡(jiǎn)化；另一方面，孿生的參數(shù)無法根據(jù)面內(nèi)單向拉伸來確定，本文作者根據(jù)沿 RD方向壓縮來確定，由于實(shí)驗(yàn)條件的限制，該實(shí)驗(yàn)值由文獻(xiàn)[19]獲得。

以沿RD方向拉伸和壓縮為例，圖3所示為實(shí)驗(yàn)測(cè)試和模擬計(jì)算的應(yīng)力-應(yīng)變曲線。由圖3可知，AZ31軋制板材除了平面各向異性外，還顯示了拉壓不對(duì)稱性，這是由于變形孿晶將引起進(jìn)一步滑移/孿生的變形阻力，該效應(yīng)通過設(shè)置潛在硬化系數(shù)(表1中“Latent”)來體現(xiàn)。圖4顯示了鎂板沿RD、45°、TD這3個(gè)方向拉伸變形過程中的r值演化，其中曲線為模擬結(jié)果，3個(gè)數(shù)據(jù)點(diǎn)為應(yīng)變0.08時(shí)對(duì)應(yīng)的實(shí)驗(yàn)結(jié)果。由圖可知，沿3個(gè)方向加載時(shí)r值均隨著應(yīng)變的增加而增大；應(yīng)變?yōu)?.08時(shí)，模擬所得的r值與實(shí)驗(yàn)測(cè)試值趨勢(shì)一致，即r值從RD到TD方向逐漸增大。表1所列為擬合得到的多晶體模型的硬化系數(shù)。

圖3 鎂板沿RD方向拉伸和壓縮的應(yīng)力-應(yīng)變曲線Fig. 3 Simulated stress-strain curves for tension or compression along RD together with experimental measurements

圖4 鎂板拉伸變形過程中r值的變化Fig. 4 Simulated evolution of r value during tension together with experimental measurements

表1 AZ31板材多晶體模型的硬化系數(shù)Table 1 Hardening coefficients for polycrystal model of AZ31 Mg alloy sheet

3 AZ31鎂板的沖壓制耳行為

3.1 制耳行為的模擬與分析

圖5所示為杯形件沖壓成形示意圖，模型幾何尺寸如表2所列。板材采用兩層C3D8R單元，沖頭、壓邊圈、凹模都定義為解析剛體，沖壓件與沖頭、壓邊圈、凹模之間的摩擦因數(shù)均為0.1，所用壓邊力為5 kN。將測(cè)得的AZ31板材的一組歐拉角(Φ1,Φ,Φ2)通過VUMAT子程序賦給各個(gè)晶粒，作為有限元模型的初始輸入條件。

圖5 杯形件沖壓成形示意圖Fig. 5 Geometric sketch for deep cup drawing tools of AZ31 Mg alloy (Dp—Punch diameter; Dd—Die diameter; Rp—Round radius of punch; Rd—Round radius of die; Db—Sheet diameter;tb—Sheet thickness)

表2 沖壓模具尺寸Table 2 Dimensions for cylindrical cup drawing tools (mm)

根據(jù)沖壓過程及沖壓件性能可知，杯底的變形幾乎可忽略，其微觀組織和織構(gòu)與初始軋板相似，本模擬中僅考慮杯壁的織構(gòu)變化情況。對(duì) AZ31板材在150 ℃條件下杯形件沖壓進(jìn)行有限元仿真，沖壓行程為30 mm。圖6所示為AZ31板材沖壓杯形件模擬結(jié)果。從圖 6(a)可以看出，杯形件呈現(xiàn)出顯著的偏離軋向45°制耳，圖6(b)顯示出顯著的(0 001)〈10 1 0〉織構(gòu)組分和拉伸孿晶組分。

不同變形機(jī)制對(duì)變形的貢獻(xiàn)按下式來計(jì)算：

圖6 AZ31板材150 ℃沖壓杯形件模擬結(jié)果Fig. 6 Simulation results for deep-drawn cup of AZ31 alloy sheet formed at 150 ℃: (a) Earing profiles; (b) Pole figures

式中：fn表示第n個(gè)晶粒的體積分?jǐn)?shù)； Δγn,α表示第n個(gè)晶粒中第α個(gè)滑移/孿生系上的剪切應(yīng)變?cè)隽?；α′表示所要?jì)算的各個(gè)變形機(jī)制；α表示考慮的所有的變形機(jī)制。計(jì)算結(jié)果如圖7所示?？梢钥闯?，沖壓初始階段拉伸孿生是主要的變形機(jī)制，柱面滑移也啟動(dòng)，隨著

孿生的進(jìn)行，晶粒轉(zhuǎn)到易滑移方向，基面滑移和錐面滑移也變得活躍。

圖7 150 ℃沖壓過程中各種變形機(jī)制的相對(duì)活動(dòng)率Fig. 7 Relative activity of each deformation mode as function of strain during deep drawing at 150 ℃

3.2 制耳行為的試驗(yàn)驗(yàn)證

對(duì)AZ31板材在150 ℃條件下進(jìn)行杯形件沖壓實(shí)驗(yàn)，所得沖壓件如圖8所示，圖中箭頭所示方向?yàn)檐堉品较?。沖壓件呈現(xiàn)出顯著的45°方向制耳。

圖8 AZ31板材150 ℃沖壓杯形件實(shí)驗(yàn)結(jié)果Fig. 8 Experimental results for deep-drawn cup of AZ31 alloy sheet formed at 150 ℃: (a) Photo of bottom; (b) Whole photo

圖9 沖壓杯形件的織構(gòu)分布實(shí)驗(yàn)結(jié)果Fig. 9 Measured texture for deep-drawn cups of AZ31 alloy sheet: (a) {0001} pole figures; (b){1010}pole figures

圖9所示為杯形件織構(gòu)的測(cè)試結(jié)果，極圖中顯示了拉伸孿晶織構(gòu)組分的大量存在，如圖 9(a)中虛線部分所示。拉伸孿生主要是由沿圓周的壓應(yīng)變引起的，因?yàn)锳Z31板材的初始織構(gòu)為大部分晶粒c軸平行于板厚方向，當(dāng)沿著圓周進(jìn)行壓縮時(shí)，板厚增加，即提供了板厚方向的拉應(yīng)力，從而晶粒c軸受拉，滿足拉伸孿生啟動(dòng)的條件。此外，圖 9同樣顯示出了(0 001)〈10 1 0〉織構(gòu)組分，該織構(gòu)組分的形成可歸結(jié)為柱面滑移的啟動(dòng)。

綜上所述，模擬所得的沖壓杯形件制耳輪廓和變形織構(gòu)均與實(shí)驗(yàn)結(jié)果吻合得較好。

4 結(jié)論

1) 建立了耦合滑移和孿生的晶體塑性有限元模型，該模型同時(shí)滿足各個(gè)方向的力學(xué)性能曲線和表征塑性各向異性特征的塑性應(yīng)變比r值，充分體現(xiàn)了鎂合金板材與織構(gòu)相關(guān)的各向異性性能。

2) 基于所開發(fā)的晶體塑性有限元模型對(duì)AZ31板材沖壓成形過程進(jìn)行有限元模擬，模擬的杯形件制耳輪廓和變形織構(gòu)均與實(shí)驗(yàn)結(jié)果具有較好的一致性。

[1] 毛衛(wèi)民. 金屬材料的晶體學(xué)織構(gòu)與各向異性[M]. 北京: 科學(xué)出版社, 2002.MAO Wei-min. Crystallographic texture and anisotropy of metals[M]. Beijing: Science Press, 2002.

[2] 鄭 瑩, 董湘懷, 李志剛. 晶體方位對(duì)鈑料成形制耳的影響[J]. 中國機(jī)械工程, 2001, 12(9): 1070-1072.ZHENG Ying, DONG Xiang-huai, LI Zhi-gang. The effect of crystal orientation on earing in sheet metal forming[J]. China Mechanical Engineering, 2001, 12(9): 1070-1072.

[3] 章海明, 董湘懷. 一種適用于率相關(guān)晶體塑性模型的準(zhǔn)隱式積分算法[J]. 精密成形工程, 2011, 3(6): 1-5.ZHANG Hai-ming, DONG Xiang-huai. An semi-implicit integration scheme for rate dependent crystal plasticity[J].Journal of Netshape Forming Engineering, 2011, 3(6): 1-5.

[4] ZHANG H M, DONG X H, WANG Q, ZENG Z. An effective semi-implicit integration scheme for rate dependent crystal plasticity using explicit finite element codes[J]. Computational Materials Science, 2012, 54: 208-218.

[5] NAKAMACHI E, XIE C L, HARIMOTO M. Drawability assessment of BCC steel sheet by using elastic/crystalline viscoplastic finite element analyses[J]. International Journal of Mechanical Sciences, 2001, 43(3): 631-652.

[6] XIE C L, NAKAMACHI E. Investigations of the formability of BCC steel sheets by using crystalline plasticity finite element analysis[J]. Materials & Design, 2002, 23(1): 59-68.

[7] KIM E Y, CHOI S H, SHIN E J, YOON J W. Simulation of earing behaviors in bake hardening steel exhibiting a strong off-γ-fiber component[J]. International Journal of Solids and Structures, 2012, 49(25): 3573-3581.

[8] ENGLER O, HIRSCH J. Polycrystal-plasticity simulation of six and eight ears in deep-drawn aluminum cups[J]. Materials Science and Engineering A, 2007, 452/453(15): 640-651.

[9] 耿小亮, 張克實(shí), 郭運(yùn)強(qiáng), 秦 亮. 非均質(zhì)材料微成形過程的晶體塑性模擬[J]. 材料科學(xué)與工藝, 2010, 18(3): 297-301.GENG Xiao-liang, ZHANG Ke-shi, GUO Yun-qiang, QIN Liang.Crystal plasticity simulation study on the micro-forming of heterogeneous materials[J]. Materials Science & Technology,2010, 18(3): 297-301.

[10] 李大永, 張少睿, 彭穎紅, 劉守榮, 仇素萍. 板材沖壓成形的晶體塑性有限元模擬[J]. 機(jī)械工程學(xué)報(bào), 2008, 44(1): 190-194.LI Da-yong, ZHANG Shao-rui, PENG Ying-hong, LIU Shou-rong, QIU Su-ping. Finite element simulation of sheet metal stamping with polycrystalline plasticity[J]. Chinese Journal of Mechanical Engineering, 2008, 44(1): 190-194.

[11] PRANTIL V C, JENKINS J T, DAWSON P R. Modeling deformation induced textures in titanium using analytical solutions for constrained single crystal response[J]. Journal of the Mechanics and Physics of Solids, 1995, 43: 1283-1302.

[12] AL-SAMMAN T, GOTTSTEIN G. Room temperature formability of a magnesium AZ31 alloy: Examining the role of texture on the deformation mechanisms[J]. Materials Science and Engineering A, 2008, 488: 406-414.

[13] HUANG X, SUZUKI K, CHINO Y, MABUCHI M. Influence of initial texture on cold deep drawability of Mg-3Al-1Zn alloy sheets[J]. Materials Science and Engineering A, 2013, 565:359-372.

[14] STYCZYNSKI A, HARTIG CH, BOHLEN J, LETZIG D. Cold rolling textures in AZ31 wrought magnesium alloy[J]. Scripta Materialia, 2004, 50(7): 943-947.

[15] STYCZYNSKI A, HARTIG CH, BOHLEN J, LETZIG D. Cold rolling textures in AZ31 wrought magnesium alloy[J]. Scripta Materialia, 2004, 50(7): 943-947.

[16] LéVESQUE J, INAL K, NEALE K W, MISHRA R K.Numerical modeling of formability of extruded magnesium alloy tubes[J]. International Journal of Plasticity, 2010, 26(1): 65-83.

[17] CLAUSEN B, TOMé C N, BROWN D W, AGNEW S R.Reorientation and stress relaxation due to twinning: Modeling and experimental characterization of Mg[J]. Acta Materialia,2008, 56(11): 2456-2468.

[18] HUANG S, ZHANG S, LI D, PENG Y. Modeling texture evolution during rolling process of AZ31 magnesium alloy with elasto-plastic self consistent model[J]. Transactions of Nonferrous Metals Society of China, 2011, 21: 1348-1354.

[19] AGNEW S R, DUYGULU ?. Plastic anisotropy and the role of non-basal slip in magnesium alloy AZ31B[J]. International Journal of Plasticity, 2005, 21(6): 1161-1193.

[20] WANG H, WU P D, TOMé C N, HUANG Y. A finite strain elastic-viscoplastic self-consistent model for polycrystalline materials[J]. Journal of the Mechanics and Physics of Solids,2010, 58(4): 594-612.

[21] POLMEAR I J. Magnesium alloys and applications[J]. Materials Science and Technology, 1994, 10(1): 1-16.

[22] 王 娜, 雷麗萍, 方 剛, 曾 攀. 鎂合金變形的晶體塑性有限元分析[J]. 稀有金屬, 2008, 32(6): 766-773.WANG Na, LEI Li-ping, FANG Gang, ZENG Pan. Deformation analysis of magnesium alloy based on crystal plasticity theory[J].Chinese Journal of Rare Metals, 2008, 32(6): 766-773.

[23] WALDE T, RIEDEL H. Simulation of earing during deep drawing of magnesium alloy AZ31[J]. Acta Materialia, 2007,55(3): 867-874.

[24] CHOI S H, KIM D H, LEE H W, SHIN E J. Simulation of texture evolution and macroscopic properties in Mg alloys using the crystal plasticity finite element method[J]. Materials Science and Engineering A, 2010, 527(4/5): 1151-1159.

[25] HUANG S, ZHANG S, LI D, PENG Y. Simulation of texture evolution during plastic deformation of FCC, BCC and HCP structured crystals with crystal plasticity based finite element method[J]. Transactions of Nonferrous Metals Society of China,2011, 21: 1817-1825.

[26] WANG X M, XU B X, YUE Z F. Micromechanical modelling of the effect of plastic deformation on the mechanical behaviour in pseudoelastic shape memory alloys[J]. International Journal of Plasticity, 2008, 24(8): 1307-1332.

[27] TAKAHASHI H, MOTOHASHI H, TOKUDA M, ABE T.Elastic-plastic finite element polycrystal model[J]. International Journal of Plasticity, 1994, 10(1): 63-80.

[28] YI S B, DAVIES C H J, BROKMEIER H G, BOLMARO R E,KAINER K U, HOMEYER J. Deformation and texture evolution in AZ31 magnesium alloy during uniaxial loading[J]. Acta Materialia, 2006, 54(2): 549-562.

[29] 黃詩堯. AZ31鎂合金擠壓成形微觀組織演化的試驗(yàn)研究與數(shù)值模擬[D]. 上海: 上海交通大學(xué), 2010.HUANG Shi-yao. Numerial and experimental research on the microstructure evolution of magnesium alloy AZ31 during extrusion process[D]. Shanghai: Shanghai Jiao Tong University,2010.

[30] DEL VALLE J A, RUANO O A. Effect of annealing treatments on the anisotropy of a magnesium alloy sheet processed by severe rolling[J]. Materials Letters, 2009, 63(17): 1551-1554.