CN112001004B - 一种解析中高频振动结构能量密度场的nurbs等几何分析方法 - Google Patents

Abstract

一种解析中高频振动结构能量密度场的NURBS等几何分析方法，先建立设计域的等几何分析模型，设计域利用NURBS曲线、NURBS曲面以及NURBS体进行描述，然后构建稳态下能量密度控制方程，再构建单元能量密度控制方程矩阵形式，然后计算基于NURBS基函数的能量密度控制方程，最后进行适应性处理，将在等几何分析框架下改造得到的能量密度控制方程代入结构的动力学分析计算中，通过求解计算得到高频激励下结构中任意位置处的能量密度响应；本发明在等几何分析的框架下，从稳态下的能量密度控制方程出发，以经过时间和空间平均的能量密度作为研究对象，采用能量有限元方法对机械结构进行动力学分析，不仅使得分析时的网格划分摆脱了振动波长的桎梏，还大大提高了精度。

Description

一种解析中高频振动结构能量密度场的NURBS等几何分析 方法

技术领域

本发明属于结构的动力学分析技术领域，具体涉及一种解析中高频振动结构能量密度场的NURBS等几何分析方法。

技术背景

复杂工况下尤其是高频激励下的结构动力学分析一直以来都是工程中备受关注的问题，以经典有限元为代表的传统方法在解决高频激励下复杂结构的动力学分析问题时的表现往往不尽人意；以航天器为例，由于航天器在飞行中受到冲击、脉动、噪声等复杂激励作用，其整流罩外气动噪声引起的振动频率可达10000Hz以上，对于这样的高频下结构动力学分析问题，结构的模态密集，且相对于整个系统的尺度，结构以短波振动为主，采用传统动力学分析方法进行分析时，为保证分析的准确性必须划分足够细小的网格，进而导致运算量增加；此外有限元分析模型只是几何模型的近似而非精确表示，高频时系统的动态响应对结构参数高度敏感，利用经典有限元方法建模时存在的细节忽略会对分析结果的准确性造成巨大影响；其几何模型与分析模型的分裂导致分析计算时在划分网格近似表达的步骤中耗费大量时间，尤其是面对类似于航天器的复杂结构时，有限元方法为了满足仿真精度必须划分非常细密的网格，导致计算的成本过大，造成分析精度与求解效率之间的矛盾，且容易产生较大误差，因此以经典有限元方法为代表的传统方法难以有效地进行高频下复杂结构的动力学分析。

发明内容

为了克服上述现有技术的缺点，本发明的目的在于提供一种解析中高频振动结构能量密度场的NURBS等几何分析方法，不仅使得分析时的网格划分摆脱了振动波长的桎梏，还大大提高了精度，为复杂工程结构的中高频响应提供计算效率更高、分析结果更加准确的动力学分析方法。

为了达到上述目标，本发明采取的技术方案为：

一种解析中高频振动结构能量密度场的NURBS等几何分析方法，包括以下步骤：

1)建立设计域的等几何分析模型：

设计域利用NURBS曲线、NURBS曲面以及NURBS体进行描述，二维结构采用NURBS曲面描述设计域形状，定义局部空间坐标系为

参数空间坐标系为ξ＝(ξ,η)，物理空间坐标系为x＝(x,y)；在参数空间中分别定义ξ和η两方向上的节点向量Ξ¹＝{ξ₁,...,ξ_n+p+1}与Ξ²＝{η₁,...,η_m+q+1}，相应方向上的NURBS基函数阶数分别为p和q，由此将设计域离散为若干NURBS单元组成的物理网格；定义P_i,j为控制点集合，参数空间中ξ和η两方向上控制点数分别为n和m，控制点总数为n_np，n_np＝m×n，由此形成了控制点组成的控制网格；

用在ξ方向p次、η方向q次的NURBS曲面描述设计域，其表达式为：

二次NURBS基函数

的表达式为：

其中，

分别表示p阶B样条基函数和q阶B样条基函数的集合，w_i,j表示权因子；

2)构建稳态下能量密度控制方程：

定义Γ表示设计域Ω的边界，从设计域中弹性波的能量流平衡原理出发得到能量密度控制微分方程：

其中，η表示阻尼系数，ω为激励频率，ε表示经过时间与空间平均的能量密度，π_in表示输入功率，c_g表示弹性波的群速度，表达式如下：

上式中，Υ表示弹性模量，h表示设计域厚度，λ表示泊松比，ρ表示材料密度；

3)构建单元能量密度控制方程矩阵形式：

以NURBS基函数作为形函数对单元能量密度ε近似表示：

其中，

表示单元控制点能量密度，n_en为每个单元上的NURBS基函数个数，由此得到式(3)的弱变分形式：

定义一个时间与空间平均的能量集中度I和能量流q分别为：

其中，

表示垂直于边界Γ的单位外法向量；

结合式(7)-式(8)，对式(6)采用一种标准的Gauss积分方法，得到单元能量密度控制方程矩阵形式：

[K^e]{E^e}＝{F^e}-{Q^e} (9)

其中，[K^e]表示单元能量矩阵，{E^e}表示单元节点能量密度向量，{F^e}表示单元输入功率向量，{Q^e}表示单元边界上的能量流向量；

4)计算基于NURBS基函数的能量密度控制方程：

4.1)NURBS基函数：

定义在参数空间上的NURBS基函数表达式为：

其中，

表示p阶B样条基函数的集合，w_i表示权因子，W(ξ)表示权函数；

NURBS基函数的一阶导通过求导法则计算得到：

4.2)通过两次映射，将物理坐标系下的NURBS单元映射到局部坐标系下的等参单元：

局部空间与参数空间之间的映射关系：

其中，(ξ_i+1,ξ_i)和(η_j+1,η_j)分别表示参数空间中NURBS单元在ξ和η方向的坐标范围；

参数空间与物理空间之间的映射关系：

dxdy＝|J_ξ|dξdη (17)

其中，

计算方法同公式(12)；

在等几何分析的框架下，几何矩阵B表达式：

其中，R_a,x表示NURBS基函数的空间导数，其计算公式如下：

上式中，

计算方法同公式(12)；

4.3)计算和组装基于NURBS基函数的单元矩阵：

在局部空间中选取高斯积分点

n_gp表示高斯积分点个数，

表示高斯积分点对应的权因子，用

表示其集合，通过高斯积分点积分得到NURBS单元上的单元能量矩阵表达式：

单元输入功率向量表达式为：

对计算得到的单元能量矩阵与单元输入功率向量采用与经典有限元相同的直接组装方法组装成全局形式[K]、{F}；

结构变化处产生不连续边界，能量密度在不连续边界上不连续；在不连续边界处新增控制点，假设两相邻单元e和单元e+1的相交边界不连续，

和

和

分别表示边界上的两单元同位控制点处能量流，在单元边界上的能量流向量用单元边界矩阵[J^{ee +1}]和单元控制点能量密度向量{E^e}表示：

按照式(24)-式(25)计算整个设计域的能量流，得到能量流向量的全局形式表达式：

{Q}＝[J]{E} (26)

其中，[J]为单元边界矩阵的全局形式，τ_ee+1和r_ee分别表示能量传递系数和反射系数；

按照上述步骤计算并组装，得到基于NURBS基函数的能量密度控制方程矩阵形式：

([K]+[J]){E}＝{F} (27)

5)适应性处理：

将在等几何分析框架下改造得到的能量密度控制方程代入结构的动力学分析计算中，通过求解计算得到高频激励下结构中任意位置处的能量密度响应。

本发明的有益效果为：

由于本发明方法是在等几何分析的框架下进行，将NURBS基函数作为形函数，可以直接利用几何信息进行分析计算，为复杂几何体分析提供更加精确的模型，无需繁琐的网格操作，这使得本发明方法的求解精度、稳定性和效率远超传统的有限元分析；同时NURBS基函数本身具有的非负性、光滑性以及在单元界面上可保持高阶连续性等特点是传统高阶有限元基函数不具备的独特优势；

本发明从稳态下的能量密度控制方程出发，以经过时间和空间平均的能量密度作为研究变量，具有能够得到结构内部任意位置响应信息的优势，可以在更宽广的频率上进行动力学分析；本发明方法为复杂机械系统的动力学分析提供了一种高效有力的工具，同时在航空航天等高端设备的振动分析控制等领域也拥有广泛的应用前景。

附图说明

图1为本发明实施例的边界条件示意图。

图2为本发明实施例中等几何分析模型示意图。

图3为本发明实施例中离散网格与NURBS基函数示意图。

图4为本发明实施例中结构能量平衡示意图。

图5为本发明实施例中单元不连续边界示意图。

图6为本发明实施例中结构能量密度响应分布示意图。

具体实施方式

下面结合附图和实施例对本发明作进一步说明，实施例设计域为四分之一圆环薄板结构，薄板尺寸参数为外圆半径R₁＝3m，内圆半径R₂＝1m，板厚h＝0.01m，定义薄板的内圆边缘固定，外圆边缘受均匀分布的单位简谐载荷，如图1所示。

一种解析中高频振动结构能量密度场的NURBS等几何分析方法，包括以下步骤：

1)建立设计域的等几何分析模型：

设计域利用NURBS曲线、NURBS曲面以及NURBS体进行描述，本实施例中设计域二维结构由一个二次NURBS曲面描述，如图2所示，定义局部空间坐标系为

参数空间坐标系为ξ＝(ξ,η)，物理空间坐标系为x＝(x,y)；在参数空间中分别定义ξ和η两方向上的节点向量Ξ¹＝{0,0,1/3,2/3,1,1}与Ξ²＝{0,0,1/3,2/3,1,1}，相应方向上的NURBS基函数阶数分别为p＝1和q＝1，由此将设计域离散为若干NURBS单元组成的物理网格；P_i,j为控制点集合，参数空间中ξ和η两方向上控制点数分别为n＝4和m＝4，控制点总数为n_np＝16，由此形成了控制点组成的控制网格；

如图3所示，用在ξ方向p次、η方向q次的NURBS曲面来描述设计域，其表达式为：

二次NURBS基函数

的表达式为：

其中，

分别表示p阶B样条基函数和q阶B样条基函数的集合，w_i,j表示权因子；

2)构建稳态下能量密度控制方程：

如图4所示，定义Γ表示设计域Ω的边界，从薄板中弹性波的能量平衡公式出发得到能量密度控制微分方程：

其中，η表示阻尼系数，取0.08，ω为激励频率，取1000Hz。ε表示经过时间与空间平均的能量密度，π_in表示输入功率，为1W，c_g表示弹性波的群速度，表达式如下：

上式中，Υ表示弹性模量，取Υ＝7.1×10¹⁰N·m²，h表示设计域厚度，λ表示泊松比，取0.35，ρ表示材料密度，取2700kg/m²；

3)构建单元能量密度控制方程矩阵形式：

以NURBS基函数作为形函数对单元能量密度ε近似表示：

其中，

表示单元控制点能量密度，n_en为每个单元上的NURBS基函数个数，由此得到式(3)的弱变分形式：

定义一个时间与空间平均的能量集中度I和能量流q分别为：

其中，

表示垂直于边界Γ的单位外法向量；

结合式(7)-式(8)，对式(6)采用一种标准的Gauss积分方法，得到单元能量密度控制方程矩阵形式：

其中，[K^e]表示单元能量矩阵，{E^e}表示单元节点能量密度向量，{F^e}表示单元输入功率向量，{Q^e}表示单元边界上的能量流向量；

4)计算基于NURBS基函数的能量密度控制方程：

4.1)NURBS基函数：

定义在参数空间上的NURBS基函数表达式为：

其中，

表示p阶B样条基函数的集合，w_i表示权因子，本实施例中取权因子全为1，W(ξ)表示权函数；

NURBS基函数的一阶导通过求导法则计算得到：

4.2)通过两次映射，将物理坐标系下的NURBS单元映射到局部坐标系下的等参单元：

局部空间与参数空间之间的映射关系：

其中，(ξ_i+1,ξ_i)和(η_j+1,η_j)分别表示参数空间中NURBS单元在ξ和η方向的坐标范围；

参数空间与物理空间之间的映射关系：

dxdy＝|J_ξ|dξdη (17)

其中

计算方法同公式(12)；

在等几何分析的框架下，几何矩阵B表达式：

其中，R_a,x表示非零NURBS基函数的空间导数，其计算公式如下：

上式中，

计算方法同公式(12)；

4.3)计算和组装基于NURBS基函数的矩阵：

在局部空间中选取高斯积分点

n_gp表示高斯积分点个数，本实施例中取4，

表示高斯积分点对应的权因子，用

表示其集合，本实施例中取权因子全为1，通过高斯积分点积分得到NURBS单元上的单元能量矩阵表达式：

单元输入功率向量表达式为：

对计算得到的单元能量矩阵与单元输入功率向量采用与经典有限元相同的直接组装方法组装成全局形式[K]、{F}；

结构变化处产生不连续边界，能量密度在不连续边界上不连续；在不连续边界处新增控制点，如图5所示，两相邻单元1和单元2的相交边界不连续，Q₂和Q₅、Q₃和Q₈分别表示边界上的两单元同位节点处能量流，在单元边界上的能量流向量用单元边界矩阵[J¹²]和单元节点能量密度向量{E^e}表示：

按照式(24)-式(25)计算整个设计域的能量流，得到能量流向量的全局形式表达式：

{Q}＝[J]{E} (26)

其中，[J]为单元边界矩阵的全局形式，τ₁₂和r₁₁分别表示能量传递系数和反射系数，具体计算见已有学术研究，此处不再赘述；

按照上述步骤计算并组装，得到基于NURBS基函数的能量密度控制方程矩阵形式：

([K]+[J]){E}＝{F} (27)

8)适应性处理：

将在等几何分析框架下改造得到的能量密度控制方程代入结构的动力学分析计算中，通过求解计算可以得到高频激励下结构中任意位置处的能量密度响应，如图6所示。

为适应不同分析需求，使用时并不局限于所述的材料相同的薄板结构，设计者可以将本发明方法运用到多相结构的动力学分析运算中，且对于更加复杂的结构模型同样能够适应并得到准确的分析结果。

Claims (1)

1.一种解析中高频振动结构能量密度场的NURBS等几何分析方法，其特征在于，包括以下步骤：

1)建立设计域的等几何分析模型：

设计域利用NURBS曲线、NURBS曲面以及NURBS体进行描述，二维结构采用NURBS曲面描述设计域形状，定义局部空间坐标系为

参数空间坐标系为ξ＝(ξ,η)，物理空间坐标系为x＝(x,y)；在参数空间中分别定义ξ和η两方向上的节点向量Ξ¹＝{ξ₁,...,ξ_n+p+1}与Ξ²＝{η₁,...,η_m+q+1}，相应方向上的NURBS基函数阶数分别为p和q，由此将设计域离散为若干NURBS单元组成的物理网格；定义P_i,j为控制点集合，参数空间中ξ和η两方向上控制点数分别为n和m，控制点总数为n_np，n_np＝m×n，由此形成了控制点组成的控制网格；

用在ξ方向p次、η方向q次的NURBS曲面描述设计域，其表达式为：

二次NURBS基函数

的表达式为：

其中，

分别表示p阶B样条基函数和q阶B样条基函数的集合，w_i,j表示权因子；

2)构建稳态下能量密度控制方程：

定义Γ表示设计域Ω的边界，从设计域中弹性波的能量流平衡原理出发得到能量密度控制微分方程：

其中，η表示阻尼系数，ω为激励频率，ε表示经过时间与空间平均的能量密度，π_in表示输入功率，c_g表示弹性波的群速度，表达式如下：

上式中，Υ表示弹性模量，h表示设计域厚度，λ表示泊松比，ρ表示材料密度；

3)构建单元能量密度控制方程矩阵形式：

以NURBS基函数作为形函数对单元能量密度ε近似表示：

其中，

表示单元控制点能量密度，n_en为每个单元上的NURBS基函数个数，由此得到式(3)的弱变分形式：

定义一个时间与空间平均的能量集中度I和能量流q分别为：

其中，

表示垂直于边界Γ的单位外法向量；

结合式(7)-式(8)，对式(6)采用一种标准的Gauss积分方法，得到单元能量密度控制方程矩阵形式：

[K^e]{E^e}＝{F^e}-{Q^e} (9)

其中，[K^e]表示单元能量矩阵，{E^e}表示单元节点能量密度向量，{F^e}表示单元输入功率向量，{Q^e}表示单元边界上的能量流向量；

4)计算基于NURBS基函数的能量密度控制方程：

4.1)NURBS基函数：

定义在参数空间上的NURBS基函数表达式为：

其中，

表示p阶B样条基函数的集合，w_i表示权因子，W(ξ)表示权函数；

NURBS基函数的一阶导通过求导法则计算得到：

4.2)通过两次映射，将物理坐标系下的NURBS单元映射到局部坐标系下的等参单元：

局部空间与参数空间之间的映射关系：

其中，(ξ_i+1,ξ_i)和(η_j+1,η_j)分别表示参数空间中NURBS单元在ξ和η方向的坐标范围；

参数空间与物理空间之间的映射关系：

dxdy＝|J_ξ|dξdη (17)

其中，

计算方法同公式(12)；

在等几何分析的框架下，几何矩阵B表达式：

其中，R_a,x表示NURBS基函数的空间导数，其计算公式如下：

上式中，

计算方法同公式(12)；

4.3)计算和组装基于NURBS基函数的单元矩阵：

在局部空间中选取高斯积分点

n_gp表示高斯积分点个数，

表示高斯积分点对应的权因子，用

表示其集合，通过高斯积分点积分得到NURBS单元上的单元能量矩阵表达式：

单元输入功率向量表达式为：

对计算得到的单元能量矩阵与单元输入功率向量采用与经典有限元相同的直接组装方法组装成全局形式[K]、{F}；

结构变化处产生不连续边界，能量密度在不连续边界上不连续；在不连续边界处新增控制点，假设两相邻单元e和单元e+1的相交边界不连续，

和

和

分别表示边界上的两单元同位控制点处能量流，在单元边界上的能量流向量用单元边界矩阵[J^ee+1]和单元控制点能量密度向量{E^e}表示：

按照式(24)-式(25)计算整个设计域的能量流，得到能量流向量的全局形式表达式：

{Q}＝[J]{E} (26)

其中，[J]为单元边界矩阵的全局形式，τ_ee+1和r_ee分别表示能量传递系数和反射系数；

按照上述步骤计算并组装，得到基于NURBS基函数的能量密度控制方程矩阵形式：

([K]+[J]){E}＝{F} (27)

5)适应性处理：

将在等几何分析框架下改造得到的能量密度控制方程代入结构的动力学分析计算中，通过求解计算得到高频激励下结构中任意位置处的能量密度响应。

Images