CN115407680A - 一种并行化的高精度时域气动弹性控制耦合仿真方法 - Google Patents

一种并行化的高精度时域气动弹性控制耦合仿真方法 Download PDF

Info

Publication number
CN115407680A
CN115407680A CN202211103424.1A CN202211103424A CN115407680A CN 115407680 A CN115407680 A CN 115407680A CN 202211103424 A CN202211103424 A CN 202211103424A CN 115407680 A CN115407680 A CN 115407680A
Authority
CN
China
Prior art keywords
thread
fluid
generalized
fluid grid
grid
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211103424.1A
Other languages
English (en)
Inventor
戴玉婷
吴优
李永昌
杨超
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beihang University
Original Assignee
Beihang University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beihang University filed Critical Beihang University
Priority to CN202211103424.1A priority Critical patent/CN115407680A/zh
Publication of CN115407680A publication Critical patent/CN115407680A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B17/00—Systems involving the use of models or simulators of said systems
    • G05B17/02—Systems involving the use of models or simulators of said systems electric
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Automation & Control Theory (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种并行化的高精度时域气动弹性控制耦合仿真方法，属于航空气动弹性领域；首先通过结构有限元软件输出飞行器结构节点坐标和结构模态，然后利用网格划分软件生成飞行器的流体网格，通过Fluent的用户自定义函数接口编写并行化径向基函数插值方法，将飞行器结构模态插值到流体网格中心点和流体网格节点上；之后计算气动力求解广义气动弹性方程组，建立并行化的气动弹性求解器。通过控制模块与气动弹性求解器的结果耦合，结构弹性变形和控制模块操纵指令引起的流体网格边界变形采用动网格技术计算，最终实现气动弹性控制的耦合时域仿真。本发明适用范围广，操作便捷，可以精确且高效的分析气动结构控制耦合问题。

Description

一种并行化的高精度时域气动弹性控制耦合仿真方法
技术领域
本发明属于航空气动弹性领域，涉及一种并行化的高精度时域气动弹性控制耦合仿真方法。
背景技术
气动力、结构弹性和控制系统的耦合是现代飞行器设计中气动伺服弹性综合问题的核心。气动伺服弹性综合是指针对给定的飞行器性能指标设计主动控制律，构成闭环反馈，从而达到预期的指标，以提高系统性能(杨超主编.飞行器气动弹性原理[M].北京航空航天大学出版社,2011)。
传统的气动伺服弹性工程分析方法一般具有如下特征：气动力的求解采用片条理论、面元法或西奥道森方法等工程气动模型；结构弹性响应的求解采用模态叠加法；控制系统表示为传递函数形式；耦合计算和分析在频域中执行。
这一套传统工程分析方法适用于经典气动伺服弹性综合问题；然而，随着市场对飞行器功能需求的多样化，各式各样的特种飞行器应运而生。特种飞行器大多应用于各种非稳态的、非线性的复杂工作环境中，由此也引出了各种非经典的气动伺服弹性综合问题。
面对这些非经典问题，传统计算分析方法在精度上无法满足现代工程需求，也无法通过频域分析研究和预测特种飞行器在实时飞行过程中遇到的各种瞬态问题，这一问题已经逐渐成为现代飞行器设计的主要瓶颈。
现阶段适用于非经典气动伺服弹性问题的分析方法较少且现有方法不足较多。因此，工业界迫切的需要一种普适性强、高效率和高精度的时域气动伺服弹性综合分析方法。
发明内容
针对上述问题，本发明提出一种并行化的高精度时域气动弹性控制耦合仿真方法，气动力的计算采用计算流体力学方法，通过用户自定义函数接口编写的结构弹性响应求解器和控制器，可以实现气动弹性控制的耦合时域仿真，适用于分析处理各种流动复杂情况的气动力。
所述的耦合仿真方法，具体步骤如下：
步骤一、使用结构有限元软件输出飞行器结构节点坐标和结构模态到txt文件中；
步骤二、使用网格划分软件生成飞行器外表面的流体网格，以并行模式启动计算流体力学软件Fluent，导入流体网格。
步骤三、针对当前时间步，通过并行化径向基函数插值方法，将飞行器的结构模态插值，得到每个分线程中流体网格的中心点模态和流体网格节点模态。
具体为：首先，在主线程中构建插值右矩阵，并传递给每个分线程；
然后，每个分线程分别各自对应飞行器一个子域的流体网格，通过读取流体网格的中心点坐标，并构建对应的第一插值左矩阵，结合插值右矩阵，将飞行器的结构模态插值得到当前流体网格中心点的模态。
同时、在每个分线程中，读取该线程内流场计算流体网格节点坐标，并构建第二插值左矩阵，结合插值右矩阵，将结构模态插值得到各流体网格节点的模态；
步骤四、每个分线程分别计算流体网格中心点的气动力，并将流体网格中心点的气动力与流体网格中心点模态相乘，对所有分线程的乘积求和，得到广义气动力传递给主线程；
步骤五、主线程利用广义气动力求解用户自定义函数中编写的广义气动弹性运动方程组，得到每阶模态对应的广义位移和广义速度，并将其传递至各分线程；
步骤六、在每个分线程中，使用流体网格节点模态和广义位移与广义速度，计算流体网格节点的真实位移与真实速度；
步骤七、主线程通过用户自定义函数调用控制模块，将流体网格节点的真实位移与真实速度传递给控制模块，分别输出对应的操纵指令给各分线程。
步骤八、在每个分线程中，构建飞行器操纵面变形的型函数，使用操纵指令与型函数计算由操纵指令引起的流体网格节点的附加位移。
步骤九、在每个分线程中，使用Fluent软件的动网格技术计算由流体网格节点的真实位移和附加位移引起的网格变形，并更新每个分线程的网格。
步骤十、判断计算是否到达终止时间步，如果是，进入步骤十一，否则，返回步骤三；
步骤十一、通过主线程输出各个时间步的广义位移和广义速度的txt文件。
本发明与现有技术相比，具有以下优势：
(1)本发明一种并行化的高精度时域气动弹性控制耦合仿真方法，通过并行化插值技术，可以支持多核并行计算，进一步提升了求解的计算速度。
(2)本发明一种并行化的高精度时域气动弹性控制耦合仿真方法，插值在每一时间步、每个线程内同时独立的执行，使得本发明可以支持网格重构、重叠网格等各种更先进的动网格技术。也可以适应各种非线性的大幅结构变形问题。对比现有技术使用了全局编号匹配搜索技术实现并行计算，在使用时必须保持网格拓扑结构和编号顺序不变，且由于仅采用了一次非并行的插值，飞行器变形偏离原始外形过大后，插值方法失效。相比之下，本发明具有更高度的普适性和更小的使用限制。
(3)本发明一种并行化的高精度时域气动弹性控制耦合仿真方法，在使用中仅需启动一次Fluent，并且仅需要一份预处理数据txt文件。对比现有技术使用了全局编号匹配搜索技术实现并行计算，在应用中需要多次启动Fluent进行数据处理，并且用户需要准备两份额外的txt文件作为预处理数据。相比之下，本发明极大地简化了用户使用流程。
附图说明
图1是本发明一种并行化的高精度时域气动弹性控制耦合仿真方法流程图；
图2是本发明实施例在一个时间步内的计算流程图。
具体实施方式
下面结合附图和示例对本发明做进一步的阐述。
本发明公开了一种并行化的高精度时域气动弹性控制耦合仿真方法，首先通过结构有限元软件输出飞行器结构节点坐标和结构模态到txt文件中，然后利用网格划分软件生成飞行器的流体网格，通过Fluent的用户自定义函数接口编写并行化径向基函数插值方法，将飞行器结构模态插值到流体网格中心点和流体网格节点上。之后通过导入计算流体力学软件Fluent计算气动力，借助用户自定义函数编写广义气动弹性方程组，建立并行化的气动弹性求解器。通过用户自定义函数编写控制模块，并与上述气动弹性求解器的结果耦合。结构弹性变形和控制模块操纵指令引起的流体网格边界变形采用动网格技术计算。本发明适用于分析处理各种流动复杂情况的气动力，适用范围广，操作便捷，通过用户自定义函数接口编写的结构弹性响应求解器和控制器，可以精确且高效的分析气动结构控制耦合问题，实现气动弹性控制的耦合时域仿真。
如图1所示，具体步骤如下：
步骤一、使用结构有限元软件输出飞行器结构节点坐标和结构模态到txt文件中；
步骤二、使用网格划分软件生成飞行器外表面的流体网格，以并行模式启动计算流体力学软件Fluent，导入流体网格，读取txt文件；
步骤三、针对当前时间步，通过Fluent的用户自定义函数接口编写并行化径向基函数插值方法，将飞行器的结构模态插值，得到每个分线程中流体网格的中心点模态和流体网格节点模态。
流体网格中心点即飞行器表面的正方形网格的中心，流体网格节点即飞行器表面的正方形网格的四个角点。
具体为：首先，在主线程中构建插值右矩阵，并传递给每个分线程；
然后，每个分线程分别各自对应飞行器一个子域的流体网格，通过读取该线程内流场计算流体网格的中心点坐标，并构建对应的第一插值左矩阵，结合插值右矩阵，将飞行器的结构模态插值得到当前流体网格中心点的模态。
同时、在每个分线程中，读取该线程内流场计算流体网格节点坐标，并构建第二插值左矩阵，结合插值右矩阵，将结构模态插值得到各流体网格节点的模态；
步骤四、每个分线程分别使用Fluent计算流体网格中心点的气动力，并将流体网格中心点的气动力与流体网格中心点模态相乘，对所有分线程的乘积求和，得到广义气动力传递给主线程；
步骤五、主线程利用广义气动力求解用户自定义函数中编写的广义气动弹性运动方程组，得到每阶模态对应的广义位移和广义速度，并将其传递至各分线程；
步骤六、在每个分线程中，使用流体网格节点模态和广义位移与广义速度，计算流体网格节点的真实位移与真实速度；
步骤七、主线程通过用户自定义函数调用控制模块，将飞行器传感器附近的流体网格节点的真实位移与真实速度传递给控制模块，分别输出对应的操纵指令给各分线程。
步骤八、在每个分线程中，构建飞行器操纵面变形的型函数，使用操纵指令与型函数计算由操纵指令引起的流体网格节点的附加位移。
步骤九、在每个分线程中，使用Fluent软件的动网格技术计算由流体网格节点的真实位移和附加位移引起的网格变形，并更新每个分线程的网格。
步骤十、判断计算是否到达终止时间步，如果是，进入步骤十一，否则，返回步骤三；
步骤十一、通过主线程输出各个时间步的广义位移和广义速度的txt文件。
实施例：
如图2所示，虚线框表示了本实施例的基本计算单元；其中有限元软件只需在程序开始前执行一次，以获得飞行器结构相关数据。主线程和所有分线程是主体，在每个时间步都要执行一次。
主线程在并行计算中只负责计算的全局调度，向用户界面或文件输出信息，接受来自用户界面或文件的输入，并且可以处理简单的计算任务。
分线程在并行计算中负责主要的计算和数据存储任务，接受主线程的指令，接受来自主线程的数据，并且可以将计算结果反馈到主线程。本实施例仅展示两个线程作示意，实际的并行计算中会用到很多的分线程。
图中x和y表示点G的坐标，
Figure BDA0003840308360000041
表示模态值，下角标j，k，p都是点或面循环中使用的编号索引，i表示模态阶数。带括号的角标(s)表示当前点是结构节点，带括号的角标(fc)表示当前点是流体网格中心点，带括号的角标(fn)表示当前点是流体网格节点。
在耦合计算中，结构节点、流体网格中心点和流体网格节点是不重合的，需要通过插值实现数据传输。
本实施例采用了分步、分线程的并行化插值方法；首先，使用如下方法在主线程构建插值右矩阵V(ss)：
Figure BDA0003840308360000051
其中，Γ(||Gj,(s)-Gj,(s)||)是径向基核函数，||·||表示二范数运算，在实际使用中，用户可以选用多种不同形式的核函数；V(ss)的角标(ss)表示该矩阵在构建时只用到了结构节点的信息。
在实际计算中，每个分线程首先对飞行器表面所有流体网格中心点进行循环，读取每个网格的中心点坐标，并计算第一插值左矩阵V(sfc)；
第k个流体网格中心点对应的插值左矩阵V(sfc)通过以下方法计算：
V(sfc)＝[1 xk,(fc) yk,(fc) Γ(||G1,(s)-Gk,(fc)||)…Γ(||Gj,(s)-Gk,(fc)||)…Γ(||Gn,(s)-Gk,(fc)||)]
V(sfc)的角标(sfc)表示该矩阵在构建时同时用到了结构节点和流体网格中心点的信息。
第k个中心点对应的第i阶模态值
Figure BDA0003840308360000052
通过以下方法计算：
Figure BDA0003840308360000053
通过Fluent软件的用户自定义函数提供的接口，读取第k个中心点处的压力矢量pk和面积矢量Ak，采用以下方法计算第k个中心点对应的第i阶模态的广义气动力fi,k：
Figure BDA0003840308360000054
其中<·>表示内积运算。
然后，对每个分线程中的所有流体网格中心点的第i阶广义气动力求和，得到
Figure BDA0003840308360000055
接着，每个分线程将各自的第i阶广义气动力求和值fi分别发送到主线程，主线程对所有接收到的广义气动力求和值再次求和，得到总的第i阶广义气动力。
广义气弹方程具有如下形式：
Figure BDA0003840308360000056
其中，I是单位矩阵；
Figure BDA0003840308360000057
是对角矩阵，其中的第i个元素ωi表示第i阶模态的角频率；f＝[f1,...fi,...fn]T是所有广义气动力组成的列向量；q＝[q1,...qi,...qn]T是所有广义位移组成的向量；
上述的广义气弹方程组使用四阶龙格库塔方法求解，并将求解得到的q传递至每个分线程。
在主线程中，假设飞行器的飞行控制系统传感器的坐标为G(sen)＝(x(sen),y(sen))，该点的第i阶模态为
Figure BDA0003840308360000058
通过以下方法可以计算传感器处的结构位移、速度和加速度：
Figure BDA0003840308360000061
控制器在实际使用时，用户可以自由选择适当的控制律进行控制系统设计。将传感器处的结构位移、速度和加速度输入到控制系统中，控制系统计算得到反馈控制指令H；将控制指令带入操纵面偏转的型函数，即可得到操纵面偏转引起的流体网格节点的附加位移
Figure BDA0003840308360000062
在每个分线程中，对飞行器表面所有流体网格节点进行循环，读取每个网格节点的坐标，并计算第p个节点的第二插值左矩阵V(sfn)，进一步计算第p个节点的第i阶模态值
Figure BDA0003840308360000063
计算方法与流体网格中心点使用的方法相同；第p个流体网格节点的真实位移是
Figure BDA0003840308360000064
得到各流体网格节点的真实位移和附加位移后，通过Fluent用户自定义函数驱动流场计算边界变形，并通过Fluent的动网格技术求解和更新计算网格。
最后进行一步CFD并行求解，在未达到终止时间步之前，图主线程和分线程中的计算过程将反复执行。

Claims (2)

1.一种并行化的高精度时域气动弹性控制耦合仿真方法，其特征在于，具体步骤如下：
首先、使用结构有限元软件输出飞行器结构节点坐标和结构模态到txt文件中；
然后、生成飞行器外表面的流体网格，以并行模式启动计算流体力学软件Fluent，导入流体网格；
针对当前时间步，通过并行化径向基函数插值方法，将飞行器的结构模态插值，得到每个分线程中流体网格的中心点模态和流体网格节点模态；
每个分线程分别计算流体网格中心点的气动力，并将流体网格中心点的气动力与流体网格中心点模态相乘，对所有分线程的乘积求和，得到广义气动力传递给主线程；主线程利用广义气动力求解广义气动弹性运动方程组，得到每阶模态对应的广义位移和广义速度，并将其传递至各分线程；每个分线程使用流体网格节点模态和广义位移与广义速度，计算流体网格节点的真实位移与真实速度；
同时、主线程通过调用控制模块，将流体网格节点的真实位移与真实速度传递给控制模块，分别输出对应的操纵指令给各分线程；在每个分线程中，构建飞行器操纵面变形的型函数，使用操纵指令与型函数计算由操纵指令引起的流体网格节点的附加位移；
最后、在每个分线程中，使用Fluent软件的动网格技术计算由流体网格节点的真实位移和附加位移引起的网格变形，并更新每个分线程的网格；直至达到终止时间步，通过主线程输出各个时间步的广义位移和广义速度的txt文件。
2.如权利要求1所述的一种并行化的高精度时域气动弹性控制耦合仿真方法，其特征在于，所述每个分线程中流体网格的中心点模态和流体网格节点模态，具体计算过程为：
具体为：首先，在主线程中构建插值右矩阵，并传递给每个分线程；
然后，每个分线程分别各自对应飞行器一个子域的流体网格，通过读取流体网格的中心点坐标，并构建对应的第一插值左矩阵，结合插值右矩阵，将飞行器的结构模态插值得到当前流体网格中心点的模态。
同时、在每个分线程中，读取该线程内流场计算流体网格节点坐标，并构建第二插值左矩阵，结合插值右矩阵，将结构模态插值得到各流体网格节点的模态。
CN202211103424.1A 2022-09-09 2022-09-09 一种并行化的高精度时域气动弹性控制耦合仿真方法 Pending CN115407680A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211103424.1A CN115407680A (zh) 2022-09-09 2022-09-09 一种并行化的高精度时域气动弹性控制耦合仿真方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211103424.1A CN115407680A (zh) 2022-09-09 2022-09-09 一种并行化的高精度时域气动弹性控制耦合仿真方法

Publications (1)

Publication Number Publication Date
CN115407680A true CN115407680A (zh) 2022-11-29

Family

ID=84165757

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211103424.1A Pending CN115407680A (zh) 2022-09-09 2022-09-09 一种并行化的高精度时域气动弹性控制耦合仿真方法

Country Status (1)

Country Link
CN (1) CN115407680A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116151084A (zh) * 2023-04-21 2023-05-23 中国空气动力研究与发展中心计算空气动力研究所 基于结构网格的模拟方法、装置、终端设备及存储介质

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116151084A (zh) * 2023-04-21 2023-05-23 中国空气动力研究与发展中心计算空气动力研究所 基于结构网格的模拟方法、装置、终端设备及存储介质

Similar Documents

Publication Publication Date Title
McMillan et al. Efficient dynamic simulation of an underwater vehicle with a robotic manipulator
CN107633117A (zh) 基于高斯过程模型的全局灵敏度分析方法
Walker et al. MATRIX x: A data analysis, system identification, control design and simulation package
CN115407680A (zh) 一种并行化的高精度时域气动弹性控制耦合仿真方法
CN116258042B (zh) 一种基于ddm的大规模传热异构并行仿真方法
Budzko et al. Symbolic-numerical analysis of equilibrium solutions in a restricted four-body problem
Balmer et al. Elastoplastic and creep analysis with the ASKA program system
CN114676522A (zh) 融合gan和迁移学习的气动形状优化设计方法及系统及设备
CN111554355B (zh) 一种基于非冯诺依曼架构的分子动力学计算方法
CN112182737B (zh) 一种基于模态法的并行化高精度颤振计算方法
CN111563297A (zh) 一种基于bp网络的超声速飞行器热环境计算方法
McDaniel et al. Efficient mesh deformation for computational stability and control analyses on unstructured viscous meshes
CN112949104A (zh) 一种协作机器人实时模态分析方法
CN114970007A (zh) 一种支持cad模型导入与网格划分的刚柔耦合多体动力学通用化仿真方法
CN115758930A (zh) 一种桨-轴耦合系统双向流固耦合数值模拟方法
Luo et al. Rigid-flexible coupled dynamics analysis of 3-revolute-prismatic-spherical parallel robot based on multi-software platform
CN113626893B (zh) 一种基于隐式几何模型的计算机力学分析数值模拟方法
Chandar et al. CU++ ET: An Object Oriented Tool for Accelerating Computational Fluid Dynamics Codes using Graphical Processing Units
CN109885968B (zh) 一种复杂运动的运动描述建模方法及系统
US20240068903A1 (en) Generation of cfd-based structurally independent aerodynamic influence coefficient matrix
WO2024103241A1 (zh) 一种基于物质点法的软体机器人仿真方法
Sistla et al. Aircraft design problem implementation under the Common Object Request Broker Architecture
Chen Programming for solving plane rigid frame based on MATLAB
Luker et al. Air Force Research Laboratory progress in fluids-structures interaction
Easterling et al. Modular Aeroelastic Sensitivities: Accuracy and Performance

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination