UMAT子程序在復(fù)合材料強(qiáng)度分析中的應(yīng)用

UMAT子程序在復(fù)合材料強(qiáng)度分析中的應(yīng)用本例使用UMAT用戶子程序進(jìn)行復(fù)合材料單層板的應(yīng)力分析和漸進(jìn)損傷壓縮強(qiáng)度分析，介紹UMAT用戶子程序編寫(xiě)方法及在Abaqus/CAE中的設(shè)置。本章使用最大應(yīng)變強(qiáng)度理論作為復(fù)合材料單層板的失效準(zhǔn)則，相應(yīng)的Fortran程序簡(jiǎn)單易讀，便于理解UAMT子程序的工作原理。知識(shí)要點(diǎn)： 強(qiáng)度分析 UMAT用戶子程序 最大應(yīng)變理論 剛度折減講師：孔祥宏版本：Abq 6.14難度：關(guān)鍵詞：強(qiáng)度分析，UMAT&.1 本章內(nèi)容簡(jiǎn)介本章通過(guò)兩個(gè)實(shí)例介紹UMAT用戶子程序在復(fù)合材料單層板的應(yīng)力分析和強(qiáng)度分析中的應(yīng)用。在第一個(gè)實(shí)例中，對(duì)一個(gè)簡(jiǎn)單的復(fù)合材料單層板進(jìn)行應(yīng)力分析，UMAT子程序主要計(jì)算應(yīng)力，不進(jìn)行強(qiáng)度分析，本例用于驗(yàn)證UMAT子程序的計(jì)算精度。在第二個(gè)實(shí)例中，對(duì)復(fù)合材料單層板進(jìn)行漸進(jìn)損傷強(qiáng)度分析，UMAT子程序用于應(yīng)力計(jì)算、強(qiáng)度分析和剛度折減。本章所用復(fù)合材料為T(mén)700/BA9916，材料屬性如表&-1所示。表&-1 T700/BA9916材料屬性參數(shù)值強(qiáng)度值E1/GPa114XT/MPa2688E2/GPa8.61XC/MPa1458E3/GPa8.61YT/MPa69.5120.3YC/MPa236130.3ZT/MPa55.5230.45ZC/MPa175G12/GPa4.16SXY/MPa136G13/GPa4.16SXZ/MPa136G23/GPa3.0SYZ/MPa95.6&.2 實(shí)例一：UMAT用戶子程序應(yīng)力分析在使用UMAT用戶子程序進(jìn)行高級(jí)應(yīng)用之前，應(yīng)該先了解UMAT子程序，熟悉UMAT子程序的工作原理，了解UMAT中的參數(shù)、變量的含義。為了便于讀者快速了解和使用UMAT，本例通過(guò)復(fù)合材料單層板的應(yīng)力分析來(lái)介紹一個(gè)簡(jiǎn)單的UMAT子程序。讀者可將本例中的單層板替換為層壓板，進(jìn)行對(duì)比分析。&.2.1問(wèn)題描述復(fù)合材料單層板幾何尺寸為15mm10mm0.15mm，纖維方向?yàn)?5，單層板的3D實(shí)體模型如圖&-1所示，X軸方向?yàn)?方向，左側(cè)面施加X(jué)軸向?qū)ΨQ邊界條件，下側(cè)面施加Y軸向?qū)ΨQ邊界條件，垂直于Z軸且Z=0的平面施加Z軸向?qū)ΨQ邊界條件，右側(cè)面施加100MPa的拉力。圖&-1 單層板邊界條件及加載情況本例中單位系統(tǒng)為mm、MPa。&.2.2 UMAT用戶子程序本例使用的UMAT用戶子程序UMAT-Stress.for的全部代碼如下，字母C及“!”之后為注釋內(nèi)容。1 SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,2 1 RPL,DDSDDT,DRPLDE,DRPLDT,3 2 STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME,4 3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT,5 4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,JSTEP,KINC)6 C 7 INCLUDE ABA_PARAM.INC8 C 9 CHARACTER*80 CMNAME10 DIMENSION STRESS(NTENS),STATEV(NSTATV),11 1 DDSDDE(NTENS,NTENS),DDSDDT(NTENS),DRPLDE(NTENS),12 2 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),13 3 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3),14 4 JSTEP(4)15 16 DIMENSION EG(6), XNU(3,3), STRAND(6), C(6,6), STRESS0(6)17 C*18 C EG.E1,E2,E3,G12,G13,G2319 C XNU.NU12，NU21,NU13，NU31，NU23，NU3220 C STRAND.STRAINT AT THE END OF THE INCREMENT21 C C.6X6 STIFFNESS MATRIX22 C STRESS0.STRESS AT THE BEGINNING OF THE INCREMENT23 C*24 C INITIALIZE XNU & C MATRIX25 XNU=026 C=027 C GET THE MATERIAL PROPERTIES-ENGINEERING CONSTANTS28 EG(1) = PROPS(1) !E1,YOUNGS MODULUS IN DIRECTION 129 EG(2) = PROPS(2) !E2,YOUNGS MODULUS IN DIRECTION 230 EG(3) = EG(2) !E3,YOUNGS MODULUS IN DIRECTION 331 XNU(1,2) = PROPS(3) !POISONS RATIO POI_1232 XNU(2,1) = XNU(1,2)*EG(2)/EG(1) !POISONS RATIO POI_2133 XNU(1,3) = XNU(1,2) !POISONS RATIO POI_1334 XNU(3,1) = XNU(1,3)*EG(3)/EG(1) !POISONS RATIO POI_3135 XNU(2,3) = PROPS(4) !POISONS RATIO POI_2336 XNU(3,2) = XNU(2,3)*EG(3)/EG(2) !POISONS RATIO POI_3237 EG(4) = PROPS(5) !G12,SHEAR MODULUS IN 12 PLANE38 EG(5) = EG(4) !G13,SHEAR MODULUS IN 13 PLANE39 EG(6) = PROPS(6) !G23,SHEAR MODULUS IN 23 PLANE40 C*41 C FILL THE 6X6 STIFFNESS MATRIX C(6,6)42 RNU = 1/(1-XNU(1,2)*XNU(2,1)-XNU(1,3)*XNU(3,1)-43 1 XNU(3,2)*XNU(2,3)-2*XNU(1,3)*XNU(2,1)*XNU(3,2)44 C STIFFNESS MATRIX C(6,6)45 C(1,1) = EG(1)*(1-XNU(2,3)*XNU(3,2)*RNU46 C(2,2) = EG(2)*(1-XNU(1,3)*XNU(3,1)*RNU47 C(3,3) = EG(3)*(1-XNU(1,2)*XNU(2,1)*RNU48 C(4,4) = EG(4)49 C(5,5) = EG(5)50 C(6,6) = EG(6)51 C(1,2) = EG(1)*(XNU(2,1)+XNU(3,1)*XNU(2,3)*RNU52 C(2,1) = C(1,2)53 C(1,3) = EG(1)*(XNU(3,1)+XNU(2,1)*XNU(3,2)*RNU54 C(3,1) = C(1,3)55 C(2,3) = EG(2)*(XNU(3,2)+XNU(1,2)*XNU(3,1)*RNU56 C(3,2) = C(2,3)57 C*58 C CALCULATE STRAIN59 DO I = 1, 6 60 STRAND(I) = STRAN(I)+DSTRAN(I)61 ENDDO62 C CALCULATE STRESS63 DO I = 1, 664 STRESS0(I) = STRESS(I)65 STRESS(I) = 066 DO J = 1, 667 STRESS(I) = STRESS(I)+C(I,J)*STRAND(J)68 ENDDO69 ENDDO70 C CALCULATE SSE71 DO I = 1, 672 SSE = SSE+0.5*(STRESS0(I)+STRESS(I)*DSTRAN(I)73 ENDDO74 C*75 C UPDATE DDSDDE76 DO I = 1, 677 DO J = 1, 678 DDSDDE(I,J) = C(I,J)79 ENDDO80 ENDDO81 RETURN82 END第1到14行及第81、82行為UMAT子程序固定格式，其中，第1到5行括號(hào)內(nèi)的變量為UMAT子程序中可以使用的變量，第10到14行定義各變量數(shù)組的維數(shù)和長(zhǎng)度。部分主要變量的含義如表&-2所示。表&-2 UMAT部分變量名及其含義STRESS增量步開(kāi)始時(shí)的應(yīng)力(S11, S22, .)，用增量步結(jié)束時(shí)的應(yīng)力計(jì)算結(jié)果對(duì)其更新STATEV(NSTATV)狀態(tài)變量（狀態(tài)變量個(gè)數(shù)），如果在材料中定義了狀態(tài)變量，則在UMAT中需要對(duì)其更新STRAN增量步開(kāi)始時(shí)的應(yīng)變(E11, E22, .)DSTRAN當(dāng)前增量步的應(yīng)變?cè)隽?E11, E22, .)NDI, NSHR, NTENS應(yīng)力、應(yīng)變的個(gè)數(shù)，NDI為正應(yīng)力或正應(yīng)變的個(gè)數(shù)，NSHR為剪應(yīng)力或剪應(yīng)變的個(gè)數(shù)，NTENS=NDI+NSHRPROPS, NPROPS材料參數(shù)、材料參數(shù)的個(gè)數(shù)DDSDDE雅克比矩陣，SSE, SPD, SCD特定的彈性應(yīng)變能、塑性耗散、蠕變耗散，只對(duì)能量輸出有影響，對(duì)其他計(jì)算結(jié)果無(wú)影響，在UMAT中需要對(duì)其更新CELENT單元特征長(zhǎng)度第15到83行為用戶自己編寫(xiě)的固定格式的Fortran程序，用于計(jì)算剛度矩陣、應(yīng)力、應(yīng)變能、雅克比矩陣。由于本例中沒(méi)有使用狀態(tài)變量，因此不需要更新STATEV，只需要更新STRESS、DDSDDE和SSE即可。第16行定義了5個(gè)數(shù)組，其中EG、STRAND、STRESS0為一維數(shù)組，XNU、C為二維數(shù)組。第18到22行為注釋部分，EG存放材料的3個(gè)彈性模量和3個(gè)剪切模量；STRAND存放當(dāng)前增量步結(jié)束時(shí)的應(yīng)變(E11, E22, .)；STRESS0存放增量步開(kāi)始時(shí)的應(yīng)力(S11, S22, .)；XNU為33的二維矩陣，存放泊松比12、21、13、31、23、32；C為66的剛度矩陣。第25、26行初始化二維數(shù)組XNU、C，使其每個(gè)元素都為0。第28到39行，讀取材料常數(shù)，計(jì)算泊松比。泊松比的計(jì)算公式如下。 (&-1)第42到56行，計(jì)算剛度矩陣。剛度矩陣的計(jì)算公式如下。 (&-2)對(duì)于本例所用材料，由于E2=E3，G12=G13，12=13，所以21=31、23=32，可以將式(&-2)化簡(jiǎn)后在UMAT中計(jì)算剛度矩陣。本例為保證UMAT-Stress.for的可讀性，剛度矩陣的計(jì)算沒(méi)有化簡(jiǎn)，直接按照式(&-2)編寫(xiě)。第59到61行，計(jì)算當(dāng)前增量步的應(yīng)變，計(jì)算公式如式(&-3)，式中上標(biāo)表示增量步序號(hào)。 (&-3)第63到69行，將當(dāng)前分析步開(kāi)始時(shí)的應(yīng)力STRESS的值賦給STRESS0，然后計(jì)算當(dāng)前增量步的應(yīng)力并賦給STRESS，當(dāng)前增量步的應(yīng)力計(jì)算如式(&-4)，式中上標(biāo)表示增量步序號(hào)，應(yīng)力i和應(yīng)變i為列向量。 (&-4)第71到73行，計(jì)算應(yīng)變能，如式(&-5)所示，式中上標(biāo)表示增量步序號(hào)；下標(biāo)表示應(yīng)力、應(yīng)變?cè)隽康姆至啃蛱?hào)，其中下標(biāo)為1、2、3表示正應(yīng)力、正應(yīng)變?cè)隽浚?、5、6表示剪應(yīng)力、剪應(yīng)變?cè)隽俊?(&-5)第76到80行為更新雅克比矩陣。由于沒(méi)有對(duì)剛度矩陣沒(méi)變化，應(yīng)力與應(yīng)變的關(guān)系如式(&-4)所示，所以應(yīng)力增量與應(yīng)變?cè)隽康年P(guān)系如式(&-6)所示，所以雅克比矩陣就等于剛度矩陣。 (&-6)第76到80行可簡(jiǎn)寫(xiě)為一行，如下所示：DDSDDE(1:6,1:6) = C(1:6,1:6)，或DDSDDE = CUMAT中的剪應(yīng)變?yōu)楣こ碳魬?yīng)變。&.2.3復(fù)合材料單層板應(yīng)力分析1、 創(chuàng)建部件及劃分網(wǎng)格l 創(chuàng)建部件在Part模塊，單擊工具區(qū)的(Create Part)，在Create Part對(duì)話框中，Name后面輸入Lam-C，Modeling Space選擇3D，Type選擇Deformable，在Base Feature區(qū)域選擇Solid、Extrusion，Approximate size使用默認(rèn)的200，單擊Continue.進(jìn)入繪圖模式。單擊工具區(qū)的(Create Lines: Rectangle (4 Lines)，在提示區(qū)輸入第1個(gè)點(diǎn)的坐標(biāo)(0,0)后按回車鍵，再輸入第2個(gè)點(diǎn)的坐標(biāo)(15,10)后按回車，再按Esc鍵或單擊鼠標(biāo)中鍵。單擊提示區(qū)的Done或鼠標(biāo)中鍵，在Edit Base Extrusion對(duì)話框Depth后面輸入0.15，單擊OK完成。l 劃分網(wǎng)格在環(huán)境欄Module后面選擇Mesh，進(jìn)入Mesh模塊。環(huán)境欄中Object選擇Part: Lam-C。單擊工具區(qū)的(Seed Part)，在Global Seeds對(duì)話框中Approximate global size后面輸入0.5，單擊OK。單擊工具區(qū)的(Mesh Part)，單擊提示區(qū)的Yes或鼠標(biāo)中鍵，完成網(wǎng)格劃分，如圖&-2所示，板的厚度方向只劃分為1層單元。圖&-2 劃分網(wǎng)格單擊工具區(qū)的(Assign Element Type)，在Element Type對(duì)話框中，選擇依次選擇Standard，Linear，3D Stress，在Hex標(biāo)簽頁(yè)中勾選Reduced integration，在Element Controls區(qū)域Hourglass control選擇Enhanced，即選擇C3D8R單元，單擊OK完成。單擊工具區(qū)的(Assign Stack Direction)，在視圖區(qū)選擇部件平行于X-Y平面的面，單擊鼠標(biāo)中鍵或提示區(qū)的Yes完成。2、 創(chuàng)建材料并給部件賦材料屬性l 創(chuàng)建材料在環(huán)境欄Module后面選擇Property，進(jìn)入Property模塊。單擊工具區(qū)的(Create Material)，在Edit Material對(duì)話框中，Name后面輸入U(xiǎn)Mat-T700；單擊GeneralUser Material，在User Material區(qū)域中Data區(qū)域的Mechanical Constants一欄依次輸入114000， 8610，0.3，0.45，4160，3000，單擊OK完成。單擊工具區(qū)的(Create Material)，在Edit Material對(duì)話框中，Name后面輸入Mat-T700；單擊MechanicalElasticityElastic，在Elastic區(qū)域中Type選擇Engineering Constants，在Data區(qū)域輸入從左到右依次輸入114000，8610，8610，0.3，0.3，0.45，4160，4160，3000，單擊OK完成。材料UMat-T700用于UMAT用戶子程序，Mat-T700用于做對(duì)比分析。輸入數(shù)據(jù)時(shí)，每輸入完一行后按回車（Enter）鍵，光標(biāo)會(huì)自動(dòng)移到下一行。也可以通過(guò)右鍵快捷菜單添加或刪除一行。本例UMAT子程序較簡(jiǎn)單，不需要使用狀態(tài)變量，因此在材料UMat-T700中沒(méi)有定義Depvar。l 給部件賦材料屬性單擊工具區(qū)的(Create Composite Layup)，在打開(kāi)的對(duì)話框中Name使用默認(rèn)名稱，Initial ply count后面輸入1，Element Type選擇Solid，單擊Continue.；在Edit Composite Layup對(duì)話框中，Layup Orientation區(qū)域的Definition選擇Coordinate system，單擊Definition下一行的(Create Datum CSYS)。在Create Datum CSYS對(duì)話框中使用默認(rèn)名稱Datum csys-1，類型選擇Rectangular，單擊Continue.，在提示區(qū)輸入原點(diǎn)坐標(biāo)(0,0,0)后按回車鍵，再輸入(1,0,0)后按回車鍵，最后輸入(0,1,0)后按回車鍵，單擊Create Datum CSYS對(duì)話框的Cancel。在Edit Composite Layup對(duì)話框中單擊(Select CSYS.)，在視圖區(qū)選擇剛創(chuàng)建的Datum csys-1，Stacking Direction選擇Element direction 3，Rotation axis選擇Axis 3。在Plies標(biāo)簽頁(yè)中，雙擊Region，在視圖區(qū)選擇部件后單擊鼠標(biāo)中鍵；右單擊Material，在快捷菜單中單擊Edit Material.，在Select Material對(duì)話框中選擇UMat-T700，單擊OK；右單擊Element Relative Thickness，在快捷菜單中單擊Edit Thickness.，在Thickness對(duì)話框中Specify Value后面輸入1后單擊OK；在Rotation Angle一欄輸入0；Integration Points使用默認(rèn)的1，單擊OK完成。在定義復(fù)合材料鋪層時(shí)，視圖區(qū)部件上會(huì)顯示鋪層方向，在Edit Composite Layup對(duì)話框中的Display標(biāo)簽頁(yè)中可以設(shè)置所需顯示的方向，在視圖區(qū)部件上白色箭頭及字母S表示Stacking Direction。3、 裝配在環(huán)境欄Module后面選擇Assembly，進(jìn)入Assembly模塊。單擊工具區(qū)的(Create Instance)，在Create Instance對(duì)話框中選擇Parts: Lam-C，單擊OK完成。4、 創(chuàng)建分析步、設(shè)置輸出變量l 創(chuàng)建分析步在環(huán)境欄Module后面選擇Step，進(jìn)入Step模塊。單擊工具區(qū)的(Create Step)，在Create Step對(duì)話框中，在Initial分析步之后插入Static, General分析步，單擊Continue.；在Edit Step對(duì)話框中使用默認(rèn)設(shè)置，單擊OK完成。l 設(shè)置輸出變量單擊工具區(qū)的(Field Output Manager)，在Field Output Requests Manager對(duì)話框中，選中F-Output-1，單擊Edit.；在Edit Field Output Request對(duì)話框中，設(shè)置如圖&-3所示，輸出整個(gè)模型最后一個(gè)增量步的S、E、U，單擊OK完成。圖&-3 場(chǎng)輸出變量設(shè)置單擊工具區(qū)的(History Output Manager)，在History Output Requests Manager對(duì)話框中，選中H-Output-1，單擊Edit.；在Edit History Output Request對(duì)話框中，設(shè)置輸出整個(gè)模型的內(nèi)能和應(yīng)變能，即ALLIE和ALLSE，單擊OK完成。5、 創(chuàng)建邊界條件及施加載荷l 創(chuàng)建邊界條件在環(huán)境欄Module后面選擇Load，進(jìn)入Load模塊。單擊工具區(qū)的（Create Boundary Condition），在Create Boundary Condition對(duì)話框中，Name后面輸入BC-X，Step選擇Initial，Category選擇Mechanical，Types for Selected Step選擇Symmetry.，單擊Continue.；在視圖區(qū)選擇裝配實(shí)例Lam-C-1左側(cè)端面，即垂直于X軸且X=0的側(cè)面，單擊鼠標(biāo)中鍵或提示區(qū)的Done，在Edit Boundary Condition對(duì)話框中選擇XSYMM，單擊OK完成。類似操作，選擇Lam-C-1的下側(cè)端面，即垂直于Y軸且Y=0的側(cè)面，創(chuàng)建邊界條件BC-Y，邊界類型為YSYMM；選擇Lam-C-1垂直于Z軸且Z=0的側(cè)面，創(chuàng)建邊界條件BC-Z，邊界類型為ZSYMM。邊界條件及加載情況見(jiàn)圖&-1。l 施加載荷單擊工具區(qū)的(Create Load)，在Create Load對(duì)話框中，Name使用默認(rèn)的Load-1，Step選擇Step-1，Category選擇Mechanical，Types for Selected Step選擇Pressure，單擊Continue.；在視圖區(qū)選擇Lam-C-1的右側(cè)端面，即垂直于X軸且X=15的側(cè)面，單擊鼠標(biāo)中鍵，在Edit Load對(duì)話框中Magnitude后面輸入-100，單擊OK完成。6、 創(chuàng)建分析作業(yè)并提交分析l 創(chuàng)建分析作業(yè)在環(huán)境欄Module后面選擇Job，進(jìn)入Job模塊。單擊工具區(qū)的（Job Manager），在Job Manager對(duì)話框中單擊Create.；在Create Job對(duì)話框中，Name后面輸入Job-Lam-Stress-Umat，Source選擇Model-1，單擊Continue.；在Edit Job對(duì)話框的General標(biāo)簽頁(yè)中，單擊User subroutine file后面的(Select.)，在相應(yīng)路徑下找到并選擇UMAT-Stress.for文件；單擊Edit Job對(duì)話框中的OK完成。在Edit Job對(duì)話框的Parallelization標(biāo)簽頁(yè)中可以設(shè)置多核并行計(jì)算。l 提交分析在Job Manager對(duì)話框中，選中Job-Lam-Stress-Umat分析作業(yè)，單擊Submit提交計(jì)算。當(dāng)Job-Lam-Stress-Umat的狀態(tài)（Status）由Running變?yōu)镃ompleted時(shí)，計(jì)算完成，單擊Results進(jìn)入可視化后處理模塊。l 保存模型單擊工具欄的File工具條中的（Save Model Database），在Save Model Database As對(duì)話框的File Name后面輸入Laminate-Umat，單擊OK完成。7、 修改材料在環(huán)境欄Module后面選擇Property，進(jìn)入Property模塊。單擊工具區(qū)的(Composite Layup Manager)，在打開(kāi)的對(duì)話框中選擇CompositeLayup-1，單擊Edit.；在Edit Composite Layup對(duì)話框的Plies標(biāo)簽頁(yè)中，右單擊Material，單擊快捷菜單中的Edit Material.；在Select Material對(duì)話框中選擇Mat-T700，單擊OK；在Edit Composite Layup對(duì)話框中單擊OK完成。8、 再次創(chuàng)建分析作業(yè)并提交分析l 創(chuàng)建分析作業(yè)在環(huán)境欄Module后面選擇Job，進(jìn)入Job模塊。單擊工具區(qū)的（Job Manager），在Job Manager對(duì)話框中單擊Create.；在Create Job對(duì)話框中，Name后面輸入Job-Lam-Stress，Source選擇Model-1，單擊Continue.；在Edit Job對(duì)話框使用默認(rèn)設(shè)置，單擊OK完成。在Edit Job對(duì)話框的Parallelization標(biāo)簽頁(yè)中可以設(shè)置多核并行計(jì)算。l 提交分析在Job Manager對(duì)話框中，選中Job-Lam-Stress分析作業(yè)，單擊Submit提交計(jì)算。當(dāng)Job-Lam-Stress的狀態(tài)（Status）由Running變?yōu)镃ompleted時(shí)，計(jì)算完成，單擊Results進(jìn)入可視化后處理模塊。l 保存模型單擊工具欄的File工具條中的(Save Model Database)保存模型。9、 可視化后處理l 顯示云圖在視圖區(qū)顯示Job-Lam-Stress-Umat.odb。長(zhǎng)按工具區(qū)的(Plot Contours on Deformed Shape)，顯示隱藏工具后單擊(Plot Contours on Undeformed Shape)，在Field Output工具條中設(shè)置輸出S11。單擊菜單欄ResultSection Points.，打開(kāi)Section Points對(duì)話框，Selection method選擇Plies，在Plies區(qū)域選擇PLY-1，單擊Apply，在視圖區(qū)顯示該鋪層的S11應(yīng)力云圖。相同操作，可以顯示各鋪層的各應(yīng)力分離的云圖。Job-Lam-Stress-Umat.odb和Job-Lam-Stress.odb的各應(yīng)力分量的云圖如圖&-4所示。(a) S11，Job-Lam-Stress-Umat(b) S11，Job-Lam-Stress(c) S22，Job-Lam-Stress-Umat(d) S22，Job-Lam-Stress(e) S12，Job-Lam-Stress-Umat(f) S12，Job-Lam-Stress圖&-4 各應(yīng)力分量對(duì)比讀者在閱讀Abaqus幫助文件Abaqus Example Problems Guide中1.4.6 Failure of blunt notched fiber metal laminates一例的exa_fml_ortho_damage_umat.f文件時(shí)注意SSE的計(jì)算公式。l 輸出應(yīng)變能在視圖區(qū)顯示Job-Lam-Stress-Umat.odb。單擊工具區(qū)的(XY Data Manager)，在打開(kāi)的對(duì)話框中單擊Create.；在Create XY Data對(duì)話框中選擇ODB history output，單擊Continue.；在History Output對(duì)話框中選擇ALLSE，單擊Save As.；在Save XY Data As對(duì)話框中Name使用默認(rèn)的XYData-1，單擊OK；在XY Data Manager對(duì)話框中選擇XYData-1，單擊Edit.；在Edit XY Data對(duì)話框中可以看到應(yīng)變能數(shù)據(jù)，分析結(jié)束時(shí)整個(gè)模型的應(yīng)變能為9.671mJ。同樣操作，觀察Job-Lam-Stress.odb中整個(gè)模型的應(yīng)變能數(shù)據(jù)，同樣為9.673mJ。通過(guò)對(duì)比，驗(yàn)證了UMAT-Stress.for中第72行SSE計(jì)算的準(zhǔn)確性。單擊菜單欄ResultHistory Output.，可以直接打開(kāi)History Output對(duì)話框。&.2.4 Inp文件解釋本例中Job-Lam-Stress-Umat.inp節(jié)選如下。*Heading* Job name: Job-Lam-Stress-Umat Model name: Model-1* Generated by: Abaqus/CAE 6.14-1*Preprint, echo=NO, model=NO, history=NO, contact=NO* PARTS* 部件Lam-C-1的節(jié)點(diǎn)、單元數(shù)據(jù)*Part, name=Lam-C*Node 1, 15., 10., 0.15 節(jié)點(diǎn)編號(hào)及坐標(biāo)*Element, type=C3D8R 1, 64, 65, 23, 22, 43, 44, 2, 1 單元編號(hào)及節(jié)點(diǎn)編號(hào)*Orientation, name=Ori-1 1., 0., 0., 0., 1., 0.3, 0.* Section: CompositeLayup-1-1* 定義復(fù)合材料鋪層*Solid Section, elset=CompositeLayup-1-1, composite, orientation=Ori-1, controls=EC-1, stack direction=3, layup=CompositeLayup-11, 1, UMat-T700, 0., Ply-1*End Part* ASSEMBLY*Assembly, name=Assembly* 使用部件Lam-C創(chuàng)建裝配實(shí)例Lam-C-1*Instance, name=Lam-C-1, part=Lam-C*End Instance*End Assembly* ELEMENT CONTROLS* 單元控制，使用沙漏控制*Section Controls, name=EC-1, hourglass=ENHANCED1., 1., 1.* * MATERIALS* 使用工程常數(shù)創(chuàng)建材料Mat-T700*Material, name=Mat-T700*Elastic, type=ENGINEERING CONSTANTS114000.,8610.,8610., 0.3, 0.3, 0.45,4160.,4160.3000.,* 使用用戶材料創(chuàng)建材料UMat-T700*Material, name=UMat-T700*User Material, constants=6114000.,8610., 0.3, 0.45,4160.,3000.* * BOUNDARY CONDITIONS* 創(chuàng)建對(duì)稱邊界條件BC-X，BC-Y和BC-Z略* Name: BC-X Type: Symmetry/Antisymmetry/Encastre*Boundary_PickedSet4, XSYMM* -* STEP: Step-1* 創(chuàng)建分析步Step-1*Step, name=Step-1, nlgeom=NO*Static1., 1., 1e-05, 1.* * LOADS* 在Step-1創(chuàng)建Pressure類型的載荷* Name: Load-1 Type: Pressure*Dsload_PickedSurf6, P, -100.* * OUTPUT REQUESTS* 設(shè)置輸出變量*Restart, write, frequency=0* * FIELD OUTPUT: F-Output-1* 場(chǎng)輸出變量U、E、S*Output, field, frequency=99999*Node OutputU, *Element Output, directions=YESE, S* * HISTORY OUTPUT: H-Output-1* 設(shè)置歷史輸出變量*Output, history*Energy OutputALLIE, ALLSE*End Step&.2.5 應(yīng)用UMAT子程序應(yīng)力分析小結(jié)本例所介紹的應(yīng)力分析是UMAT用戶子程序最簡(jiǎn)單的應(yīng)用，讀者在了解UMAT子程序的基礎(chǔ)上，可以結(jié)合相應(yīng)的強(qiáng)度理論、剛度折減方法對(duì)UMAT子程序中的剛度矩陣、雅克比矩陣進(jìn)行計(jì)算和更新，從而達(dá)到剛度折減和漸進(jìn)損傷強(qiáng)度分析的目的。&.3 實(shí)例二：UMAT用戶子程序漸進(jìn)損傷強(qiáng)度分析本例在實(shí)例一使用的UMAT子程序的基礎(chǔ)上，通過(guò)增加失效判定，對(duì)UMAT子程序中的剛度矩陣進(jìn)行折減，達(dá)到漸進(jìn)損傷強(qiáng)度分析的目的。為了使本例的UMAT子程序簡(jiǎn)單易讀，本例對(duì)復(fù)合材料單層板進(jìn)行纖維方向的拉伸強(qiáng)度分析，僅考慮纖維方向拉伸破壞。本例使用本章實(shí)例一的模型，稍做修改，用于本例的強(qiáng)度分析。&.3.1問(wèn)題描述復(fù)合材料單層板的3D實(shí)體模型的幾何尺寸及邊界條件同本章的實(shí)例一的模型。將右側(cè)面100MPa的拉力替換為X軸正方向0.4mm的位移載荷。本例中單位系統(tǒng)為mm、MPa。&.3.2 UMAT用戶子程序本例使用的UMAT用戶子程序UMAT-Strength.for的全部代碼如下，字母C及“!”之后為注釋內(nèi)容。1 SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,2 1 RPL,DDSDDT,DRPLDE,DRPLDT,3 2 STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME,4 3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT,5 4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,JSTEP,KINC)6 C 7 INCLUDE ABA_PARAM.INC8 C 9 CHARACTER*80 CMNAME10 DIMENSION STRESS(NTENS),STATEV(NSTATV),11 1 DDSDDE(NTENS,NTENS),DDSDDT(NTENS),DRPLDE(NTENS),12 2 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),13 3 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3),14 4 JSTEP(4)1516 DIMENSION STRESS0(6),STRAND(6),C(6,6),CD(6,6),DCDE(6,6)17 PARAMETER (ALPHA = 1000.0, LAMBDA = 0.0, DMAX = 0.9999, DRND = 3)18 C*19 C STRESS0.STRESS AT THE BEGINNING OF THE INCREMENT20 C STRAND.STRAIN AT THE END OF THE INCREMENT21 C C.6X6 STIFFNESS MATRIX22 C CD.6X6 DAMAGED STIFFNESS MATRIX23 C DCDE.D CD/D E 24 C STATEV(1).DAMAGE VARIABLE D25 C*26 C GET THE MATERIAL PROPERTIES27 E1 = PROPS(1) !E1,YOUNGS MODULUS IN DIRECTION 128 E2 = PROPS(2) !E2 = E3,YOUNGS MODULUS IN DIRECTION 2 & 329 XNU12 = PROPS(3) !POISONS RATIO POI_12,XNU13 = XNU1230 XNU21 = XNU12*E2/E1 !POISONS RATIO POI_21,XNU31 = XNU21 31 XNU23 = PROPS(4) !POISONS RATIO POI_23,XNU32 = XNU2332 G12 = PROPS(5) !G12 = G13,SHEAR MODULUS IN 12 & 13 PLANE33 G23 = PROPS(6) !G23,SHEAR MODULUS IN 23 PLANE34 STH = PROPS(7) !FAILURE STRESS IN 1 DIRECTION IN TENSION35 C*36 C STIFFNESS MATRIX C(6,6)37 RNU = 1/(1-2*XNU12*XNU21-XNU23*2-2*XNU12*XNU21*XNU23)38 C = 039 C(1,1) = E1*(1-XNU23*2)*RNU40 C(2,2) = E2*(1-XNU12*XNU21)*RNU41 C(3,3) = C(2,2)42 C(4,4) = G1243 C(5,5) = G12 44 C(6,6) = G2345 C(1,2) = E1*(XNU21+XNU21*XNU23)*RNU46 C(2,1) = C(1,2)47 C(1,3) = C(1,2)48 C(3,1) = C(1,2)49 C(2,3) = E2*(XNU23+XNU12*XNU21)*RNU50 C(3,2) = C(2,3)51 C CALCULATE THE STRAIN AT THE END OF THE INCREMENT52 DO I = 1, 653 STRAND(I) = STRAN(I) + DSTRAN(I)54 ENDDO55 C*56 C CALCULATE THE FAILURE COEFFICIENT 57 STRANF = STH/E158 IF (STRAND(1) 0) THEN59 F = STRAND(1)/STRANF60 ELSE61 F = 062 ENDIF63 C CALCULATE D,DAMAGE V