1 引言
当前广泛应用于航空航天领域以及民用结构的纤维增强树脂基复合材料,一般是由增强纤维和树脂基体构成的单向板按一定的叠层顺序铺设组成。而复合材料的宏观力学性能主要由其各组分的力学特性以及细观力学性能决定[1]。对于纤维增强复合材料，由于其结构特点存在三相组分材料：增强相、基体相和界面相。3种组分由于材料性质不同，也分别呈现了不同的力学特性。一般增强相为高模量、高强度的纤维材料，常见的有碳纤维和玻璃纤维等；基体相常见的有聚合物和金属等模量和强度较低而延展性好的材料；增强相与基体相的结合面被称为界面相，界面相的性能与基体和纤维的力学性能不同，可以在2种组分材料间传递应力并影响复合材料内部细观的应力分布。由于复合材料的细观损伤与破坏较为复杂，其损伤包括了各组分材料的细观损伤和相互作用引起的损伤累积。基于宏观力学方法分析单向板的力学特性时往往将复合材料单向板看作各向异性材料，无法综合考虑材料细观的结构特征和组分性能。而细观力学方法从微观尺度上将各组分材料看作均质各向同性材料，分析各组分结合后的应力应变关系和破坏过程，因此细观力学方法常被用来分析复合材料单项板的力学响应。在研究纤维增强复合材料的损伤过程和失效机理时，有必要从细观尺度，结合材料的细观结构特点和力学性能，综合考虑各组分的力学响应以及三者之间的相互影响进行分析。
Voigt[1]最早基于混合率模型通过细观力学方法来预报复合材料的宏观材料特性。但该解析模型计算精度较低，无法准确反映复合材料在各方向上的力学特性。此后有学者在该混合率模型的基础上进一步提出修正的混合率模型，但该修正模型仍无法准确预报材料的宏观性能。尽管该解析模型形式较为简单且计算结果与试验一定程度上吻合，但是对于复合材料的各向异性特性以及各组分材料的力学特性缺少分析。
为了更准确地分析复合材料细观组分性能和材料宏观力学性能的关系，研究人员提出了通过有限元方法，结合细观力学理论分析复合材料的力学行为和失效的方法。该方法通过在复合材料内选取具有能够反映复合材料各组分构成的代表性典型结构，建立能够反映复合材料细观结构特征的代表性体积单元(representative volume element，RVE)，并赋予各组分材料相应的本构关系研究复合材料在细观的应力应变响应和失效过程中的破坏机理。Brockenbrough 等[2]通过细观力学有限元的方法研究了金属基复合材料的宏观力学特性以及细观结构特征和各组分体积分数对分析结果的影响。Sun 和 Vaidya[3]首先提出了在 RVE 模型中应用周期性边界条件，并结合平均应力平均应变理论，预测了单向板复合材料的弹性性能。Xia 等[4]提出了用于周期性代表性体积单元有限元分析的统一形式的周期性边界条件，但该统一周期性边界条件存在过度约束的问题。 Xia 等[5]进一步通过施加基于位移的周期性边界条件进行了复合材料的细观有限元分析，并认为该边界条件满足连续周期性的条件。Zhang 等[6]建立了包含一根纤维的复合材料单胞模型(unit cell model)，研究了单向板复合材料的偏轴力学响应。该模型中针对基体材料的非线性特性应用了粘弹性本构关系，模拟了单向板在 45° 偏轴拉/剪混合应力状态下的力学响应，但是模型中未考虑基体/纤维界面的影响。Paley 和 Aboudi[7]提出了基于细观力学的通用单胞方法(generalized method of cells，GMC)，该方法以单胞结构中的各子胞应变为求解量，使其适用于复合材料的弹塑性力学分析。唐占文[8]利用通用单胞方法研究了界面的力学特性对材料宏观模量的影响，认为界面模量对单向板复合材料的纵向模量没有明显影响。同时也研究了界面层厚度对材料横向力学性能的影响，认为材料两方向上的宏观弹性模量随界面厚度增加有一定的提高。Zhang 等[9]研究了碳纤维/环氧树脂和玻璃纤维/环氧树脂复合材料的界面形态，认为界面层的厚度应为纤维半径的 0～0.1 倍。Asp 等[10]研究了树脂基体的屈服和失效过程，并认为该材料对静水压力敏感，其拉伸和压缩性能呈较明显不对称现象。Fiedler 等[11]研究了环氧树脂材料的静态拉伸和压缩性能，实验结果表明该材料在单向载荷下具有较明显非线性力学响应。这些研究多集中于分析复合材料各组分在载荷作用下的力学响应，且一般认为树脂基体材料损伤后发生脆性破坏，对树脂基体的塑性变形较少考虑。
本工作主要研究纤维增强复合材料在纵向和横向载荷下各组分的力学特性及其相互影响，基于 Abaqus 有限元软件建立该复合材料的代表性体积单元模型，通过引入合适的本构模型和失效准则实现该材料的宏观性能预测以及细观损伤失效机理研究。
2 复合材料细观有限元模型
2.1 复合材料细观力学模型
2.1.1 纤维模型
本研究所针对的材料为 AS4/8552纤维增强环氧树脂基复合材料。由于基体和界面相的强度远小于纤维强度，且在关于碳纤维增强复合材料的横向力学测试中尚没有关于纤维横向破坏的研究报道[12]，可认为纤维增强复合材料的横向破坏模式中无纤维破坏的发生。因此本文中的数值模拟不考虑纤维的横向拉伸和压缩破坏，采用横观各向同性线弹性本构关系表征 AS4 碳纤维的弹性力学响应。根据Naya等[13]的研究，AS4碳纤维的纵向和横向弹性模量、剪切模量以及泊松比如表1所示[13]。
表1 AS4 碳纤维材料性能参数
Table 1 Properties of AS4 carbon fibre
E1 /GPaE2 /GPaG12 /GPaG23 /GPaυ12υ2323112.9711.284.450.30.46
2.1.2 基体弹塑性模型
针对聚合物基复合材料的研究表明，静水压力对基体的力学性能具有显著的影响。静水压力会导致材料的弹性模量增加，且产生拉伸和压缩屈服应力不对称现象[14]。Vyas等[15]研究了树脂基复合材料在横向压缩过程中的静水压力依存现象，指出在压缩载荷作用下基体表现出更明显的塑性特性且具有更高的屈服强度。Vaughan 和 McCarthy[16]应用Mohr-Coulomb屈服准则来描述树脂基体的静水压力依存现象。由于该准则无法分析基体的损伤破坏过程，且Mohr-Coulomb 条件在主应力平面的屈服面的不等边六边形，其不连续点容易造成数值计算不收敛。因此，利用扩展 Drucker-Prager屈服准则结合金属延性损伤准则来描述该树脂基体分别在拉伸与压缩载荷下的非对称力学响应和渐进损伤过程。
扩展 Drucker-Prager 屈服准则实质上是 Mises 屈服准则的一种改进，该准则在主应力平面的屈服面为处处连续的椭圆，可表示为：
F=t-ptanβ-d=0
(1)
(2)
式中：p为静水压力；q为Mises等效应力；r为第三偏应力不变量； β为t-p应力面上屈服面的斜率，即内摩擦角；k为材料的三周拉伸和三轴压缩屈服应力比。当k取值为k=1时，该准则退化为Mises屈服条件，拉伸和压缩屈服应力相同，当k取值为0.778≤k≤1时，可使屈服面为外凸曲面。d为材料的内聚力，可由单向压缩屈服应力σc确定：
(3)
由于扩展Drucker-Prager准则是于Mohr-Coulomb准则的基础上提出，其材料参数β和d可由Mohr-Coulomb准则中转换得到：
(4)
(5)
Quinson等[17]研究了多种非晶态聚合物在单向拉伸、压缩和剪切载荷下的屈服强度和屈服点，并指该种材料的摩擦角取值为 7°～15°，对于 8552 树脂基体，其摩擦角为Φ=15°。因此，扩展 Drucker-Prager 准则的摩擦角取值为β=23.9°。材料的粘聚强度c决定了材料在纯剪切状态下的破坏强度，在缺少试验标定该参数时，Vogler和Kyriakides[18]认为粘聚强度可取材料的剪切强度，即 78.8 MPa。由式(5)可计算得到扩展 Drucker-Prager 准则的粘聚强度d=130.5MPa。
为反映基体材料的拉伸压缩不对称性，基体的损伤起始和损伤演化准则使用金属材料的延性损伤准则来描述。该准则将损伤起始时的等效塑性应变表示为应力三轴度η的函数，在纯剪切应力状态下η=0，在单向拉伸或压缩应力状态下应力三轴度可取值
基体损伤出现在当塑性应变达到损伤起始塑性应变
时，损伤变量可表示为：
(6)
式中：L为网格尺寸的特征长度；
为破坏时的等效塑性位移，且有
为破坏应力，Gf为单位面积断裂能。
(7)
其中，
和
分别为损伤起始和最终失效时的等效塑性应变。根据 Fiedler等[11]的研究，本文中拉伸和压缩损伤起始等效塑性应变分别取值为 0.05 和 0.5。参考 González 和 Llorca[12]的测试结果，取基体材料的断裂能为 0.5 J/m2。本文中数值模拟时使用的8552环氧树脂基体的力学性能参数如表2所示[12]。
表2 8552环氧树脂基体力学性能参数
Table 2 Properties of 8552 Epoxy matrix
Em/GPaυmσmyt/MPaσmyc/MPa)5.070.35121176β/(°)d/MPaG/(J·m-2)k23.9130.50.50.8
表2中：σmyt为拉伸强度； σmyc为压缩强度； β为内摩擦角。
2.1.3 纤维/基体界面内聚力模型
复合材料界面相有2种破坏模式：正向拉伸破坏和切向剪切破坏。Rodríguez等[19]采用纤维推出法测量了碳纤维增强环氧树脂基复合材料的界面剪切强度，测得的界面剪切强度为64 MPa。复合材料界面的正向强度无法通过试验获得，Rodríguez等[19]通过纤维推出试验结合有限元分析，认为界面的正向强度与切向强度相等。
采用内聚力模型(cohesive zone model)描述纤维与基体界面的损伤与失效过程。该模型由 Dugdale[20]在研究金属的塑性屈服时提出。常见的内聚力本构模型有双线型、三角函数型、抛物线型和指数型。采用双线型内聚力本构模型，由图1所示牵引-分离准则(traction-separation law)描述界面层单元的损伤起始、扩展和失效，图中横坐标为位移δ，纵坐标为应力σ，根据界面法向本构关系中的区别，内聚力单元的σ和δ这2个参数可定义为：
(8)
式中：Tn为法向牵引力；σmax为法向强度; Δn为法向张开位移； δn为达到损伤时法向张开位移。
图1 界面层内聚力模型的牵引-分离准则示意图
Fig.1 Traction-separation criterion of the cohesive model
由图示双线性内聚力模型可以看出，材料在达到强度极限前，界面单元呈完整状态，其材料性能保持良好，呈线弹性应力-应变关系，达到强度极限后刚度呈线性降低软化，线弹性阶段曲线斜率为刚度K，曲线所围成三角形面积为材料断裂释放的能量。
内聚力模型中材料的线弹性阶段本构关系可表示为：
σn=Kδn, σs=Kδs, σt=Kδt
(9)
式中：σn为界面的法向正应力；σs和σt为两个切向应力；δn、δs和δt分别为界面的法向位移和2个切向位移; K为界面层材料的初始刚度。界面初始刚度是内聚力模型中为满足2种单元之间保持刚性连接而引入的非物理量，K的取值需要足够大以消除纤维和界面间的应力不连续，本文中研究选取K=105 GPa/m[12]。
当界面处的牵引力达到极限强度时，单元开始产生损伤，界面不再呈线弹性力学响应。内聚力模型的损伤模型包括损伤起始判据和损伤演化判据。损伤起始判据采用二次名义应力准则描述，损伤演化过程为基于能量的Benzaggagh-Kenane (BK) 退化准则[21]。二次名义应力准则表达式如下所示：
(10)
其中，Nmax，Smax，Tmax为界面不产生损伤的最大法向正应力和2个切向应力，〈〉为Macaulay符号，有：
(11)
该式表示当界面层法向受压时，其压缩刚度保持不变，不随损伤演化而退化；当界面层法向受拉伸载荷并达到最大拉伸应力时，产生损伤，采用刚度折减法控制拉伸刚度的损伤退化：
(12)
式中：
和
分别为损伤后界面层的法向正应力和两个切向应力；d为界面损伤因子。
引入等效位移
描述损伤因子的演化：
(13)
则损伤因子d可表示为：
(14)
式中：
和
表示界面单元未损伤和完全失效时的等效位移；
为界面单元节点上的当前最大等效位移。
基于能量的BK混合损伤扩展判据为：
(15)
式中：
为Ⅰ型断裂韧性；
为Ⅱ型断裂韧性；GS为切应力做的功；GT为3个应力的总功。本文中模拟界面的断裂韧性的取值参考González和Llorca[12]的工作，认为其Ⅰ型和Ⅱ型断裂韧性均为 100 J/m2。AS4/8552纤维/基体界面特性参数如表3所示。
表3 AS4/8552 纤维/基体界面特性参数
Table 3 Properties of the interface of AS4/8552 composite
K/(GPa·m-1)σN/MPaτT/MPaτL/MPaG/(J·m-1)105646464100
2.2 代表性体积单元(RVE)
为准确分析材料的横向及纵向力学响应，本文中采用三维实体 RVE 模型来分析纤维增强复合材料的纵向与横向细观力学响应，纤维排布方式为菱形排列，如图2所示。
图2 纤维排布方式示意图
Fig.2 Illustration of the fibre alignment of the RVE
建立的复合材料RVE模型宽度w为32 μm，高度l为16 μm,厚度为1 μm,纤维直径7 μm，其体积分数为 50%，界面层厚度0.1 μm，为纤维半径的0.03倍。
采用8节点六面体线性减缩积分单元C3D8R对纤维和基体进行单元划分，对纤维/基体界面应用8节点三维cohesive单元COH3D8，单元长度为0.1 μm，这样可在界面处划分一层内聚力单元。所建立代表性体积单元模型如图3所示。
图3 RVE 模型示意图
Fig.3 Illustration of the RVE model
2.3 周期性边界条件
代表性体积单元是反映材料宏观特性的最小结构单元，其排布具有周期性的特点，对应的边界处会呈现周期性变化，RVE 单元在对应的边界处应满足连续性条件：应力连续和位移连续。因此，在利用有限元方法进行细观力学模拟时应施加周期性边界条件。根据 Li[22]提出的周期性边界条件，本文中所建立的三维RVE模型施加的周期性边界条件控制方程为：
(16)
式(13)中j+和j-分别表示单胞的第j组相对平行面，而已知在周期性单胞相对平行面的对应点上，
是一致的，由式(1)可以得到：
(17)
式(14)中，当单胞是具有周期性特点的RVE单胞时，
为一个常量，当给定一个确切的
值时，就可以确定式(14)中表达出的位移差的值，其为常值。
3 结果分析
由于纤维增强复合材料的纵向拉伸和压缩破坏模式较为复杂，纤维强度分布的随机性使得断裂过程很难用细观力学有限元的方法进行模拟。因此本文中主要分析复合材料的横向破坏过程。
通过给所建立纤维增强复合材料代表性体积单元模型施加不同位移载荷，采用 Abaqus/Explicit 显式有限元法计算了模型的横向力学响应；并使用均匀化方法得到材料的宏观模量性能，如表4所示[13]。由表中可以看出，模型分析结果与Naya等[13]的试验结果对比误差较小，证明本文中所建立的细观模型可以用于复合材料的横向渐进损伤和破坏分析。
表4 复合材料弹性模量模拟结果与试验结果
Table 4 Comparison of elastic modulus between simulation results and experimental results
项目E22/GPaE33/GPaG23/GPa模拟9.239.233.01试验9.69.63.19误差/%3.93.95.6
3.1 横向拉伸
图4为该RVE模型的横向拉伸应力-应变曲线。由图中可以看出，复合材料单向板在横向拉伸载荷下呈较明显的非线性应力-应变响应。这是由于树脂基体和界面的非线性特性所导致的。由均匀化方法得到材料的横向拉伸强度为62.9 MPa，对比试验结果[13]误差为1.5%。
图4 横向拉伸应力-应变曲线
Fig.4 Transvers tension stress-strain curve
图5、图6所示为该复合材料的横向拉伸损伤演化与失效过程。由图可见，在横向拉伸载荷作用下，基体出现较明显塑性变形，界面处的应力集中导致界面层损伤起始，随着载荷增加，界面处损伤进一步累积，引起界面层破坏，材料出现脱粘，基体应力随之增大，产生塑性变形并导致损伤累积，产生裂纹并逐渐扩展至界面处，直至与界面裂纹连通，最终导致材料拉伸失效。
图5 复合材料的横向拉伸损伤演化 与失效过程云图
Fig.5 Transverse tension damage process and failure
图6 复合材料的基体与界面层横向拉伸损伤演化 与失效过程云图
Fig.6 Transverse tension damage process and failure of matrix and interface
3.2 横向压缩
该单向板复合材料横向压缩应力-应变曲线如图7。由横向压缩损伤演化图与应力-应变曲线可以看出，相较于横向拉伸，复合材料的横向压缩具有更为明显的非线性特征，且具有更高的屈服强度，呈明显的拉/压不对称现象。由均匀化方法得到材料的横向压缩强度为 249.6 MPa，对比试验结果[13]误差为 6.9%。由图8 、图9所示模型的损伤和破坏过程可知，在压缩载荷作用下，基体/纤维界面处的剪应力引起界面初始损伤，同时基体的塑性变形增加，随着载荷进一步增加，基体多处出现剪切破坏形成裂纹，基体裂纹在剪应力作用下扩展直至连通形成塑性剪切带，并与界面裂纹贯通，材料发生最终破坏。
图7 横向压缩应力-应变曲线
Fig.7 Transverse compression stress-strain curve
图8 复合材料的横向压缩损伤演化 与失效过程云图
Fig.8 Transverse compression damage process and failure
图9 复合材料的基体与界面层横向压缩损伤 演化与失效过程云图
Fig.9 The transverse compression damage process and failure of matrix and interface
在考虑基体和界面损伤的条件下，AS4/8552复合材料的横向拉伸和压缩强度预测值与文献中试验结果[13]对比如表5所示[13]，其中，Yt为材料的横向拉伸破坏强度，Yc为横向压缩破坏强度。可以看出，材料的横向拉伸和压缩破坏强度预测值与文献中试验结果对比误差较小，说明本文中所建立的基于基体延性损伤模型和界面内聚力损伤模型的复合材料RVE模型可以用于树脂基复合材料的横向力学性能预测。
表5 AS4/8552复合材料横向拉伸和压缩强度 模拟结果与试验结果
Table 5 Comparison of transverse tension and compression strength between simulation results and experimental results
项目Yt/MPaYc/MPa模拟62.9249.6试验63.9268.0误差/%1.56.9
4 结论
1) 建立了纤维增强树脂基复合材料的细观代表性体积单元模型，结合扩展 Drucker-Prager 屈服准则和金属韧性准则描述树脂基体的损伤过程和失效；使用内聚力模型模拟纤维/基体界面的损伤和破坏；
2) 采用显式有限元方法分析复合材料的横向拉伸和压缩渐进损伤过程和破坏机理。计算结果与试验结果对比显示复合材料横向模量误差在 6% 以内，横向拉伸和压缩强度误差在 7% 以内，准确性较高；
3) 纤维增强复合材料横向拉伸破坏主要是基体开裂和界面脱粘；横向压缩破坏主要是基体压缩剪切破坏和界面剪切破坏。
[1] VOIGT W.Ueber die Beziehung zwischen den beiden Elasticitätsconstanten isotroper Körper[J/OL].Annalen der Physik,1889,274(12):573-587.DOI:10.1002/andp.18892741206.
[2] BROCKENBROUGH J R,SURESH S,WIENECKE H A.Deformation of metal-matrix composites with continuous fibers:geometrical effects of fiber distribution and shape[J/OL].Acta Metallurgica et Materialia,1991,39:735-752.DOI:10.1016/0956-7151(91)90274-5.
[3] SUN C T,VAIDYA R S.Prediction of composite properties from a representative volume element[J/OL].Composites Science and Technology,1996,56:171-179.DOI:10.1016/0266-3538(95)00141-7.
[4] XIA Z,ZHANG Y,ELLYIN F.A unified periodical boundary conditions for representative volume elements of composites and applications[J/OL].International Journal of Solids and Structures,2003,40(8):1907-1921.DOI:10.1016/S0020-7683(03)00024-6.
[5] XIA Z,ZHOU C,YONG Q,et al.On selection of repeated unit cell model and application of unified periodic boundary conditions in micro-mechanical analysis of composites[J/OL].International Journal of Solids and Structures,2006,43(02):266-278.DOI:10.1016/j.ijsolstr.2005.03.055.
[6] ZHANG Y,XIA Z,ELLYIN F.Nonlinear viscoelastic micromechanical analysis of fibre-reinforced polymer laminates with damage evolution[J/OL].Micromechanics of Materials,2005,42(02):591-604.DOI:10.1016/j.ijsolstr.2004.06.021.
[7] PALEY M,ABOUDI J.Micromechanical analysis of composites by the generalized cells model[J/OL].Mechanics of Materials,1992,14(02):127-139.DOI:10.1016/0167-6636(92)90010-B.
[8] TANG Zhanwen.Study on the analytical model of Macro-Meso-Progressive damage of composite materials considering interfacial phases[D].China:Harbin Institute of Technology,2013.
[9] ZHANG J,GUO Q,HUSON M,et al.Interphase study of thermoplastic modified epoxy matrix composites:Phase behaviour around a single fibre influenced by heating rate and surface treatment[J/OL].Composites Part A:Applied Science and Manufacturing,2010,41(6):787-794.DOI:10.1016/j.compositesa.2010.02.016.
[10] ASP L E,BERGLUND L A,TALREJA R.A criterion for crack initiation in glassy polymers subjected to a composite-like stress state[J/OL].Composites Science and Technology,1996,56(11):1291-1301.DOI:10.1016/S0266-3538(96)00090-5.
[11] FIEDLER B,HOJO M,OCHIAI S,et al.Failure behavior of an epoxy matrix under different kinds of static loading[J/OL].Composites Science and Technology,2001,61(11):1615-1624.DOI:10.1016/S0266-3538(01)00057-4.
[12] GONZ
LEZ C,LLORCA J.Mechanical behavior of unidirectional fiber-reinforced polymers under transverse compression:Microscopic mechanisms and modeling[J/OL].Composites Science and Technology,2007,67(13):2795-2806.DOI:10.1016/j.compscitech.2007.02.001.
[13] NAYA F,GONZLEZ C,LOPES C S,et al.Computational micromechanics of the transverse and shear behavior of unidirectional fiber reinforced polymers including environmental effects[J/OL].Composites Part A:Applied Science and Manufacturing,2017,92:146-157.DOI:10.1016/j.compositesa.2016.06.018.
[14] RHEE K Y,LEE S M,PARK S J.Effect of hydrostatic pressure on the mechanical behavior of seawater-absorbed carbon/epoxy composite[J/OL].Materials Science and Engineering:A,2004,384(01):308-313.DOI:10.1016/j.msea.2004.06.028.
[15] VYAS G M,PINHO S T,ROBINSON P.Constitutive modelling of fibre-reinforced composites with unidirectional plies using a plasticity-based approach[J/OL].Composites Science and Technology,2011,71(08):1068-1074.DOI:10.1016/j.compscitech.2011.03.009.
[16] VAUGHAN T J,MCCARTHY C T.Micromechanical modelling of the transverse damage behaviour in fibre reinforced composites[J/OL].Composites Science and Technology,2011,71(3):388-396.DOI:10.1016/j.compscitech.2010.12.006.
[17] QUINSON R,PEREZ J,RINK M,et al.Yield criteria for amorphous glassy polymers[J/OL].Journal of Materials Science,1997,32(05):1371-1379.DOI:10.1023/A:1018525127466.
[18] VOGLER T J,KYRIAKIDES S.Inelastic behavior of an AS4/PEEK composite under combined transverse compression and shear.Part I:experiments[J/OL].International Journal of Plasticity,1999,15(8):783-806.DOI:10.1016/S0749-6419(99)00011-X.
[19] RODR
GUEZ M,MOLINA-ALDAREGUA J M,GONZLEZ C,et al.A methodology to measure the interface shear strength by means of the fiber push-in test[J/OL].Composites Science and Technology,2012,72(15):1924-1932.DOI:10.1016/j.compscitech.2012.08.011.
[20] DUGDALE D S.Yielding of steel sheets containing slits[J/OL].Journal of the Mechanics and Physics of Solids,1960,8(2):100-104.DOI:10.1016/0022-5096(60)90013-2.
[21] BENZEGGAGH M L,KENANE M.Measurement of mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites with mixed-mode bending apparatus[J/OL].Composites Science and Technology,1996,56(4):439-449.DOI:10.1016/0266-3538(96)00005-X.
[22] LI S.Boundary conditions for unit cells from periodic microstructures and their implications[J/OL].Composites Science and Technology,2008,68(09):1962-1974.DOI:10.1016/j.compscitech.2007.03.035.