一种含孔复合材料层合板孔周应力分布的计算方法

技术领域

本发明属于复合材料结构力学与结构设计领域，具体涉及一种含孔复合材料层合 板孔周应力分布的计算方法，用于精确分析和计算平面应力状态下对称铺层的含孔复 合材料层合板孔周的平均应力与各单层主方向应力分布状态。

背景技术

碳纤维增强复合材料(CarbonFiberReinforcedPolymer，CFRP)由于其高比强度、 高比刚度、耐疲劳和易设计性等特点，作为零部件能有效地降低飞机、运载火箭、导 弹和卫星的结构重量，在航空航天领域中得到了广泛的应用。复合材料在航空器上的 应用部位，经历了由次承力结构向主承力结构过渡的过程，其用量的多少成为衡量技 术是否先进的重要标准之一。

在实际生产中，由于减重、连接、安装、管路铺设、检查维护等功能需要，复合 材料层合板上不可避免地含有开孔结构，如机翼上便于安装和检修的口盖，复合材料 翼梁上通过电气液压管路的开口等。开口打断了纤维和基体在该区域的连续性，影响 了结构的完整性，使得孔结构成为零部件的薄弱部位。在外载荷的作用下，孔周区域 极易发生应力集中现象，使得层合板的刚度降低、承载能力减小，严重时还会导致复 合材料零部件功能失效，轻者缩短零部件使用寿命，重者造成安全事故。

复合材料作为各向异性材料，其纤维与基体两相结合带来的非均匀性和各向异性 特性使得含孔结构应力分析较为复杂，层合板的材料属性、几何尺寸、开孔形状、加 载方式、铺层方向、铺层顺序等对孔周应力分布均有不同程度的影响，且层合板在厚 度方向上应分布不均匀，各单层应力分布不连续，针对传统金属材料的应力分析方法 不再适用。

众多学者针对含孔复合材料层合板结构孔周应力分布开展了许多研究。德国的 Kirsch首先得出圆孔附近应力集中的结果(G.Kirsch,Dietheoriederunddie bedürfnissederfestigkeitslehre,Springer,1898)，俄国的Klossowski求出椭圆孔附近应力 集中的公式(Timoshenko,S.,Woinowsky-Krieger,S.,Theoryofplatesandshells, McGrew-HillNewYork,1959)，这两项研究主要适用于各项同性材料，但为各向异性 材料研究做了铺垫；20世纪20年代末，苏联Muskhelishvili等人把复变函数引入弹性 力学，用保角映射把一个由不规则分段光滑的曲线构成的孔变换到单位圆上，导出复 变函数的应力表达式及其边界条件，进而获得一批应力集中的精确解(N.I. Muskhelishvili,Somebasicproblemsofthemathematicaltheoryofelasticity；fundamental equations,planetheoryofelasticity,torsion,andbending,P.Noordhoff,Groningen,1963)； Lekhnitskii对各向异性弹性力学进行了系统性的研究，详细论述了对于各向异性材料 在不同载荷下的应力分析方法(S.G.Lekhnitskii,Anisotropicplates,GordonandBreach, NewYork,1968.)；Kirchhoff-Love基于“直法线假设”、“法线长度保持不变”、“z向应 力可以忽略”这三个假设提出了板壳理论(Reddy,J.N.,Theoryandanalysisofelastic platesandshells,CRCPressandFrancis,2007)，得到薄板中面变形方程，成为研究复合 材料的经典层合板理论基础。这些都构成了含孔复合材料层合板应力分析和后续研究 的理论基础。

在实际工程应用中，复合材料层合板根据功能设定需要开口，而开口后层合板的 孔周应力集中水平与承载能力通常由实验或是有限元数值模拟计算，从而存在一些问 题：1)需要做多组对比实验，时间和财力耗费较大；2)数值计算时间成本大，而且 计算时间随计算精度的提高而成指数级增大；3)精确度浮动性比较大，建模水平、边 界条件和载荷工况的真实性直接影响计算准确度；3)试验和有限元数值模拟均未深入 研究各项因素对孔周各单层应力分布的影响规律。

发明内容

要解决的技术问题

本发明旨在至少在一定程度上解决上述技术问题之一或至少提供一种有用的商业 选择。为此，本发明的一个目的在于通过发明人在探索开口椭圆程度、铺层顺序、铺 层方向、铺层厚度以及载荷形式对孔周应力分布的影响规律的研究基础上，提出一种 含孔复合材料层合板孔周应力分布的计算方法，准确衡量复合材料层合薄板开口孔周 各单层的应力水平，继而能够为寻求减弱开口边缘应力集中程度的方法服务，增强含 孔薄板类零部件承载能力。

技术方案

本发明提出的一种含孔复合材料层合板孔周应力分布的计算方法，基于 Muskhelishvili的弹性力学复势理论与复合材料断裂力学分析方法，构建层合板孔周应 力分布模型，实现各单层孔周主方向应力的精确计算，为复合材料层合板开口设计和 后续的孔周损伤与疲劳分析奠定基础。

所述一种含孔复合材料层合板孔周应力分布的计算方法，其特征在于：包括以下 步骤：

步骤1：将含孔复合材料层合板的面内受力转换为在坐标系x′oy′下的面内应力载 荷，使得在坐标系x′oy′下的剪切应力载荷得到ox′方向面内应力载荷oy′方向面内应力载荷坐标系x′oy′与层合板全局坐标系xoy的偏转角为β；

步骤2：将坐标系x′oy′下的面内应力载荷转换到全局坐标系xoy下，得到全局坐标 系xoy下的等效面内应力载荷：

${\sigma }_{x,y}^{\infty }=\left(\begin{array}{c}{\sigma }_x^{\infty }\\ {\sigma }_y^{\infty }\\ {\tau }_{xy}^{\infty }\end{array}\right)=\left(\begin{array}{c}\frac{\sigma }{2}[(\lambda +1)+(\lambda -1)cos2\beta ]\\ \frac{\sigma }{2}[(\lambda +1)-(\lambda -1)\cos 2\beta ]\\ \frac{\sigma }{2}\left[(\lambda -1)sin2\beta \right]\end{array}\right);$

步骤3：根据步骤2得到的全局坐标系xoy下的等效面内应力载荷，通过以下公式 得到全局坐标系xoy下的层合板厚度方向孔周各点平均应力

$\left(\begin{array}{c}{\bar{\sigma }}_x\\ {\bar{\sigma }}_y\\ {\bar{\tau }}_{xy}\end{array}\right)=\left(\begin{array}{c}{\sigma }_x^{\infty }+2\mathrm{Re}\left[\frac{{\mu }_1^2{A}_1}{p_1-q_1(cos2\theta +i·\sin 2\theta )}+\frac{{\mu }_2^2{B}_1}{q_2-p_2(\cos 2\theta +i·\sin 2\theta )}\right]\\ {\sigma }_y^{\infty }+2\mathrm{Re}\left[\frac{A_1}{p_1-q_1(\cos 2\theta +i·\sin 2\theta )}+\frac{B_1}{q_2-p_2(\cos 2\theta +i·\sin 2\theta )}\right]\\ {\tau }_{xy}^{\infty }-2\mathrm{Re}\left[\frac{{\mu }_1{A}_1}{p_1-q_1(cos2\theta +i·\sin 2\theta )}+\frac{{\mu }_2{B}_1}{q_2-p_2(\cos 2\theta +i·\sin 2\theta )}\right]\end{array}\right)$

其中Re表示取实部，i为虚数标记，θ为孔周位置角；

μ₁、μ₂为层合板变形协调特征方程的特征根，其中μ₁＝α₁+i·β₁，μ₂＝α₂+i·β₂，且 β₁＞0,β₂＞0，层合板变形协调特征方程为

a₁₁μ⁴-2a₁₆μ³+(2a₁₂+a₆₆)μ²-2a₂₆μ+a₂₂＝0

层合板变形协调特征方程中的系数为层合板整体等效柔度矩阵的元素，层合板整体等 效柔度矩阵为

$\left(\begin{array}{ccc}{a}_{11}& a_{12}& a_{16}\\ a_{12}& a_{22}& a_{26}\\ a_{16}& a_{26}& a_{66}\end{array}\right);$

参数A₁、B₁为

$A_1=\frac{[{\mu }_2(K_1+{\bar{K}}_2)-(K_3+{\bar{K}}_4)]+m[{\mu }_2({\bar{K}}_1+K_2)-({\bar{K}}_3+K_4)]}{{\mu }_1-{\mu }_2},$

$B_1=\frac{[{\mu }_1(K_1+{\bar{K}}_2)-(K_3+{\bar{K}}_4)]+m[{\mu }_1({\bar{K}}_1+K_2)-({\bar{K}}_3+K_4)]}{{\mu }_1-{\mu }_2}$

其中m为层合板孔的椭圆程度，a、b分别为层合板孔的长半轴与短 半轴；K₁、K₂、K₃、K₄为：

$K_1=\frac{1}{2}[B(1+i·{\mu }_1)+(B^{\prime }+i·C^{\prime })(1+i·{\mu }_2)]$

$K_2=\frac{1}{2}[B(1-i·{\mu }_1)+(B^{\prime }+i·C^{\prime })(1-i·{\mu }_2)],$

$K_3=\frac{1}{2}[{\mu }_{1}B(1+i·{\mu }_1)+{\mu }_2(B^{\prime }+i·C^{\prime })(1+i·{\mu }_2)]$

$K_4=\frac{1}{2}[{\mu }_{1}B(1-i·{\mu }_1)+{\mu }_2(B^{\prime }+i·C^{\prime })(1-i·{\mu }_2)]$

对应分别为K₁、K₂、K₃、K₄的共轭；参数B,B′,C′为：

$B=\frac{{\sigma }_x^{\infty }+({\alpha }_2^2+{\beta }_2^2){\sigma }_y^{\infty }+2{\alpha }_2{\tau }_{xy}^{\infty }}{2[{({\alpha }_2-{\alpha }_1)}^2+({\beta }_2^2-{\beta }_1^2)]}$

$B^{\prime }=\frac{({\alpha }_1^2-{\beta }_1^2-2{\alpha }_1{\alpha }_2){\sigma }_y^{\infty }-{\sigma }_x^{\infty }-2{\alpha }_2{\tau }_{xy}^{\infty }}{2[{({\alpha }_2-{\alpha }_1)}^2+({\beta }_2^2-{\beta }_1^2)]};$

$C^{\prime }=\frac{({\alpha }_1-{\alpha }_2){\sigma }_x^{\infty }+\left[{\alpha }_2\right({\alpha }_1^2-{\beta }_1^2)-{\alpha }_1({\alpha }_2^2-{\beta }_2^2\left)\right]{\sigma }_y^{\infty }+\left[\right({\alpha }_1^2-{\beta }_1^2)-({\alpha }_2^2-{\beta }_2^2\left)\right]{\tau }_{xy}^{\infty }}{2{\beta }_2[{({\alpha }_2-{\alpha }_1)}^2+({\beta }_2^2-{\beta }_1^2)]}$

参数p₁、p₂、q₁、q₂为：

$p_1=\frac{(1+m)+i·{\mu }_1(1-m)}{2},p_2=\frac{(1+m)+i·{\mu }_2(1-m)}{2},$

$q_1=\frac{(1+m)-i·{\mu }_1(1-m)}{2},q_2=\frac{(1+m)-i·{\mu }_2(1-m)}{2},$

步骤4：根据步骤3得到的全局坐标系xoy下的层合板厚度方向孔周各点平均应力， 得到复合材料层合板各层孔周主方向上的主应力，其中第k层孔周主方向上的主应力 为${\left(\begin{array}{c}{\sigma }_L\\ {\sigma }_T\\ {\tau }_{LT}\end{array}\right)}_k:$

${\left(\begin{array}{c}{\sigma }_L\\ {\sigma }_T\\ {\tau }_{LT}\end{array}\right)}_k={\left(\begin{array}{ccc}{\cos }^2\alpha & {\sin }^2\alpha & 2\sin \alpha \cos \alpha \\ {\sin }^2\alpha & {\cos }^2\alpha & -2\sin \alpha \cos \alpha \\ -\sin \alpha \cos \alpha & \sin \alpha \cos \alpha & {\cos }^2\alpha -{\sin }^2\alpha \end{array}\right)}_k{\left(\begin{array}{ccc}{\bar{Q}}_{11}& {\bar{Q}}_{12}& {\bar{Q}}_{16}\\ {\bar{Q}}_{12}& {\bar{Q}}_{22}& {\bar{Q}}_{26}\\ {\bar{Q}}_{16}& {\bar{Q}}_{26}& {\bar{Q}}_{66}\end{array}\right)}_k\left(\begin{array}{ccc}{a}_{11}& a_{12}& a_{16}\\ a_{12}& a_{22}& a_{26}\\ a_{16}& a_{26}& a_{66}\end{array}\right)\left(\begin{array}{c}{\bar{\sigma }}_x\\ {\bar{\sigma }}_y\\ {\bar{\tau }}_{xy}\end{array}\right),$其中矩阵${\left(\begin{array}{ccc}{\cos }^2\alpha & {\sin }^2\alpha & 2\sin \alpha cos\alpha \\ {\sin }^2\alpha & {\cos }^2\alpha & -2sin\alpha cos\alpha \\ -\sin \alpha \cos \alpha & \sin \alpha \cos \alpha & {\cos }^2\alpha -{\sin }^2\alpha \end{array}\right)}_k$为由全局坐标系到复合材料层合板第 k层局部坐标系的应力转轴矩阵，其中复合材料层合板第k层局部坐标系与全局坐标系 的偏转角为复合材料层合板第k层的铺层方向角α；矩阵${\left(\begin{array}{ccc}{\bar{Q}}_{11}& {\bar{Q}}_{12}& {\bar{Q}}_{16}\\ {\bar{Q}}_{12}& {\bar{Q}}_{22}& {\bar{Q}}_{26}\\ {\bar{Q}}_{16}& {\bar{Q}}_{26}& {\bar{Q}}_{66}\end{array}\right)}_k$为复合材料层 合板第k层的等效刚度矩阵。

根据本发明的实施例，该方法基于断裂力学的复变方法建立，方法计算量小、成 本低廉、计算精度高、成本低、速度快；而且该方法不仅仅能够计算复合材料层合板 上的圆孔和椭圆孔的孔周平均应力与各单层主方向上的应力，而且在面对复合材料层 合板上的不规则孔时，只需要将不规则孔边界投影到复平面的单位圆上，即可通过本 方法计算相应的孔周平均应力与各单层主方向上的应力，适用范围广泛，能够满足实 际复合材料结构设计与应用过程中的实际需要。

附图说明

本发明的上述和/或附加的方面和优点从结合下面附图对实施例的描述中将变得 明显和容易理解，其中：

图1含孔复合材料层合板规格及标记；

图2含孔复合材料层合板平面应力状态；

图3层合板各单层面内孔周各点L向主应力曲线；

图4层合板各单层面内孔周各点T向主应力曲线；

图5层合板各单层面内孔周各点切应力曲线。

具体实施方式

下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终 相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参 考附图描述的实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限 制。

在本发明的描述中，需要理解的是，术语“中心”、“纵向”、“横向”、“长度”、“宽 度”、“厚度”、“上”、“下”、“前”、“后”、“左”、“右”、“竖直”、“水平”、“顶”、“底”、 “内”、“外”、“顺时针”、“逆时针”等指示的方位或位置关系为基于附图所示的方位 或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或 元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限 制。

此外、术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重 要性或者隐含指明所指示的技术特征的数量。因此，限定有“第一”、“第二”的特征 可以明示或者隐含地包括一个或者更多个该特征。在本发明的描述中，“多个”的含义 是两个或两个以上，除非另有明确具体的限定。

实际生产中常用的复合材料层合板如图1所示。为避免铺层不对称造成的翘曲， 层合板为正交对称铺层结构，板料最小边长为w，单层厚度为t_k，层合板整体厚度为h。 层合板中心含有功能性的椭圆形开口，且边长w与开口包络圆直径d满足w/d≥4。开 口的椭圆程度由系数m表示，定义为：

$m=\frac{a-b}{a+b}(a\ge b)$

参数a、b分别为椭圆形开口的长半轴与短半轴。特别地，a＝b时开口为圆，a＞＞b时 开口为狭缝。

复合材料层合板面内尺寸远大于厚度尺寸，由此层合板可看成是薄板，在面内的 受力状态可简化为平面应力状态。所以可以将含孔复合材料层合板的面内受力转换为 在某一坐标系x′oy′下的面内应力载荷，使得在坐标系x′oy′下的剪切应力载荷得到ox′方向面内应力载荷oy′方向面内应力载荷坐标系x′oy′与层合 板全局坐标系xoy的偏转角为β，如图2所示。

并可将坐标系x′oy′下的面内应力载荷转换到全局坐标系xoy下，得到全局坐标系 xoy下的等效面内应力载荷：

${\sigma }_{x,y}^{\infty }=\left(\begin{array}{c}{\sigma }_x^{\infty }\\ {\sigma }_y^{\infty }\\ {\tau }_{xy}^{\infty }\end{array}\right)=\left(\begin{array}{c}\frac{\sigma }{2}[(\lambda +1)+(\lambda -1)cos2\beta ]\\ \frac{\sigma }{2}[(\lambda +1)-(\lambda -1)\cos 2\beta ]\\ \frac{\sigma }{2}\left[(\lambda -1)sin2\beta \right]\end{array}\right);$

参数λ和β决定应力载荷形式。特别地，当λ＝0时，为单轴载荷，当λ≠0,β≠0时 为双轴载荷；当λ＝-1,β＝π/4或3π/4为剪应力载荷。

本实施例中碳纤维增强复合材料层合板铺层顺序为[+45/0/-60/90]_s，所用的碳纤维 环氧树脂力学性能参数如表1所示，单层厚度t_k＝0.25mm，总厚度h＝2mm，板中心含 有圆孔m＝0，孔径比为w/d≥4，受到全局坐标系下x轴方向上的单轴拉伸应力载荷， 大小为${\sigma }_{x^{\prime }}^{\infty }=10Mpa.$

表1碳纤维/环氧树脂(T700/3501-6)力学性能

层合板受到x轴方向上的单轴拉伸应力载荷取λ＝0，β＝π/2，得到全局坐 标系下的等效面内应力载荷为：

${\sigma }_{x,y}^{\infty }=\left(\begin{array}{c}{\sigma }_x^{\infty }\\ {\sigma }_y^{\infty }\\ {\tau }_{xy}^{\infty }\end{array}\right)=\left(\begin{array}{c}10\\ 0\\ 0\end{array}\right)Mpa.$

在得到全局坐标系xoy下的等效面内应力载荷后，通过以下公式得到全局坐标系 xoy下的层合板厚度方向孔周各点平均应力

$\left(\begin{array}{c}{\bar{\sigma }}_x\\ {\bar{\sigma }}_y\\ {\bar{\tau }}_{xy}\end{array}\right)=\left(\begin{array}{c}{\sigma }_x^{\infty }+2\mathrm{Re}\left[\frac{{\mu }_1^2{A}_1}{p_1-q_1(cos2\theta +i·\sin 2\theta )}+\frac{{\mu }_2^2{B}_1}{q_2-p_2(\cos 2\theta +i·\sin 2\theta )}\right]\\ {\sigma }_y^{\infty }+2\mathrm{Re}\left[\frac{A_1}{p_1-q_1(\cos 2\theta +i·\sin 2\theta )}+\frac{B_1}{q_2-p_2(\cos 2\theta +i·\sin 2\theta )}\right]\\ {\tau }_{xy}^{\infty }-2\mathrm{Re}\left[\frac{{\mu }_1{A}_1}{p_1-q_1(cos2\theta +i·\sin 2\theta )}+\frac{{\mu }_2{B}_1}{q_2-p_2(\cos 2\theta +i·\sin 2\theta )}\right]\end{array}\right)$

其中Re表示取实部，i为虚数标记，θ为孔周位置角。

μ₁、μ₂为层合板变形协调特征方程的特征根。

层合板变形协调特征方程为a₁₁μ⁴-2a₁₆μ³+(2a₁₂+a₆₆)μ²-2a₂₆μ+a₂₂＝0，通过以下 过程推导得到并求解：

含孔复合材料层合板受到面内载荷时，其变形在弹性限度范围内，层合板平均应 力应变关系为：

$\left(\begin{array}{c}{\overline{ϵ}}_x\\ {\overline{ϵ}}_y\\ {\bar{\gamma }}_{xy}\end{array}\right)=\left(\begin{array}{ccc}{a}_{11}& a_{12}& a_{16}\\ a_{12}& a_{22}& a_{26}\\ a_{16}& a_{26}& a_{66}\end{array}\right)\left(\begin{array}{c}{\bar{\sigma }}_x\\ {\bar{\sigma }}_y\\ {\bar{\tau }}_{xy}\end{array}\right)$

其中矩阵$\left(\begin{array}{ccc}{a}_{11}& a_{12}& a_{16}\\ a_{12}& a_{22}& a_{26}\\ a_{16}& a_{26}& a_{66}\end{array}\right)$为层合板整体等效柔度矩阵。

应用弹性平面力学复变函数方法，层合板厚度方向上的平均应力以Airy应力函数 U(x,y)的二阶倒数表示：

${\bar{\sigma }}_x=\frac{{\partial }^{2}U}{\partial y^2},{\bar{\sigma }}_y=\frac{{\partial }^{2}U}{\partial x^2},{\bar{\tau }}_{xy}=\frac{{\partial }^{2}U}{\partial x\partial y}$

层合板在外载荷作用下会产生变形，为保证变形后仍保持整体性和连续性，其平 均应力分量应满足应变协调方程：

$\frac{{\partial }^2{\overline{ϵ}}_x}{\partial y^2}+\frac{{\partial }^2{\overline{ϵ}}_y}{\partial x^2}=\frac{{\partial }^2{\bar{\gamma }}_{xy}}{\partial x\partial y}$

由此得到一个四阶齐次偏微分方程：

$a_{11}\frac{{\partial }^{4}U}{\partial x^4}-2a_{16}\frac{{\partial }^{4}U}{\partial x^3\partial y}+(2a_{12}+a_{66})\frac{{\partial }^{4}U}{\partial x^2\partial y^2}-2a_{26}\frac{{\partial }^{4}U}{\partial x\partial y^3}+a_{22}\frac{{\partial }^{4}U}{\partial y^4}=0$

该方程通解依赖于层合板变形协调特征方程的根，将本实施例中复合材料层合板的几 何及力学性能参数带入层合板变形协调特征方程，得到特征根：

μ₁＝0.58i+0.15，μ₂＝1.38i+0.43

${\bar{\mu }}_1=-0.58i+0.15,{\bar{\mu }}_2=-1.38i+0.43$

由于接下来还需要用到特定特征根的实部和虚部，即取其中其中μ₁＝α₁+i·β₁， μ₂＝α₂+i·β₂，且β₁＞0,β₂＞0，对应本实施例，α₁、α₂、β₁、β₂按照约定取值为： α₁＝0.15，α₂＝0.43，β₁＝0.58，β₂＝1.38，用以表征复合材料层合板的各向异性程度。

表征应力函数内边界条件的复参数A₁、B₁为

$A_1=\frac{[{\mu }_2(K_1+{\bar{K}}_2)-(K_3+{\bar{K}}_4)]+m[{\mu }_2({\bar{K}}_1+K_2)-({\bar{K}}_3+K_4)]}{{\mu }_1-{\mu }_2},$

$B_1=\frac{[{\mu }_1(K_1+{\bar{K}}_2)-(K_3+{\bar{K}}_4)]+m[{\mu }_1({\bar{K}}_1+K_2)-({\bar{K}}_3+K_4)]}{{\mu }_1-{\mu }_2}$

其中K₁、K₂、K₃、K₄为：

$K_1=\frac{1}{2}[B(1+i·{\mu }_1)+(B^{\prime }+i·C^{\prime })(1+i·{\mu }_2)]$

$K_2=\frac{1}{2}[B(1-i·{\mu }_1)+(B^{\prime }+i·C^{\prime })(1-i·{\mu }_2)],$

$K_3=\frac{1}{2}[{\mu }_{1}B(1+i·{\mu }_1)+{\mu }_2(B^{\prime }+i·C^{\prime })(1+i·{\mu }_2)]$

$K_4=\frac{1}{2}[{\mu }_{1}B(1-i·{\mu }_1)+{\mu }_2(B^{\prime }+i·C^{\prime })(1-i·{\mu }_2)]$

对应分别为K₁、K₂、K₃、K₄的共轭；参数B,B′,C′为：

$B=\frac{{\sigma }_x^{\infty }+({\alpha }_2^2+{\beta }_2^2){\sigma }_y^{\infty }+2{\alpha }_2{\tau }_{xy}^{\infty }}{2[{({\alpha }_2-{\alpha }_1)}^2+({\beta }_2^2-{\beta }_1^2)]}$

$B^{\prime }=\frac{({\alpha }_1^2-{\beta }_1^2-2{\alpha }_1{\alpha }_2){\sigma }_y^{\infty }-{\sigma }_x^{\infty }-2{\alpha }_2{\tau }_{xy}^{\infty }}{2[{({\alpha }_2-{\alpha }_1)}^2+({\beta }_2^2-{\beta }_1^2)]};$

$C^{\prime }=\frac{({\alpha }_1-{\alpha }_2){\sigma }_x^{\infty }+\left[{\alpha }_2\right({\alpha }_1^2-{\beta }_1^2)-{\alpha }_1({\alpha }_2^2-{\beta }_2^2\left)\right]{\sigma }_y^{\infty }+\left[\right({\alpha }_1^2-{\beta }_1^2)-({\alpha }_2^2-{\beta }_2^2\left)\right]{\tau }_{xy}^{\infty }}{2{\beta }_2[{({\alpha }_2-{\alpha }_1)}^2+({\beta }_2^2-{\beta }_1^2)]}$

参数p₁、p₂、q₁、q₂为：

$p_1=\frac{(1+m)+i·{\mu }_1(1-m)}{2},p_2=\frac{(1+m)+i·{\mu }_2(1-m)}{2},.$

$q_1=\frac{(1+m)-i·{\mu }_1(1-m)}{2},q_2=\frac{(1+m)-i·{\mu }_2(1-m)}{2},$

上面直接给出了化简之后的相关公式，发明人实际推导过程为：

根据Lekhnitskii的各向异性板理论，层合板厚度方向上的平均应力表示为：

其中的应力解析函数ψ(z₂)可以表示成如下形式：

$\psi \left(z_2\right)={\gamma }^{\prime }{\mathrm{lnz}}_2+(B^{\prime }+{\mathrm{iC}}^{\prime })z_2+\underset{n=1}{\overset{\infty }{\Sigma }}\frac{B_n}{z_2^n}$

在层合板外边界上，平均应力分量与所施加的应力载荷平衡，即有：

$\underset{x\to \infty }{\lim }{\bar{\sigma }}_x={\sigma }_x^{\infty };\underset{x\to \infty }{\lim }{\bar{\sigma }}_y={\sigma }_x^{\infty };\underset{x\to \infty ,y\to \infty }{\lim }{\bar{\tau }}_y={\tau }_{xy}^{\infty }$

从得到层合板厚度方向上的平均应力为：

${\bar{\sigma }}_x{|}_{\infty }=2\mathrm{Re}[{\mu }_1^{2}B+{\mu }_2^2(B^{\prime }+{\mathrm{iC}}^{\prime })]$

${\bar{\sigma }}_y{|}_{\infty }=2\mathrm{Re}[B+(B^{\prime }+{\mathrm{iC}}^{\prime })]$

${\bar{\tau }}_{xy}{|}_{\infty }=-2\mathrm{Re}[{\mu }_{1}B+{\mu }_2(B^{\prime }+{\mathrm{iC}}^{\prime })]$

该方程即对应参数B,B′,C′的表达式，将本实施例中的外边界载荷条件与特征根参数带 入后得到参数B、B′和C′的值：B＝2.61，B′＝-2.61，C′＝1.10。

内边界上，通过保角映射方法，将xoy平面上层合板内椭圆边界投影到复平面z上 的单位圆求解，各向同性情形下的映射函数为：

z＝ω(ζ)

复合材料层合板为各向异性材料，需要将各向同性映射函数改写为复平面内的各向异 性映射函数：

$z_j={\omega }_j\left(\zeta \right)=\frac{1}{2}(1-i·{\mu }_j)\omega \left(\zeta \right)+\frac{1}{2}(1+i·{\mu }_j)\omega \left(\frac{1}{\zeta }\right),(j=1,2)$

其中ζ＝e^iθ，θ∈[0°，360°)，同样的应力解析函数在复平面内的表达式为：

$\psi \left(\zeta \right)=\psi \left(z_2\right)={\gamma }^{\prime }\ln \left[{\omega }_2\left(\zeta \right)\right]+(B^{\prime }+{\mathrm{iC}}^{\prime })\left[{\omega }_2\left(\zeta \right)\right]+\underset{k=1}{\overset{n}{\Sigma }}{B}_k{\zeta }^{-k}$

在弹性力学中，平面问题的内边界条件为：

本例涉及的含孔层合板孔周为自由边界，孔周无牵引力作用，则有：

$\gamma =0,{\gamma }^{\prime }=0;\bar{X}=0,\bar{Y}=0$

所以得到内边界条件方程有

$2\mathrm{Re}[K_2(\zeta +\frac{m_1}{\zeta })+K_1(\frac{1}{\zeta }+m_1\zeta )+\underset{k=1}{\overset{n}{\Sigma }}{A}_k{\zeta }^{-k}+\underset{k=1}{\overset{n}{\Sigma }}{B}_k{\zeta }^{-k}]=0$

$2\mathrm{Re}[K_4(\zeta +\frac{m_1}{\zeta })+K_3(\frac{1}{\zeta }+m_1\zeta )+{\mu }_1\underset{k=1}{\overset{n}{\Sigma }}{B}_k{\zeta }^{-k}+{\mu }_2\underset{k=1}{\overset{n}{\Sigma }}{B}_k{\zeta }^{-k}]=0$

求解该方程，得到参数K₁、K₂、K₃、K₄的表达式，对应本实施例中：

K₁＝1.28+5.48×10^-1i，K₂＝-1.28+5.48×10^-1i

K₃＝-8.36×10^-1+1.22i，K₄＝8.36×10^-1-3.78i

并得到本实施例中，表征应力函数内边界条件的复参数A₁、B₁为

A₁＝4.0820-2.9532i，B₁＝-4.0820+2.9532i

求解过程中的辅助复参数p₁、p₂、q₁、q₂的值为

p₁＝2.10×10^-1+7.52×10^-2i，p₂＝1.19+2.16×10^-1i

q₁＝7.90×10^-1-7.52×10^-2i，q₂＝-1.19-2.16×10^-1i

将上述参数值带入层合板孔周厚度方向平均应力表达式，得到

${\bar{\sigma }}_x=2\mathrm{Re}\left(\begin{array}{c}\frac{-7.67\times {10}^{-1}+1.64i}{(2.10\times {10}^{-1}+7.52\times {10}^{-2}i)-(7.90\times {10}^{-1}-7.52\times {10}^{-2}i)·(cos2\theta +i·\sin 2\theta )}\\ +\frac{1.06\times {10}^1-2.35\times {10}^{-1}i}{(1.92\times {10}^{-1}+2.16\times {10}^{-1}i)-(-1.19-2.16\times {10}^{-1}i)·(cos2\theta +i·\sin 2\theta )}+5\end{array}\right)$

${\bar{\sigma }}_y=2\mathrm{Re}\left(\begin{array}{c}\frac{4.08-2.95i}{(2.10\times {10}^{-1}+7.52\times {10}^{-2}i)-(7.90\times {10}^{-1}-7.52\times {10}^{-2}i)·(cos2\theta +i·\sin 2\theta )}\\ +\frac{-4.08+2.95i}{(1.92\times {10}^{-1}+2.16\times {10}^{-1}i)-(-1.19-2.16\times {10}^{-1}i)·(cos2\theta +i·\sin 2\theta )}\end{array}\right).$

${\bar{\tau }}_{xy}=-2\mathrm{Re}\left(\begin{array}{c}\frac{2.33+1.93i}{(2.10\times {10}^{-1}+7.52\times {10}^{-2}i)-(7.90\times {10}^{-1}-7.52\times {10}^{-2}i)·(cos2\theta +i·\sin 2\theta )}\\ +\frac{-2.33-6.93i}{(1.92\times {10}^{-1}+2.16\times {10}^{-1}i)-(-1.19-2.16\times {10}^{-1}i)·(cos2\theta +i·\sin 2\theta )}\end{array}\right)$

根据得到的全局坐标系xoy下的层合板厚度方向孔周各点平均应力，得到复合材 料层合板各层孔周主方向上的主应力，其中第k层孔周主方向上的主应力为${\left(\begin{array}{c}{\sigma }_L\\ {\sigma }_T\\ {\tau }_{LT}\end{array}\right)}_k:$

${\left(\begin{array}{c}{\sigma }_L\\ {\sigma }_T\\ {\tau }_{LT}\end{array}\right)}_k={\left(\begin{array}{ccc}{\cos }^2\alpha & {\sin }^2\alpha & 2\sin \alpha \cos \alpha \\ {\sin }^2\alpha & {\cos }^2\alpha & -2\sin \alpha \cos \alpha \\ -\sin \alpha \cos \alpha & \sin \alpha \cos \alpha & {\cos }^2\alpha -{\sin }^2\alpha \end{array}\right)}_k{\left(\begin{array}{ccc}{\bar{Q}}_{11}& {\bar{Q}}_{12}& {\bar{Q}}_{16}\\ {\bar{Q}}_{12}& {\bar{Q}}_{22}& {\bar{Q}}_{26}\\ {\bar{Q}}_{16}& {\bar{Q}}_{26}& {\bar{Q}}_{66}\end{array}\right)}_k\left(\begin{array}{ccc}{a}_{11}& a_{12}& a_{16}\\ a_{12}& a_{22}& a_{26}\\ a_{16}& a_{26}& a_{66}\end{array}\right)\left(\begin{array}{c}{\bar{\sigma }}_x\\ {\bar{\sigma }}_y\\ {\bar{\tau }}_{xy}\end{array}\right),$其中矩阵${\left(\begin{array}{ccc}{\cos }^2\alpha & {\sin }^2\alpha & 2\sin \alpha cos\alpha \\ {\sin }^2\alpha & {\cos }^2\alpha & -2sin\alpha cos\alpha \\ -\sin \alpha \cos \alpha & \sin \alpha \cos \alpha & {\cos }^2\alpha -{\sin }^2\alpha \end{array}\right)}_k$为由全局坐标系到复合材料层合板第 k层局部坐标系的应力转轴矩阵，其中复合材料层合板第k层局部坐标系与全局坐标系 的偏转角为复合材料层合板第k层的铺层方向角α；矩阵${\left(\begin{array}{ccc}{\bar{Q}}_{11}& {\bar{Q}}_{12}& {\bar{Q}}_{16}\\ {\bar{Q}}_{12}& {\bar{Q}}_{22}& {\bar{Q}}_{26}\\ {\bar{Q}}_{16}& {\bar{Q}}_{26}& {\bar{Q}}_{66}\end{array}\right)}_k$为复合材料层 合板第k层的等效刚度矩阵。

本实施例中层合板45°层内主应力表达式为：

层合板0°层内主应力表达式为：

层合板-60°层内主应力表达式为：

层合板90°层内主应力表达式为：

由于层合板各单层在孔周位置角θ∈[0，360)内有两个周期，因此应力结果取 θ∈[0，180)内值，且以正则化的应力σ_L,T/σ来表示，并与有限元模拟结果进行对比。各 单层L向主应力σ_L、T向主应力σ_T、面内切应力τ_LT及其相应的有限元数值模拟结果分 别如图3、图4、图5所示。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例 性的，不能理解为对本发明的限制，本领域的普通技术人员在不脱离本发明的原理和 宗旨的情况下在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。