MATLAB編程在層合板強(qiáng)度分析中的應(yīng)用1)

董紀(jì)偉 陳培見 程紅梅 羅寧 張桂民

(中國礦業(yè)大學(xué)力學(xué)與土木工程學(xué)院，江蘇徐州221116)

強(qiáng)度問題是復(fù)合材料板力學(xué)性能研究和結(jié)構(gòu)設(shè)計(jì)的核心問題之一。從分析尺度來看，板的強(qiáng)度分析包括宏觀力學(xué)、細(xì)觀力學(xué)和微納米力學(xué)分析[1]。其中，層合板強(qiáng)度的宏觀力學(xué)分析是復(fù)合材料力學(xué)課程的教學(xué)重點(diǎn)和難點(diǎn)。由于該問題理論性強(qiáng)、邏輯關(guān)系復(fù)雜、數(shù)學(xué)運(yùn)算繁瑣、手工計(jì)算量大，學(xué)生在學(xué)習(xí)時(shí)往往會(huì)產(chǎn)生畏難情緒，接受知識(shí)困難。另一方面，MATLAB以其強(qiáng)大的數(shù)學(xué)計(jì)算、符號(hào)運(yùn)算、圖形處理等功能和簡單易學(xué)的編程語言，作為輔助手段可促進(jìn)學(xué)生對知識(shí)的理解和應(yīng)用，提升教學(xué)效果。目前，在力學(xué)基礎(chǔ)課和專業(yè)課中已有許多應(yīng)用MATLAB的成功案例[2-5]。本文將探討一種基于MATLAB編程的層合板強(qiáng)度分析方法，通過該方法不僅可使問題的求解大為簡化，還便于討論主要參數(shù)對層合板強(qiáng)度的影響，為層合板的強(qiáng)度設(shè)計(jì)提供依據(jù)。

1 問題的提出[6]

如圖1所示的三層正交鋪設(shè)對稱層合板，其正交鋪設(shè)比m=0.2，總厚度為t，承受軸向拉力Nx作用。各單層材料為玻璃/環(huán)氧，其性能為：E1=54 GPa，E2=18 GPa，ν21=0.25，G12=8.8 GPa，Xt=Xc=1050 MPa，Yt=28 MPa，Yc=140 MPa，S=42 MPa。試用MATLAB求層合板的極限載荷Nx/t，并討論正交鋪設(shè)比對極限載荷的影響。

圖1 三層正交鋪設(shè)對稱層合板示意圖

正交鋪設(shè)層合板單層剛度形式簡單，強(qiáng)度分析的運(yùn)算量相對較小，因此常作為《復(fù)合材料力學(xué)》教材中的典型例題[6]。然而，由于層合板破壞過程復(fù)雜，破壞形式多樣，手工計(jì)算過程仍較為繁雜，且計(jì)算效率低，每次只能預(yù)測一種參數(shù)下的層合板強(qiáng)度。事實(shí)上，對于正交鋪設(shè)層合板，正交鋪設(shè)比m是影響其力學(xué)性能的主要參數(shù)，它表示0°層總厚度與90°層總厚度的比值。那么，既要研究正交鋪設(shè)比對層合板強(qiáng)度及破壞規(guī)律的影響，從而對該類層合板進(jìn)行合理的厚度設(shè)計(jì)，又要保證計(jì)算效率，手工計(jì)算顯然是不可行的，而應(yīng)采用參數(shù)化編程的方法對該問題進(jìn)行分析。

2 層合板的強(qiáng)度分析及MATLAB編程

要將復(fù)雜的層合板強(qiáng)度分析通過計(jì)算機(jī)編程實(shí)現(xiàn)，首先應(yīng)對強(qiáng)度分析過程中涉及到的基本理論和強(qiáng)度分析流程有一個(gè)清晰的認(rèn)識(shí)，再將其轉(zhuǎn)化為計(jì)算機(jī)語言進(jìn)行計(jì)算分析。

2.1 層合板的單層應(yīng)力計(jì)算

層合板的破壞是其單層板逐層破壞的過程，而大多數(shù)關(guān)于單層板的強(qiáng)度理論都是基于應(yīng)力的理論，因此，單層板的應(yīng)力計(jì)算是層合板強(qiáng)度分析的基礎(chǔ)。基于經(jīng)典層合板理論，層合板單層應(yīng)力分析流程如圖2所示。

首先，由各單層的正軸模量計(jì)算出層合板的正軸柔度矩陣S，然后對其求逆得到正軸剛度矩陣Q，如下所示

對于例題中的正交鋪設(shè)層合板，0°層的偏軸剛度和正軸剛度相同，90°層的偏軸剛度只需將正軸剛度系數(shù)Q11與Q22互換，即

其次，計(jì)算層合板的剛度和柔度。對于正交鋪設(shè)對稱層合板，不存在拉彎或拉剪的耦合，因此單向受拉時(shí)只需計(jì)算其拉伸剛度A，對其求逆得到層合板柔度A′，計(jì)算公式為

其中t1=mt/(m+1)，t2=t/(m+1)。

圖2 層合板單層應(yīng)力計(jì)算流程圖

再次，計(jì)算層合板的中面應(yīng)變。在軸向拉力Nx下，該層合板的內(nèi)力-應(yīng)變關(guān)系為

最后，計(jì)算各單層的應(yīng)力。由單層的應(yīng)力-應(yīng)變關(guān)系得0°和90°層的應(yīng)力分別為

單層板的剛度及層合板應(yīng)力的計(jì)算分別用兩個(gè)函數(shù)文件表示，程序如下：

2.2 單層板的強(qiáng)度準(zhǔn)則

進(jìn)行強(qiáng)度分析時(shí)，需選擇合適的強(qiáng)度理論判斷各單層是否發(fā)生強(qiáng)度失效?？紤]到單層板拉、壓強(qiáng)度的不同，本文采用目前最流行的Tsai-Wu張量強(qiáng)度理論作為單層板的強(qiáng)度失效判據(jù)，其表達(dá)式為

式中

其中Xt和Xc分別表示纖維方向(縱向)的抗拉、抗壓強(qiáng)度；Yt和Yc分別表示垂直纖維方向(橫向)的抗拉、抗壓強(qiáng)度；S表示面內(nèi)剪切強(qiáng)度。

計(jì)算Tsai-Wu準(zhǔn)則方程左邊項(xiàng)的MATLAB程序如下：

2.3 強(qiáng)度分析流程

單層的破壞并不意味著整個(gè)層合板的破壞，但會(huì)使該單層的剛度降低，通常情況下采用剛度折減法對其進(jìn)行處理。由于單層板的破壞具有各向異性，在對某破壞層進(jìn)行剛度折減時(shí)，并不是將其所有的剛度系數(shù)都降為零，而應(yīng)根據(jù)破壞類型對其剛度系數(shù)進(jìn)行相應(yīng)的折減，具體折減方案見表1。

表1 單層板的剛度系數(shù)折減方案

表1 中，各單層的破壞類型可由最大應(yīng)力理論確定。根據(jù)式(7)和式(8)計(jì)算出發(fā)生破壞單層主方向的應(yīng)力，若某個(gè)應(yīng)力分量接近相應(yīng)的基本強(qiáng)度值時(shí)，單層板就沿該方向發(fā)生破壞。例如：若σ1＞0，σ2＞0，σ1＜＜Xt，σ2≈Yt，|τ12|＜＜S，則該層就沿2方向發(fā)生破壞。

對破壞的單層進(jìn)行剛度折減后，應(yīng)重新計(jì)算層合板的剛度和柔度，再檢驗(yàn)是否有新的單層發(fā)生連鎖破壞，如果有，則重復(fù)剛度折減的步驟，直到?jīng)]有連鎖破壞發(fā)生為止。接下來，由于層合板剛度的降低，使其載荷-變形關(guān)系變?yōu)榉蔷€性，應(yīng)采用載荷增量法進(jìn)行下一步的計(jì)算，即：假設(shè)載荷增量ΔNx，由式(6)～式(8)計(jì)算中面應(yīng)變增量Δε0和各層應(yīng)力增量Δσ，將此時(shí)各層的總應(yīng)力σ再代入強(qiáng)度準(zhǔn)則中確定下一次破壞時(shí)的載荷。以此類推，直到所有層沿加載方向均不能承載，則層合板最終破壞，此時(shí)的載荷即為層合板的極限載荷。

基于上述理論，例題中層合板強(qiáng)度分析的流程如圖3所示。

圖3 層合板強(qiáng)度分析流程圖

求解過程中，需根據(jù)Tsai-Wu準(zhǔn)則確定單層破壞時(shí)的載荷Nx/t或載荷增量ΔNx/t，這就需要求解一元二次方程，手算比較困難，而利用MATLAB符號(hào)運(yùn)算的功能比較容易實(shí)現(xiàn)，以求解首次破壞時(shí)的Nx/t為例，其程序如下

按照圖3所示的分析流程，編制MATLAB主程序求解例題，函數(shù)文件為“function[ex，Nxdt]=OrSy TFA(E1，E2，miu21，G12，Xt，Xc，Yt，Yc，S，m)”，該函數(shù)返回的參數(shù)為層合板的極限載荷和極限應(yīng)變，而計(jì)算過程中的鋪層失效信息、每次破壞時(shí)的載荷和應(yīng)變通過fprint語句輸出到結(jié)果文件“fail.dat”中。

3 結(jié)果及討論

編制完M文件后，在MATLAB的命令行窗口對該函數(shù)進(jìn)行調(diào)用，語句為[ex，Nxdt]=OrSy TFA(54000，18000，0.25，8800，1050，1050，28，140，42，0.2)，執(zhí)行完該語句可得到本例的計(jì)算結(jié)果為：拉伸極限載荷(Nx/t)u=198.33 MPa，極限應(yīng)變?chǔ)?x=0.019 06，與文獻(xiàn)[6]的結(jié)果基本吻合，存在微小差別的原因在于文獻(xiàn)[6]采用的是Hill-Tsai強(qiáng)度準(zhǔn)則，而本文考慮了單層板拉、壓強(qiáng)度的不同，采用了Tsai-Wu強(qiáng)度準(zhǔn)則。

圖4 為拉伸極限載荷隨正交鋪設(shè)比的變化規(guī)律。從圖中可以看出，極限載荷Nx/t隨m的增大而增大，這是因?yàn)閤方向?yàn)?°層的縱向，但對90°層卻是橫向，則0°層x方向的承載能力遠(yuǎn)大于90°層，m增大即為0°層的比例增加，因此Nx/t也隨之增加。

圖4 拉伸極限載荷隨正交鋪設(shè)比的變化規(guī)律

但不難發(fā)現(xiàn)，當(dāng)m從0.1增加到0.15時(shí)，Nx/t增幅突然提高，隨后逐漸減小，這是為什么呢？可以通過在“fail.dat”文件中分析層合板的破壞過程尋找答案。圖5顯示了m分別為0.05，0.1和0.15時(shí)層合板的破壞過程。從中可以看出：當(dāng)m≤0.1時(shí)，90°層的破壞隨之引起0°層y方向的連鎖破壞，但0°層在x方向還可以繼續(xù)承載，需要用增量法進(jìn)行第二次計(jì)算，層合板經(jīng)歷兩次加載過程即發(fā)生破壞；但當(dāng)m＞0.1時(shí)，90°層的破壞不會(huì)引起0°層的連鎖破壞，層合板的破壞經(jīng)歷三次加載過程才最終破壞。由此可見，當(dāng)單層破壞時(shí)，檢驗(yàn)有無其他層發(fā)生連鎖破壞這一步是十分必要的。

圖5 層合板破壞過程的記錄文件

4 結(jié)論

本文以復(fù)合材料力學(xué)課程中的難點(diǎn)--層合板的強(qiáng)度分析為例，提出一種將MATLAB編程引入復(fù)合材料力學(xué)課程中的方法。相比傳統(tǒng)教學(xué)中繁瑣的公式計(jì)算和推導(dǎo)相比，該方法的優(yōu)點(diǎn)在于：(1)化繁為簡，提高計(jì)算效率。層合板的強(qiáng)度分析往往需要多步計(jì)算才能最終得到極限載荷，計(jì)算工作量龐大。但事實(shí)上，每步計(jì)算都包含許多重復(fù)的部分，如：單層板的應(yīng)力計(jì)算、單層破壞載荷的計(jì)算等。利用MATLAB編程，既可將這些需重復(fù)計(jì)算的部分編制為子函數(shù)進(jìn)行調(diào)用，使程序更簡潔，還可解決許多復(fù)雜的數(shù)學(xué)計(jì)算問題，使得分析過程大為簡化。(2)便于對該問題進(jìn)行參數(shù)化討論。采用參數(shù)化編程，不只局限于對一種參數(shù)下的層合板進(jìn)行強(qiáng)度分析，還可研究層合板的強(qiáng)度及破壞規(guī)律隨某些參數(shù)的變化規(guī)律，為層合板的設(shè)計(jì)提供可靠的理論依據(jù)。