CN106874611A - 一种基于超体积迭代策略的含区间参数结构响应区间的分析方法 - Google Patents

Abstract

本发明公开了一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，本发明将结构响应量作为一个目标函数，通过优化算法求解响应量的最大/小值，从而确定响应量的上/下界。本发明采用的优化算法是一种基于超体积迭代策略全局优化算法，采用辛普森积分的策略在不确定参数范围内自适应地布置积分点，寻找这些积分点中对应目标函数的最大/小值，然后以目标函数最大/小值所对应的积分点为起始点，基于牛顿法得到全局的最大/小值，即得到了不确定性参数在非线性系统中传播的上/下界。该方法可以在不确定输入条件下，准确地确定结构的区间响应传播边界。

Description

一种基于超体积迭代策略的含区间参数结构响应区间的分析 方法

技术领域

本发明涉及含区间参数结构响应的不确定性传播领域，特别涉及一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，该方法是基于超体积迭代策略全局优化算法而提出的。

背景技术

现代化工程结构复杂而精细，但是由于制造工艺的限制、认知的匮乏、测量的误差以及环境因素的变化，不确定性问题广泛存在于工程问题中，因此在分析和设计时需要更加注重先进数值算法和设计技术的研究。对复杂系统而言，不确定因素会因为系统的非线性而使响应的区间上/下界较难确定，现在解决该问题的方法主要有区间参数摄动法、顶点组合法等。前者往往基于一阶Taylor级数展开理论，将问题在很小的范围内线性化以此达到简化问题的目的，但是其不足也显而易见，即只能针对非线性不太强的问题在较小不确定度的情况，对于机翼等中复合材料具有较大分散性的区间响应分析问题就会有较大误差。顶点组合法针对单调问题可以精确地得到响应区间上/下界，如果不是单调的问题，则该方法将不再适用。以上都是基于区间数学的理论方法，蒙特卡洛法也可以解决该问题并得到较为精确的区间上/下界，甚至可以可作为对比其他方法精确度的一个参考基准，该方法结果越精确意味着随机抽样次数越多，对于复杂结构将带来较大的计算负担。因此针对非线性较强且具有较大不确定度的问题，研究一种在保证一定计算效率的条件下可以精确地确定响应量区间上下界的方法，十分有意义。

本发明基于优化方法的思路将响应量作为优化的目标函数，通过求解其最大/小值即可得到区间的上/下界。这样将一个区间响应的不确定性传播转化为一个优化问题，不仅使问题简化，无需考虑区间数学的理论方法，同时也克服上述诸多问题。在诸多优化方法中，牛顿法由于其较快的收敛速度广泛应用于各优化问题中，但是牛顿法针对的是局部优化问题，当针对全局优化问题时，牛顿法以及其他基于梯度信息的算法将陷入局部最优解区域无法跳出，从而无法得到全局最优解。启发式算法是现在发展较好的一类全局优化算法，包括模拟退火法、遗传算法、粒子群算法等，该类方法通过模拟或揭示某些自然现象或过程而得到发展，其思想和内容涉及数学、物理学、生物进化、人工智能、统计力学等方面，为解决复杂问题提供新的思路和手段，但是该类方法属于随机算法，没有明确的收敛准则。本发明基于超体积迭代策略全局优化算法，该方法以目标函数在可行域的超体积作为迭代的收敛准则，可以快速准确地收敛到全局最优解。

发明内容

本发明要解决的技术问题是：克服现有方法和技术的不足，提供一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，本发明充分考虑了不确定性因素对结构响应带来的影响，建立了一种基于优化思路的求解不确定性传播的响应区间上/下界的方法，相较于传统的方法拓展了不确定性传播分析的方法，并可以快速准确地确定非线性较强、具有较大不确定度系统中结构响应区间的上/下界。

本发明采用的技术方案为：一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，具体流程见图1，其实现步骤如下：

步骤(1)基于MSC.Patran构建结构的有限元模型，包括建立几何模型，赋予材料属性，施加荷载和边界条件。

步骤(2)确定由步骤(1)构建的模型中的不确定参数x及其不确定区间。由于制作工艺、测量误差等因素导致结构中的材料和荷载具有一定的不确定性，结构中材料的弹性模量作用在结构上的荷载这里E和P分别表示弹性模量和荷载的下界，和分别表示弹性模量和荷载的上界，然后在MSC.Patran中修改相应不确定参数的数值，并且该修改过程需要用到MSC.Patran的宏录制功能最后形成记录该修改过程的ses文件。其中区间不确定性参数向量x可以表示为：

其中，x^c＝(E^c,P^c,…)为区间中点，x^r＝(E^r,P^r,…)为区间半径，e∈Ξⁿ，Ξⁿ定义为所有元素包含在[-1,1]内的n维向量集合，符号“×”定义为两个向量各对应元素乘积仍为n维向量。

步骤(3)修改结构有限元模型中对应的不确定参数。在步骤(2)中形成的ses文件中找到修改模型中不确定参数的关键性命令，利用MATLAB的读写功能修改ses文件中修改模型的命令，使得可以在MATLAB中赋值来修改有限元模型中对应的不确定参数x(如弹性模量、荷载)，完成有限元模型可被MATLAB调用并被修改的功能。

需要说明的是：第一次执行该步骤时，首先取不确定参数区间上/下界及区间中点修改结构有限元模型，第二次执行及之后将按照步骤(9)中增加的积分点来修改不确定参数。

步骤(4)生成新的有限元模型。通过步骤(3)生成新的结构有限元模型(db文件)，新的有限元模型除了上述修改的部分改变，其它都保持原来属性。

步骤(5)利用MATLAB调用有限元求解器MSC.Nasran读取步骤(4)中db文件，计算步骤(3)中结构的响应，形成xdb文件。

步骤(6)利用MATLAB调用MSC.Patran读取步骤(5)中的xdb文件，对计算结果进行后处理生成结果报告(rpt文件)，以MATLAB读取相关的响应量。

步骤(7)删除计算的过程文件，因为结算过程中产生的过程文件对下一次计算没有意义，甚至可能会和下次循环中的文件难以区分，而影响下次循环中对相关文件的调用、读取或修改，同时为了节约计算机储存空间同时不影响后续计算，删除与后续计算无关的过程文件，如步骤(5)中的xdb文件和步骤(6)中的rpt文件。

步骤(8)以结构的响应为被积函数，在不确定参数x的区间范围内进行自适应辛普森积分，所谓“自适应”即积分点的分布和被积函数在可行域内的变化趋势相关，如果在某区间内被积函数变化平缓，则积分点分布较为稀疏；如果某段区间被积函数变化剧烈，则该区间内积分点则较密，如图2所示，以一维问题来作说明，首先在积分域[a,b]的两端和中点c布置3个积分点，对函数f(·)进行辛普森积分，如下式所示：

为了验证3个积分点计算结果是否满足精度，分别在[a,c]和[c,b]域内中点增加积分点d和e，分别在这两个区间采用上式积分并求和为Q_s，如果Q_s和Q的相对误差δ满足容许误差ε，即3个积分点计算结果满足计算精度；如果Q_s和Q的相对误差δ不满足容许误差ε，则在在[a,c]和[c,b]域内分别按照上述方法验证各自区间内是否满足计算精度，如此循环，直至满足精度为止。这样就可以保证被积函数的最大/小值两侧的积分点之间的函数为凸函数，以便进一步进行牛顿法迭代得到最优解，如图3所示。

步骤(9)判断步骤(8)中积分的结果是否满足精度，判断是否满足精度的方法是判断增加积分点前后计算得到的积分值的相对误差δ是否满足容许误差ε要求。如果积分结果精度满足容许误差ε进入下一步；如果不满足计算精度则返回到步骤(3)进行新一轮的循环。其中，容许误差ε一般设置为10^-3，而相对误差δ是迭代过程中前后两次计算的积分值只差与前一次积分值的比值，即为：

其中，Q和Q_s分别是迭代前后的积分值。

步骤(10)选取现有的所有积分点集合中所对应目标函数的最大/小值及对应的积分点，如图3所示，图3(a)为全局的积分点分布情况，选取这些积分点中目标函数的最大/小值的积分点，如图3(b)所示。

步骤(11)以(11)选取的积分点作为起点，由于该点到最优解处之间为凸函数甚至是单调函数，故可以采用采用收敛速度较快的牛顿法进行局部寻优，得到的最优解即为结构响应的上/下界。

步骤(12)结束。

本发明与现有技术相比的有点在于：

本发明以优化思路为切入点，提出了基于超体积迭代策略全局优化算法的结构区间响应的不确定性传播分析方法，拓展了解决结构区间响应不确定性传播问题的思路。基于该思路可以很好地解决非线性较强系统在不确定度较大时的区间响应传播问题，可以较为精确地确定区间响应的上下界，相较于区间参数摄动法精度更高，相较于蒙特卡洛法计算效率更高，相较于遗传算法等优化算法具有明确的收敛准则。

附图说明

图1是本发明一种基于超体积迭代策略的含区间参数结构响应区间的分析方法流程图；

图2是一维自适应辛普森积分示意图；

图3是一维问题自适应积分点配置结果，其中，图3(a)为全局的积分点分布情况，图3(b)为选取这些积分点中目标函数的最大/小值的积分点；

图4是桁架结构及相关物理参数。

具体实施方式

下面结合附图以及具体实施方式进一步说明本发明。

如图3所示，本发明提出了一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，包括以下步骤：

步骤(1)基于MSC.Patran构建结构的有限元模型，包括建立几何模型，赋予材料属性，施加荷载和边界条件。

步骤(2)确定由步骤(1)构建的模型中的不确定参数x及其不确定区间。由于制作工艺、测量误差等因素导致结构中的材料和荷载具有一定的不确定性，结构中材料的弹性模量作用在结构上的荷载这里E和P分别表示弹性模量和荷载的下界，和分别表示弹性模量和荷载的上界，然后在MSC.Patran中修改相应不确定参数的数值，并且该修改过程需要用到MSC.Patran的宏录制功能最后形成记录该修改过程的ses文件。其中区间不确定性参数向量x可以表示为：

其中，x^c＝(E^c,P^c,…)为区间中点，x^r＝(E^r,P^r,…)为区间半径，e∈Ξⁿ，Ξⁿ定义为所有元素包含在[-1,1]内的n维向量集合，符号“×”定义为两个向量各对应元素乘积仍为n维向量。

步骤(3)修改结构有限元模型中对应的不确定参数。在步骤(2)中形成的ses文件中找到修改模型中不确定参数的关键性命令，利用MATLAB的读写功能修改ses文件中修改模型的命令，使得可以在MATLAB中赋值来修改有限元模型中对应的不确定参数x(如弹性模量、荷载)，完成有限元模型可被MATLAB调用并被修改的功能。

需要说明的是：第一次执行该步骤时，首先取不确定参数区间上/下界及区间中点修改结构有限元模型，第二次执行及之后将按照步骤(9)中增加的积分点来修改不确定参数。

步骤(4)生成新的有限元模型。通过步骤(3)生成新的结构有限元模型(db文件)，新的有限元模型除了上述修改的部分改变，其它都保持原来属性。

步骤(5)利用MATLAB调用有限元求解器MSC.Nasran读取步骤(4)中db文件，计算步骤(3)中结构的响应，形成xdb文件。

步骤(6)利用MATLAB调用MSC.Patran读取步骤(5)中的xdb文件，对计算结果进行后处理生成结果报告(rpt文件)，以MATLAB读取相关的响应量。

步骤(7)删除计算的过程文件，因为结算过程中产生的过程文件对下一次计算没有意义，甚至可能会和下次循环中的文件难以区分，而影响下次循环中对相关文件的调用、读取或修改，同时为了节约计算机储存空间同时不影响后续计算，删除与后续计算无关的过程文件，如步骤(5)中的xdb文件和步骤(6)中的rpt文件。

步骤(8)以结构的响应为被积函数，在不确定参数x的区间范围内进行自适应辛普森积分，所谓“自适应”即积分点的分布和被积函数在可行域内的变化趋势相关，如果在某区间内被积函数变化平缓，则积分点分布较为稀疏；如果某段区间被积函数变化剧烈，则该区间内积分点则较密，如图2所示，以一维问题来作说明，首先在积分域[a,b]的两端和中点c布置3个积分点，对函数f(·)进行辛普森积分，如下式所示：

为了验证3个积分点计算结果是否满足精度，分别在[a,c]和[c,b]域内中点增加积分点d和e，分别在这两个区间采用上式积分并求和为Q_s，如果Q_s和Q的相对误差δ满足容许误差ε，即3个积分点计算结果满足计算精度；如果Q_s和Q的相对误差δ不满足容许误差ε，则在在[a,c]和[c,b]域内分别按照上述方法验证各自区间内是否满足计算精度，如此循环，直至满足精度为止。这样就可以保证被积函数的最大/小值两侧的积分点之间的函数为凸函数，以便进一步进行牛顿法迭代得到最优解，如图3所示。

步骤(9)判断步骤(8)中积分的结果是否满足精度，判断是否满足精度的方法是判断增加积分点前后计算得到的积分值的相对误差δ是否满足容许误差ε要求。如果积分结果精度满足容许误差ε进入下一步；如果不满足计算精度则返回到步骤(3)进行新一轮的循环。其中，容许误差ε一般设置为10^-3，而相对误差δ是迭代过程中前后两次计算的积分值只差与前一次积分值的比值，即为：

其中，Q和Q_s分别是迭代前后的积分值。

步骤(10)选取现有的所有积分点集合中所对应目标函数的最大/小值及对应的积分点，如图3所示，图3(a)为全局的积分点分布情况，选取这些积分点中目标函数的最大/小值的积分点，如图3(b)所示。

步骤(11)以步骤(11)选取的积分点作为起点，由于该点到最优解处之间为凸函数甚至是单调函数，故可以采用采用收敛速度较快的牛顿法进行局部寻优，得到的最优解即为结构响应的上/下界。

步骤(12)结束。

实施例：

为了充分地了解本发明的特点及其对工程实际的适用性，本发明建立如图所示的包含弹性模量和荷载不确定性的桁架结构，并对结构区间位移传播进行分析。该桁架结构及相关物理参数如图4所示，其中，杆①、⑤和⑥的弹性模量的不确定度为恒荷载的不确定度为活荷载的不确定度为本发明的方法和蒙特卡洛法、区间参数摄动法的结果作对比，计算结果如表1所示。通过对比发现三者的计算结果相差无几，并且从表3看出，区间参数摄动法具有一定优势。但是这是针对较小的不确定度的工况，系统的非线性较弱。当结构中不确定度较大时，就会使系统的非线性增强，区间参数摄动法对小参数的强烈依赖性就体现出来了。当杆①、⑤和⑥的弹性模量的不确定度为恒荷载的不确定度为活荷载的不确定度为时，计算的结果如表2所示。本发明方法较蒙特卡洛法调用计算有限元计算次数更少，计算量较少具有很大优势。

表1小不确定度工况下桁架3、4节点位移(单位：×10^-5m)

表2大不确定度工况下桁架3、4节点位移(单位：×10^-5m)

表3各方法调用有限元程序次数对比

本发明未详细阐述部分属于本领域技术人员的公知技术。

以上仅是本发明的具体步骤，对本发明的保护范围不构成任何限制；其可扩展应用于非线性系统区间相应的不确定性传播分析中，凡采用等同变换或者等效替换而形成的技术方案，均落在本发明权利保护范围之内。

Claims (6)

1.一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，其特征在于，该方法实现步骤如下：

步骤(1)基于MSC.Patran构建结构的有限元模型，包括建立几何模型，赋予材料属性，施加荷载和边界条件；

步骤(2)确定由步骤(1)构建的模型中的不确定参数x及其不确定区间，然后在MSC.Patran中修改相应不确定参数的数值，并且该修改过程需要用到MSC.Patran的宏录制功能最后形成记录该修改过程的ses文件；

步骤(3)在步骤(2)中形成的ses文件中找到修改模型中不确定参数的关键性命令，利用MATLAB的读写功能修改ses文件中修改模型的命令，使得可以在MATLAB中赋值来修改有限元模型中对应的不确定参数x，完成有限元模型可被MATLAB调用并被修改的功能；

步骤(4)通过步骤(3)生成新的结构有限元模型，存为db文件，新的有限元模型除了上述修改的部分改变，其它都保持原来属性；

步骤(5)利用MATLAB调用有限元求解器MSC.Nastran读取步骤(4)中db文件，计算步骤(3)中结构的响应，形成xdb文件；

步骤(6)利用MATLAB调用MSC.Patran读取步骤(5)中的xdb文件，对计算结果进行后处理生成结果报告，存为rpt文件，以MATLAB读取相关的响应量；

步骤(7)删除计算的过程文件，因为结算过程中产生的过程文件对下一次计算没有意义，甚至可能会和下次循环中的文件难以区分，而影响下次循环中对相关文件的调用、读取或修改，同时为了节约计算机储存空间且不影响后续计算，删除与后续计算无关的过程文件，如步骤(5)中的xdb文件和步骤(6)中的rpt文件；

步骤(8)以结构的响应为被积函数，在不确定参数x的区间范围内进行自适应辛普森积分，所谓“自适应”即积分点的分布和被积函数在可行域内的变化趋势相关，如果在某区间内被积函数变化平缓，则积分点分布较为稀疏；如果某段区间被积函数变化剧烈，则该区间内积分点则较密，这样就可以保证被积函数的最大/小值两侧的积分点之间的函数为凸函数，以便进一步进行牛顿法迭代得到最优解；

步骤(9)判断步骤(8)中积分的结果是否满足精度，判断是否满足精度的方法是判断增加积分点前后计算得到的积分值的相对误差δ是否满足容许误差ε要求；如果积分结果精度满足容差进入下一步，如果不满足计算精度则返回到步骤(3)进行新一轮的循环；

步骤(10)选取现有的所有积分点集合中所对应目标函数的最大/小值及对应的积分点，选取这些积分点中目标函数的最大/小值的积分点；

步骤(11)以步骤(11)选取的积分点作为起点，由于该点到最优解处之间为凸函数甚至是单调函数，故可以采用采用收敛速度较快的牛顿法进行局部寻优，得到的最优解即为结构响应的上/下界；

步骤(12)结束。

2.根据权利要求1所述的一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，其特征在于：在整个实施流程中以MATLAB编写代码为主程序，首先基于MSC.Patran构建结构的有限元模型，并录制修改不确定参数的宏文件，并用MATLAB修改模型文件中相应的参数，然后调用MSC.Nastran计算结构有限元模型，以MATLAB调用MSC.Patran生成计算结果报告，运用MATLAB读取结果报告中响应量并进行自适应辛普森积分，循环执行直至收敛。

3.根据权利要求1所述的一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，其特征在于：如步骤(2)中区间不确定性参数向量x可以表示为：

$\begin{array}{c}x=\[\underset{\‾}{x},\stackrel{\‾}{x}\]=\[x^c-x^r,x^c+x^r\]\\ =x^c+x^r\[-1,1\]\\ =x^c+x^r\×e\end{array}$

其中，x^c＝(E^c,P^c,…)为区间中点，x^r＝(E^r,P^r,…)为区间半径，e∈Ξⁿ，Ξⁿ定义为所有元素包含在[-1,1]内的n维向量集合，符号“×”定义为两个向量各对应元素乘积仍为n维向量。

4.根据权利要求1所述的一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，其特征在于：步骤(8)中对响应量积分时采用的是自适应辛普森积分，首先在积分域[a,b]的两端和中点c布置3个积分点，对函数f(·)进行辛普森积分，如下式所示：

$Q=\frac{b-a}{6}\left(f\right(a)+4f(\frac{a+b}{2})+f(b\left)\right)$

为了验证3个积分点计算结果是否满足精度，分别在[a,c]和[c,b]域内中点处增加积分点d和e，分别在这两个区间采用上式积分，并将两区间积分值求和为Q_s，如果Q_s和Q的相对误差δ满足容许误差，即3个积分点计算结果满足计算精度；如果Q_s和Q的相对误差δ不满足容许误差，则在在[a,c]和[c,b]域内分别按照上述方法验证各自区间内是否满足计算精度，如此循环，直至满足精度为止。

5.根据权利要求4所述的一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，其特征在于：在步骤(3)中首先计算不确定参数取区间上/下界及区间中点处时的结构位移响应，如果步骤(9)中需要返回步骤(3)中时，则增加积分点，修改结构有限元模型的不确定参数。

6.根据权利要求1所述的一种基于超体积迭代策略的含区间参数结构响应区间的分析方法，其特征在于：为了保证计算的精确性，同时也为了使计算可以较快地收敛，步骤(9)中的容许误差ε一般设置为10^-3，而相对误差δ是指迭代过程中前后两次计算的积分值只差与前一次积分值的比值，即：

$\mathrm{\δ}=\frac{|Q-Q_s|}{Q}\×100\%$

其中，Q和Q_s分别是迭代前后的积分值。