1.本发明公开一种获取热环境双向功能梯度梁固有频率的方法与算法，涉及热环境双向功能梯度梁固有频率的获取方法技术领域。
背景技术：
2.现有技术中，求解变系数偏微分方程问题的常用方法包括同伦法，里兹法，谐波微分求积法,里兹法结合newmark法等方法。这些方法多为数值解法，计算量大。结合传递函数法和摄动法是一种有效的计算方法，可以将解变系数偏微分方程问题转化为常系数偏微分方程问题；但是当变系数过多时或阶数过多时方程系数复杂程度增加，其求解精度有所下降，甚至不能满足工程需求。这种问题在解析求解热环境下双向梯度功能梯度材料梁振动问题尤为显著。
3.本
技术实现要素：
4.本发明目的在于，提供一种获取热环境双向功能梯度梁固有频率的方法与算法，解决求解热环境下双向梯度功能梯度材料梁振动问题精度不佳的问题。
5.为实现上述技术目的，达到上述技术效果，发明是通过以下技术方案实现：
6.一种获取热环境双向功能梯度梁固有频率的方法，包括如下步骤：
7.步骤a、定义无量纲小参数ε，令
[0008][0009]
步骤b、利用摄动法将所述状态方程展开，得到一系列常系数偏微分方程；
[0010]
步骤c、求解零系数偏微分方程；
[0011]
步骤d骤b，计算方程系数表达式，并基于该表达式和零阶摄动结果，迭代计算摄动量，并代回步骤b的频率展开式中，获得摄动后的固有频率。
[0012]
进一步的，对ε的具体应用为：
[0013]
步骤a、将振动控制方程变形为带有ε的状态方程
[0014]
材料参数为
[0015][0016]
引入归一化参数，
[0017][0018]
以及：
[0019][0020]
其中，系数为
[0021][0022]
根据传递函数法定义状态向量
[0023][0024]
利用状态向量将控制方程写为状态方程形式，为
[0025][0026]
其中，
[0027][0028]
其中元素为
[0029][0030]
其对应边界条件为
[0031]
m(ω)η(0，ω) n(ω)η(1，ω)＝0；
[0032]
其中，m(ω)和n(ω)分别表示梁左端和右端的边界条件选择矩阵。
[0033]
步骤b、摄动展开
[0034]
基于一阶摄动法，利用ε分解相关参数，将式转化为常系数微分方程。
[0035]
ω＝ω0 εω1；
[0036]
[0037][0038][0039][0040]
式中，ωj(j＝0，1，2
…
)为与固有频率ω有关的参数；
[0041]
略去ε的高阶项(即ε2)，得到摄动后的两组控制方程组，即零阶摄动控制方程，和一阶摄动控制方程，为
[0042][0043]
m0(ω0)η0(0，ω0) n0(ω0)η0(1，ω0)＝0；
[0044]
和
[0045][0046]
步骤c、求解零阶摄动控制方程
[0047]
首先计算系数φ0，根据摄动法，
[0048]
φ0(ω0)＝[φ
mn
|
ε＝0
]6×6；
[0049]
其中，φ
mn
为矩阵的元素。
[0050]
进一步的，当ε＝0时，梯度材料退化为均匀材料，材料参数为左下角材料参数，为
[0051][0052]
带入计算获得：
[0053][0054]
得到：
[0055][0056]
带入边界条件，得到
[0057][0058]
当特征向量η0(0，ω0)存在非零解时，需满足
[0059][0060]
由此，可计算控制方程的零阶摄动解ω0，也是均匀材料梁的固有频率，其对应的固有振型为
[0061]
[0062]
其中，v0为满足时，矩阵
[0063][0064]
的特征向量；
[0065]
所述步骤d为，根据摄动法计算φ1，迭代计算目标值；
[0066][0067]
只需迭代计算当
[0068]
o(ε2)＝min(o(ε2))；
[0069]
时，所对应的ω1，再代回ω＝ω0 εω1，即可得到ω迭代后的结果。
[0070]
本发明的另一目的在于，公开一种获取热环境双向功能梯度梁固有频率的方法的实现算法，其特征在于，所述步骤具体包括：
[0071]
步骤a：输入原始问题的变系数偏微分方程；
[0072]
l(μ(i，t))＝a(x，t)μi(x，t) b(x，t)μi(x，t) c(x，t)μi(x，t) ...＝f(x，t)；i∈(1，2，...，n)；
[0073]
步骤b：选择合适的无量纲化摄动小参数ε，并引入原始方程；
[0074]
a(x，t)＝a0(x，t) εa1(x，t) ... εa
i-1
(x，t)；
[0075]
将其代入原始方程，得到：
[0076]
l
′
(μ
′
(x，t)，ε)＝(a0(x，t) εa1(x，t)μ
′i(x，t) b(x，t)μ
′
i-1
(x，t)...)＝f(x，t)；
[0077]
步骤c：调用摄动展开函数，将带有摄动小参数的方程展开为一系列常系数偏微分方程；
[0078]
μ
′i(x，t)＝μ0(x，t) εμ1(x，t) ε2μ2(x，t) ...；
[0079]
将其代入带有摄动参数的方程l
′
(μ
′
(x，t)，ε)，逐阶收集关于ε的项，得到一系列常系数偏微分方程；
[0080]
步骤d：求解零阶摄动方程的函数；
[0081]
l0(μ0(x，t))＝f0(x，t)；
[0082]
步骤e：通过迭代计算固有频率的函数，根据求解出的零阶摄动方程，得到一阶、二阶等高阶摄动方程：
[0083]
a(x，t)μ
″i(x，t) b(x，t)μ
′i(x，t) c(x，t)μi(x，t)＝fi(x，μ
i-1
)；i∈(1，2，...，n)；
[0084]
步骤f：调用迭代函数，根据摄动法计算得到的固有频率进行迭代求解依次求解，μ1(x)，μ2(x)...直到满足收敛条件；
[0085]
步骤g：输出原始问题的更精确的固有频率，
[0086]
μ(x，ε)＝μ0(x) μ1(x) ε2μ2(x) ...。
[0087]
有益效果：
[0088]
通过摄动法，考虑梯度材料特点构造小参数，将非线性问题转化为线性问题。结合传递函数获得结果，相比其他的方法复杂度低而且保持了良好的精度。此外基于摄动法与传递函数法，针对其摄动展开误差较大的环节，提出了一种改进的迭代方法，大大简化了计算步骤，也大大提高了计算精度。在计算热环境下双向梯度功能梯度梁的振动固有频率时，简化了计算步骤，且最大精度提高了15.67％。
[0089]
当然，实施本发明的任一产品并不一定需要同时达到以上所述的所有优点。
附图说明
[0090]
图1为本发明实施例所述一种求解变系数偏微分方程问题的方法逻辑图；
[0091]
图2为本发明实施例所述一种求解变系数偏微分方程问题的方法的应用逻辑图；
[0092]
图3为本发明实施例所述一种求解变系数偏微分方程问题的方法的实现算法逻辑图；
[0093]
图4为本发明实施例所述一种2d fg梁的尺寸和作标定义图；
具体实施方式
[0094]
为了更清楚地说明本发明实施例的技术方案，下面将结合附图对实施例对本发明进行详细说明。
[0095]
本发明开发了一种改进的半解析半数值的求解方法，既简化了计算步骤，又能得到较高精度的结果。该方法能够便利、精确、半解析半数值地求解双向梯度功能梯度梁在热环境下的振动固有频率，其改进思路新颖，可扩展应用于类似的数学和力学领域。
[0096]
具体而言，一种获取热环境双向功能梯度梁固有频率的方法，包括如下步骤：步骤a、定义无量纲小参数ε，令
[0097][0098]
所述无量纲化摄动小参数通常表示为ε，用ε修改原始问题的方程，获得带有ε的状态方程；
[0099]
步骤b、利用摄动法将所述状态方程展开，得到一系列常系数偏微分方程；
[0100]
步骤c、求解零系数偏微分方程；
[0101]
步骤d计算方程系数表达式，并基于该表达式和零阶摄动结果，迭代计算摄动量，并代回步骤b的频率展开式中，获得摄动后的固有频率。
[0102]
可见本实施例基于摄动法与传递函数法，针对其摄动展开误差较大的环节，提出了一种改进的迭代方法，大大简化了计算步骤，也大大提高了计算精度。
[0103]
实施例1
[0104]
本实施例中，双向梯度功能梯度梁的坐标定义、材料分布和尺寸如图4所示，其中l和h分别代表梁的轴向长度和截面高度。有效材料参数沿截面高度和轴向梯度变化。梁处于均匀/线性/非线性温度场下，其上表面温度为tc，下表面温度为tm。如图4所示。
[0105]
基于幂律假设，假设材料属性随温度变化，则梁的有效材料参数可表示为弹性模量e(x，z，t)，密度ρ(x，z，e)和泊松比v(x，z，t)为
[0106][0107]
其中，p代表材料的有效弹性模量e(x，z，t)，密度ρ(x，z，e)，泊松比v(x，z，e)，热膨胀系数α(x，z，e)和热传导系数k(x，z，e)，p
x
，pz为在x和z方向上的幂律指数，下标m代表左下角的材料，下标c代表右上角的材料。
[0108]
假设的线性温度场t(z)
lin
表达式为：
[0109][0110]
求解一维热传导方程得到的非线性温度场t(x，z，t)
non
表达式为
[0111][0112]
根据timosheno梁和哈密顿原理，建立梁的振动控制方程组，为
[0113][0114]
其中，u(x)，w(x)，φ(x)为待求位移场，ω为固有频率。式中各项系数为
[0115][0116]
基于前述内容，所述步骤具体包括：
[0117]
步骤a、将振动控制方程变形为带有ε的状态方程
[0118]
定义无量纲小参数ε，令
[0119][0120]
材料参数为
[0121][0122]
引入归一化参数，
[0123][0124]
以及：
[0125][0126]
其中，系数为
[0127][0128]
根据传递函数法定义状态向量
[0129][0130]
利用状态向量将控制方程写为状态方程形式，为
[0131][0132]
其中，
[0133][0134]
其中元素为
[0135][0136]
其对应边界条件为
[0137]
m(ω)η(0，ω) n(ω)η(1，ω)＝0；
[0138]
其中，m(ω)和.n(ω).分别表示梁左端()和右端()的边界条件选择矩阵。
[0139]
步骤b、摄动展开
[0140]
基于一阶摄动法，利用ε分解相关参数，将式转化为常系数微分方程。
[0141]
ω＝ω0 εω1；
[0142][0143][0144][0145][0146]
式中，ωj(j＝0，1，2
…
)为与固有频率ω有关的参数；
[0147]
略去ε的高阶项(即ε2)，得到摄动后的两组控制方程组，即零阶摄动控制方程，和一阶摄动控制方程，为
[0148][0149]
m0(ω0)η0(0，ω0) n0(ω0)η0(1，ω0)＝0；
[0150]
和
[0151][0152]
步骤c、求解零阶摄动控制方程
[0153]
首先计算系数φ0，根据摄动法，
[0154]
φ0(ω0)＝[φ
mn
|
ε＝0
]6×6；
[0155]
其中，φ
mn
为矩阵的元素；
[0156]
当ε＝0时，梯度材料退化为均匀材料，材料参数为左下角材料参数，为
[0157][0158]
带入计算获得：
[0159][0160]
得到：
[0161][0162]
带入边界条件，得到
[0163][0164]
当特征向量η0(0，ω0)存在非零解时，需满足
[0165][0166]
由此，可计算控制方程的零阶摄动解ω0，也是均匀材料梁的固有频率，其对应的固有振型为
[0167][0168]
其中，v0为满足时，矩阵
[0169][0170]
的特征向量；
[0171]
所述步骤d为，根据摄动法计算φ1，迭代计算目标值；
[0172]
[0173]
只需迭代计算当
[0174]
o(ε2)＝min(o(ε2))；
[0175]
时，所对应的ω1，再代回ω＝ω0 εω1，即可得到ω迭代后的结果。
[0176]
本实施例中，对比了非线性温度场时以及gdqm计算结果(用已有成熟方法计算的结果)和文献(wei-ren chen，heng chang，vibration analysis of bi-directional functionally graded timoshenko beams using chebyshev collocation method，international journal of structural stability and dynamics，doi：10.1142/s0219455421500097)，本方法计算出的固有频率，未改进的方法计算出的固有频率，和gdqm计算出的固有频率(用已有成熟方法计算的结果)，如表1所示。
[0177]
表1改变p
x
和pz时，本发明方法，未改进的方法和gdqm计算的固有频率对比
[0178][0179][0180]
可见，本发明计算结果与gdqm方法非常接近，相较于为改进的方法，计算精度有了很大提高，具有很高实用价值。
[0181]
实施例2
[0182]
根据前述实施例内容，本实施例公开一种获取热环境双向功能梯度梁固有频率的方法的实现算法，所述步骤具体包括：
[0183]
步骤a：输入原始问题的变系数偏微分方程；
[0184]
l(μ(i，t))＝a(x，t)μi(x，t) b(x，t)μi(x，t) c(x，t)μi(x，t) ...＝f(x，t)；i∈(1，2，...，n)；
[0185]
步骤b：选择合适的无量纲化摄动小参数ε，并引入原始方程；
[0186]
a(x，t)＝a0(x，t) εa1(x，t) ... εa
i-1
(x，t)；
[0187]
将其代入原始方程，得到：
[0188]
l
′
(μ
′
(x，t)，ε)＝(a0(x，t) εa1(x，t)μ
′i(x，t) b(x，t)μ
′
i-1
(x，t)...)＝f(x，t)；
[0189]
步骤c：调用摄动展开函数，将带有摄动小参数的方程展开为一系列常系数偏微分方程；
[0190]
μ
′i(x，t)＝μ0(x，t) εμ1(x，t) ε2μ2(x，t) ...；
[0191]
将其代入带有摄动参数的方程l
′
(μ
′
(x，t)，ε)，逐阶收集关于ε的项，得到一系列常系数偏微分方程；
[0192]
步骤d：求解零阶摄动方程的函数；
[0193]
l0(μ0(x，t))＝f0(x，t)；
[0194]
步骤e：通过迭代计算固有频率的函数，根据求解出的零阶摄动方程，得到一阶、二阶等高阶摄动方程：
[0195]
a(x，t)μ
″i(x，t) b(x，t)μ
′i(x，t) c(x，t)μi(x，t)＝fi(x，μ
i-1
)；i∈(1，2，...，n)；
[0196]
步骤f：调用迭代函数，根据摄动法计算得到的固有频率进行迭代求解依次求解，μ1(x)，μ2(x)...直到满足收敛条件；
[0197]
步骤g：输出原始问题的更精确的固有频率，
[0198]
μ(x，ε)＝μ0(x) μ1(x) ε2μ2(x) ...。
[0199]
实际上，本实施例所公开的算法可以通过任意一种电子设备或者计算机实现，可以应用任意一种包含数学库的编程语言。优选的可以使用matlab。当然基于算法转化为代码语言属于智力活动范畴，本实施例不再赘述。
[0200]
具体的来说，将本发明中涉及的方法转化为可实现的算法逻辑实际上才体现了本发明的实际目的。本算法可以应用于多种情况的计算，例如变系数偏微分方程在其他领域中的应用，如流体力学、结构力学、热力学等；或者其他类型的功能梯度材料结构问题，如板、壳等。本发明由于其简单、解析法计算的特征，适用于计算机计算，并且不需要过高的浮点精度。
[0201]
本方法具有很好的扩展性，因此可以根据输入的内容，将不同的变系数偏微分方程转化为常系数偏微分方程，对于复杂的非线性问题，本方法可以作为一种有效的近似解析方法，降低问题的求解难度。为工程实践带来便利。
[0202]
以上仅是该申请的实施例部分，并非对该申请做任何形式上的限制。对以上实施例所做的任何简单的修改、等同变化及修饰，仍属于该申请技术方案保护的范围内。