一种预测孔边多裂纹扩展的方法
技术领域
本发明属于安全验证领域，尤其是涉及一种预测孔边多裂纹扩展的方法。
背景技术
飞机作为民航运输的承载主体，在交通运输系统中发挥着至关重要的作用，当飞机运营过程中发生事故时，不仅对航空公司等企业造成巨大的经济损失，更会威胁到公众的生命安全，损害公众利益。为提高民航运输的安全水平，保持飞机全寿命周期的固有适航性，对民航飞机进行持续、有效的安全性评估是一项必要工作。民用飞机安全性评估贯穿飞机全寿命周期，设计阶段需要开展系统安全性分析(System Safety Analysis，SSA)，确定飞机固有安全性水平；在运营阶段，民机面临着腐蚀、疲劳损伤、重着陆、雷击等风险事件，随着服役时间增加，飞机结构的风险水平呈现上升趋势。由各类风险事件导致的结构损伤一旦发生，可能会造成结构失效，引发灾难事故。
我国的民航制造业起步较晚，目前ARJ-21飞机进入持续载客运营阶段，C919飞机进入试飞取证阶段，国产宽体机项目开始有序推进，但是国产民机的研发和制造仍然处于起步阶段，还没有建立起体系化的运营风险评估方法和流程。目前国内尚未展开对民机运营阶段的持续安全性分析和评估的研究，随着越来越多的国产民机投入运营，这部分空白急需填补。
发明内容
有鉴于此，本发明旨在提出一种预测孔边多裂纹扩展的方法，以解决上述问题中的不足之处。
为达到上述目的，本发明的技术方案是这样实现的：
一种预测孔边多裂纹扩展的方法，包括以下步骤：
S1、设定多个孔，同时每个孔的两侧水平位置取两个点，形成以下的测试点：
A、B、C、D、E、F、···、N；
S2、利用蒙特卡洛仿真方法对上述的测试点裂纹扩展初始寿命进行模拟实验，生成各个测试点的裂纹扩展初始寿命：
N00,N01,N02,N03,N04,N06…Nn；
S3、定义初始裂纹寿命是裂纹从0发展到Xmm时的循环次数，基于实验测得初始裂纹对数寿命满足N0～N(μ,σ)的对数正态分布；
S4、设置指定循环次数Nd，Nd也称作为指定寿命，各个测试点的循环次数到达Nd时无论裂纹扩展情况如何，结束模拟循环；
S5、判断到达Nd时测试点的裂纹生成的情况，将测试点的裂纹扩展初始寿命与Nd做比较，如果小于Nd认定该位置会产生裂纹，如果大于Nd认定到指定寿命时，该位置不会产生裂纹；
S6、将测试点的裂纹初始寿命从小到大进行排序，裂纹扩展初始寿命设置为多个测试点的裂纹扩展初始寿命中初始寿命的最小值Nmin，设定该位置开始时长度为Xmm，如果此时有其余位置寿命相同，则同时设置为Xmm，其余寿命大于Nmin的位置记裂纹尺寸为0；
S7、将测试点的裂纹扩展初始寿命从小到大进行排序，依次以各个测试点的裂纹扩展初始寿命为基准，对其余的测试点进行步骤S6中的对比；
S8、对步骤S6和S7中得到的每个裂纹大于0的位置进行裂纹扩展计算，设置循环次数每增加Z，需要计算每个测试点裂纹增长长度da，利用下式：
其中C＝2.34×10-8，m＝3.427，dN即为Z；
S9、当每次内层循环得到各个裂纹的长度后，加入判断：
G＝305*(97.74-15-每个测试点的裂纹长度总和)/97.74；
S10、针对步骤S9中求得的G值，设定基于G值的G判断规则，根据G判断规则将标记结果标记为安全/不安全；
S11、统计安全与不安全的次数，显示每种不安全原因的次数，计算试件的失效概率，作为该程序的最终结果，失效概率的计算公式为：
其中N为总的外循环的次数。
进一步的，步骤S8中的ΔK的计算方法为：
当到达循环次数Nd时裂纹长度为0，则ΔK为0；
当到达循环次数Nd时裂纹长度大于0时：
ΔK＝β*108*(π*(2.5+a))0.5 (2)
其中a为孔边的测试点的循环次数到达循环次数Nd时的裂纹长度，β因子为B1、B2、B3、B4连乘形式，B1、B2、B3、B4对不同形式的裂纹，其系数不同。
进一步的，B1的计算方法：
其中，x＝(a+2.5)/2.5，a为当前裂纹长度，2.5表示孔半径，Fh2是孔边单裂纹时使用的系数B1；
B2的计算方法：
其中，x＝(a+2.5)/2.5，a为当前裂纹长度，2.5表示孔半径，Fh1是孔边双裂纹时使用的系数B2；
B3的计算方法：
其中，x＝5/(15-1.5a)，FohA是相邻孔的影响因子B3；
B4的计算方法：
其中，x＝2*a1/(12-a1-a2)，a1为待求位置处的裂纹长度，a2为裂纹发展方向相邻的孔边裂纹长度；FOCA是相邻裂纹影响因子B4。
进一步的，在对B3进行计算时，引入B3判断方法：a超过规定长度时，在到达指定循环次数之前，提前结束本次内循环，结果标记为不安全，原因记为单裂纹长度超限；
在对B4进行计算时，引入B4判断方法：12-a1-a2低于规定长度时，在到达指定循环次数之前，提前结束本次内循环，结果标记为不安全，原因记为单裂纹长度超限。
相对于现有技术，本发明所述的一种预测孔边多裂纹扩展的方法具有以下有益效果：
变量的分布不需要近似，可以通过抽样手段完成变量模拟分布；可模拟相关性以及各个变量之间的依赖关系；可以通过增加循环次数来确保精确度；可迅速更改模型，进行变参数分析。在安全性评估领域，Monte-Carlo方法相比于故障树和马尔可夫等分析方法上的优势在于可以模拟更复杂的系统行为，例如多裂纹的随机扩展和结构失效行为，使用一般性的概率性方法，很难进行准确评估。所以采用Monte-Carlo方法进行孔边多裂纹的裂纹扩展失效的预测时，需对模拟实验的运算过程进行优化，在保证精度的情况下，减少运算量，提高运算效率。
附图说明
构成本发明的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
图1为实施例中取用的三个孔以及三个孔的测试点ABCDEF示意图；
图2为B1计算过程中的散点数据进行线性拟合过程示意图；
图3、4为B2计算过程中边界修正因子与几何参数的特定曲线函数关系示意图；
图5是一种预测孔边多裂纹扩展的方法流程图。
具体实施方式
需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。
下面将参考附图并结合实施例来详细说明本发明。
有内外两层循环，外层循环是蒙特卡洛模拟实验，主要用来生成初始的裂纹扩展寿命，本次模拟实验设置外层循环的循环次数为10万次(根据模拟条件、模拟精度可作调整)，具体流程如下：
如图1所示，设定3个孔，A B C D E F六个位置(每个孔的两侧水平位置)均可能产生裂纹，为了模拟多裂纹试样，采取蒙特卡洛仿真方法，每次外层循环，生成6个裂纹初始寿命(N00,N01,N02,N03,N04,N06)。初始裂纹寿命是裂纹从0发展到1mm时的循环次数，基于实验测得初始裂纹对数寿命满足N0～N(μ,σ)(实验测得μ＝4.969，σ＝0.326)的对数正态分布。该方法可用来探究给定循环次数下，多裂纹的扩展以及剩余强度情况，试件出现初始裂纹寿命到试件失效，寿命区间大概在7万到12万次，给定循环Nd范围大概在此区间即可(该给定循环参数可调)。
设置指定循环次数Nd＝86000，到达该次数时，无论裂纹扩展情况如何，结束该次内层循环。
判断到达指定寿命Nd时，6个位置裂纹生成的情况，用第一步MC方法抽取到的6个初始寿命与86000次做比较，如果小于86000次，认定该位置会产生裂纹，如果大于86000次，认定到指定寿命时，该位置不会产生裂纹，每个位置的裂纹长度始终为0。例如A＝110000、B＝80000、C＝90000、D＝110000、E＝85000、F＝92000，A、C、D、F四个位置认定为无裂纹产生，进行标记，后续不再计算(不计算裂纹增长长度da)。
将ABCDEF六个寿命从小到大进行排序，裂纹扩展初始寿命设置为6个裂纹中初始寿命的最小值80000，该位置开始时长度为1mm，如果此时有其余位置寿命相同，则同时设置为1mm，其余寿命大于80000的位置记裂纹尺寸为0。
往后循环次数增加，每次需要判断，如果循环次数大于等于该位置的初始寿命，则将该位置处的裂纹尺寸改为1。
进入裂纹扩展环节，对于裂纹大于0的位置(判断)，设置循环次数每增加1000，需要计算每个位置裂纹增长长度da，利用下式：
其中C＝2.34×10-8，m＝3.427，由于上方设置循环次数每增加1000就需要计算每个位置裂纹增长长度da，因此dN为1000；
ΔK的计算较为复杂，需要大量判断。具体的计算过程为：
(1)当裂纹长度为0时，ΔK为0；
(2)当裂纹长度大于0时：
ΔK＝β*108*(π*(2.5+a))0.5 (2)
其中a为孔边到设定循环次数的裂纹总长度，β因子为B1、B2、B3、B4连乘形式，B1、B2、B3、B4对不同形式的裂纹，其系数不同。具体计算如下：
1)B1的计算：
通过对图2中的散点数据进行非线性拟合，得到Fh2(x)关于a/R的函数(3)式：
其中，x＝(a+2.5)/2.5，a为当前裂纹长度，2.5表示孔半径，Fh2是孔边单裂纹时使用的系数B1。
2)B2的计算：
从图3可看出，边界修正因子Fh1与几何参数a/R具有特定的曲线函数关系。利用MATLAB对图3中的数据进行非线性拟合，得到关于Fh1的曲线函数表达式见公式(4)：
其中，x＝(a+2.5)/2.5，a为当前裂纹长度，2.5表示孔半径。Fh1是孔边双裂纹时使用的系数B2；
3)B3的计算：
其中，x＝5/(15-1.5a)，FohA是相邻孔的影响因子B3(裂纹扩展方向相临孔)，此处加入判断，当a超过7mm时，在到达本实施例中的指定循环次数(86000)之前，提前结束本次内循环，结果标记为不安全，原因记为单裂纹长度超限；
4)B4的计算：
其中，x＝2*a1/(12-a1-a2)，a1为待求位置处的裂纹长度，a2为裂纹发展方向相邻的孔边裂纹长度(由于A点和B点为同一个孔的边缘部位的两个测试点，例如图中的a1为A点处的裂纹长度，则a2为B点处的裂纹长度)。FOCA是相邻裂纹影响因子B4(裂纹扩展方向的孔边裂纹)。
此处加入判断，12-a1-a2低于2mm时，在到达本实施例中指定循环次数(86000)之前，提前结束本次内循环，结果标记为不安全，原因记为单裂纹长度超限。
β因子的具体计算方式：
孔边某位置如果此时裂纹为0，则不计算，或令k＝0；
1)对于两侧的A、F两个位置有裂纹，即两侧的裂纹(该位置裂纹长度大于0)，如果孔的另一侧，裂纹长度等于0，即与AF对应的BE侧无裂纹，则使用以下公式：
β＝B1 (7)
如果孔的另一侧，裂纹长度大于0(也就是说如图2中的AF有裂纹，BE也有裂纹)，即已有萌生的裂纹，则：
β＝B2 (8)
2)对于B、E两个边孔的内侧位置初裂纹，首先判断对应的A、F位置如果无裂纹，则使用B1；有裂纹则使用B2，无论何种情况，B3此时一定要有，如果对应的C或D处无裂纹，则不乘B4：即β＝B1*B3或β＝B2*B3，如果对应的C、D位置有裂纹则需要再乘以B4，即β＝B1*B3*B4或β＝B2*B3*B4。
3)对于C、D处的裂纹，同B、E相同，考虑对应C、D位置的裂纹情况，β＝B1*B3或β＝B2*B3或β＝B1*B3*B4或β＝B2*B3*B4。
将计算得到的参数β代入下式：
ΔK＝β*108*(π*(2.5+a))0.5 (9)
得到ΔK，将各个位置处的ΔK记录下来，并最终代入裂纹扩展公式(1)。
每次循环次数增长1000，即dN＝1000计算出各个裂对应的da，与已有的裂纹长度相加，得到每次循环后各个裂纹的aj＝aj+daj。
当每次内层循环得到各个裂纹的长度后，加入判断：
G＝305*(97.74-15-所有裂纹的长度总和(从A到F位置))/97.74 (10)
如果G大于150则将结果标记为安全；
如果到达指定循环次数86000后，G仍大于150，则该次内循环的结果为安全；
如果在到达指定循环次数(86000)之前，G值已经低于150，结束本次内循环计算，标记结果为不安全.
如果达到指定循环次数标记结果为不安全，记录为一次失效ni＝1；如果标记结果为安全，则记录ni＝0。
关于G的说明：共有三个判断，如果不满足任何一个判断条件，结果均标记为不安全，原因记为剩余强度过低；只有三个条件都满足时，到达指定循环次数86000后，本次内循环结果为安全。
经过G判断结束，当外层循环次数到达时，统计每次内层循环的结果，统计安全与不安全的次数，显示每种不安全原因的次数，计算试件的失效概率，作为该程序的最终结果。
失效概率的计算公式为：
其中N为总的外循环的次数。
本领域普通技术人员可以意识到，结合本文中所公开的实施例描述的各示例的单元及方法步骤，能够以电子硬件、计算机软件或者二者的结合来实现，为了清楚地说明硬件和软件的可互换性，在上述说明中已经按照功能一般性地描述了各示例的组成及步骤。这些功能究竟以硬件还是软件方式来执行，取决于技术方案的特定应用和设计约束条件。专业技术人员可以对每个特定的应用来使用不同方法来实现所描述的功能，但是这种实现不应认为超出本发明的范围。
在本申请所提供的几个实施例中，应该理解到，所揭露的方法和系统，可以通过其它的方式实现。例如，以上所述单元的划分，仅仅为一种逻辑功能划分，实际实现时可以有另外的划分方式，例如多个单元或组件可以结合或者可以集成到另一个系统，或一些特征可以忽略，或不执行。上述单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部单元来实现本发明实施例方案的目的。
最后应说明的是：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围，其均应涵盖在本发明的权利要求和说明书的范围当中。
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
- 上一篇：石墨接头机器人自动装卡簧、装栓机
- 下一篇：一种基于AIS数据的船舶水污染预测方法
