CN109960865A - 薄板无网格Galerkin结构动力响应分析的GPU加速方法 - Google Patents

薄板无网格Galerkin结构动力响应分析的GPU加速方法 Download PDF

Info

Publication number
CN109960865A
CN109960865A CN201910212828.6A CN201910212828A CN109960865A CN 109960865 A CN109960865 A CN 109960865A CN 201910212828 A CN201910212828 A CN 201910212828A CN 109960865 A CN109960865 A CN 109960865A
Authority
CN
China
Prior art keywords
matrix
gpu
vector
thin plate
time
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.)
Granted
Application number
CN201910212828.6A
Other languages
English (en)
Other versions
CN109960865B (zh
CN109960865B9 (zh
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.)
Xiangtan University
Original Assignee
Xiangtan 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 Xiangtan University filed Critical Xiangtan University
Priority to CN201910212828.6A priority Critical patent/CN109960865B9/zh
Publication of CN109960865A publication Critical patent/CN109960865A/zh
Publication of CN109960865B publication Critical patent/CN109960865B/zh
Application granted granted Critical
Publication of CN109960865B9 publication Critical patent/CN109960865B9/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2117/00—Details relating to the type or aim of the circuit design
    • G06F2117/08—HW-SW co-design, e.g. HW-SW partitioning
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/06—Power analysis or power optimisation
    • 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

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Complex Calculations (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种薄板无网格Galerkin结构动力响应分析的GPU加速方法。它通过CPU将薄板数据读入到主机内存，根据输入数据计算GPU并行所需相关数据，并将得到的数据复制到GPU全局存储器中；在GPU上设置线程块与线程数量，并行加速计算节点相应的形函数值，然后建立GPU线程块与交叉节点对一一映射模式，并行加速组装薄板的总体刚度矩阵与质量矩阵；在GPU上根据边界条件修正总体刚度矩阵，然后由质量矩阵和修正后的刚度矩阵求得阻尼矩阵、有效刚度矩阵，并将有效刚度矩阵采用三角分解；在GPU中采用Newmark法对薄板进行动力响应分析，最后输出求解的位移、速度和加速度的结果。本发明极大地提高了动力响应分析的求解效率。

Description

薄板无网格Galerkin结构动力响应分析的GPU加速方法
技术领域
本发明属于计算机辅助工程中薄板的结构动力仿真分析技术领域，具体涉及一种薄板无网格Galerkin(伽辽金)结构动力响应分析的GPU(Graphic Processing Unit，图形处理器)加速方法。
背景技术
薄板动力响应分析是研究薄板在动载荷的激励下所产生的位移、速度和加速度等力学参量的动态响应信息，是薄板结构动力设计与灵敏度分析及优化的重要基础之一。在工程中，许多领域都需要采用板壳结构，如车辆、船舶、飞机和压力容器等，然而，这些设备经常处于各种激励的作用，组成部件中的薄板不可避免地产生各种各样的振动。当薄板振动量超出允许范围之后，将会影响设备的工作性能，缩短设备的使用寿命，甚至造成严重的事故。因此，现代工程设计对薄板振动问题的解决提出了更严格的要求。
随着有限元法在工程领域的应用，薄板动力响应可以通过有限元法来模拟分析计算。但有限元法是基于网格的数值方法，网格的质量和尺寸严重影响计算的精度和计算工作量，而且复杂结构的网格生成也极具挑战性。近年来，无网格法得到迅速发展，它不需要划分网格，克服了有限元法对网格的依懒性，在涉及网格畸变时显示出明显优势。目前，无网格Galerkin法是应用最广的无网格计算方法，具有计算精度高、前后处理简便、自适应性好等优点，在结构动力学问题的分析中得到广泛应用。然而薄板的刚度矩阵和质量矩阵的组装耗时长、阻尼矩阵和总体离散系统方程求解速度慢等缺点，严重影响无网格Galerkin法在薄板结构动力分析中的应用。
近年来，随着计算技术和计算机技术的迅猛发展，CPU(Central ProcessingUnit)具有体积小、重量轻结构灵活、价格低廉等优点，但随着晶体管集成工艺进入瓶颈，晶体管密度增长达到物理极限。尽管Intel公司提出了多核CPU，但对CPU提高计算能力帮助有限，与此同时CPU芯片的时钟速率也已经达到了极限，很难再提高。然而，GPU在数据并行处理能力和存储器带宽上逐渐优于CPU，且GPU具有超高的浮点计算能力、性价比高、功耗小等优点，在通用计算领域取得了很大的进展。再者，GPU统一计算设备架构(CUDA)的常用编程模型是直接使用C语言进行开发，大大减少了编程工作。利用CUDA C调用GPU对薄板的刚度矩阵和质量矩阵并行组装，及薄板的阻尼矩阵和总体离散系统方程的求解，极大地提高了结构动力响应分析的计算效率，因此研究基于GPU并行加速的薄板无网格Galerkin法结构动力响应分析具有重要意义。
发明内容
本发明的目的在于针对薄板无网格Galerkin法结构动力响应分析中存在刚度矩阵与质量矩阵组装耗时长、阻尼矩阵和总体离散系统方程求解效率低等问题，提供一种薄板无网格Galerkin结构动力响应分析的GPU加速方法。
本发明的薄板无网格Galerkin结构动力响应分析的GPU加速方法，包括如下顺序的步骤：
(1)通过CPU读取薄板的节点坐标、背景积分网格、材料物理性能、约束条件和动态载荷数据，以及求解设置参数与动力响应的计算时间步长，并存储至主机内存中；在薄板内部和位移边界布置积分点，并求解薄板内部和位移边界上的积分点数据；通过对薄板的节点、积分点循环搜索，确定每个节点影响域内的积分点数据，以及每个积分点定义域内的节点数据；再对薄板的节点循环搜索，确定交叉节点对信息；然后，将薄板内部和位移边界上的积分点数据、交叉节点对信息存储至主机内存中；
(2)在CPU上将薄板的节点坐标、积分点数据、积分点定义域内的节点数据、交叉节点对信息由主机内存复制到GPU全局存储器中，设置GPU中的线程块与线程的数量，并行加速计算节点的形函数值、形函数一阶导数值及形函数二阶导数值；并行加速组装总体刚度矩阵K、总体惩罚刚度矩阵Kα和总体质量矩阵M，且以按行压缩存储CSR格式存储于GPU全局存储器中；
(3)在CPU中根据设置的时间步长与总时间，读取薄板的节点坐标、作用力以及力作用的时间，依次求得Newmark法每一时间步长下的总体外力向量，并对应存储于总体外力向量矩阵F中，最后由主机内存复制至GPU全局存储器中；
(4)读取GPU全局存储器中总体刚度矩阵K和总体惩罚刚度矩阵Kα数据，根据边界条件对总体刚度矩阵进行修正其中，为修正后的总体刚度矩阵；
(5)读取GPU全局存储器中总体质量矩阵M、修正后的总体刚度矩阵设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
(6)设置位移向量矩阵a、速度向量矩阵加速度向量矩阵以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、σ＝1/2、 c6＝Δt(1-δ)、c7＝δΔt，Δt为时间步长；
(7)读取GPU全局存储器中修正后的总体刚度矩阵总体质量矩阵M、阻尼矩阵C，在GPU上形成有效刚度矩阵即：并对有效刚度矩阵进行三角分解，即：其中，为有效刚度矩阵，L为三角分解后的下三角矩阵，LT为三角分解后的上三角矩阵；
(8)求解每一时间步长(t＝0,Δt,2Δt,…)的位移、速度、加速度；其中，t为计算时间，Δt为时间步长；
(9)将主机内存中存储的所求节点每一时间步长的位移、速度、加速度的计算结果写入到输出文件中。
具体的，所述步骤(8)包括如下步骤：
(a)在CPU中读取位移向量矩阵a、速度向量矩阵加速度向量矩阵计算质量矩阵系数向量Mx与阻尼矩阵系数向量Cx： 并复制到GPU的全局存储器中；其中，at为t时刻位移向量，为t时刻速度向量，为t时刻加速度向量；
(b)读取GPU全局存储器中总体力向量矩阵F、总体质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷向量即：其中，Ft+Δt为t+Δt时刻的总体力向量；
(c)在GPU上采用回代的方法求解线性方程组，计算t+Δt时刻的位移向量at+Δt，即其中，at+Δt为t+Δt时刻的位移向量，时刻的有效载荷向量，L为三角分解后的下三角矩阵，LT为三角分解后的上三角矩阵；
(d)在GPU中计算t+Δt时刻的加速度向量速度向量即
(e)读取t+Δt时刻的位移向量at+Δt、速度向量加速度向量并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵加速度向量矩阵所对应的位置；
(f)判断t+Δt时刻是否达到所设置的总时间，若达到结束步骤(8)，若没达到，转到步骤(a)。
本发明与现有技术相比，具有如下有益效果：
(1)本发明在GPU上同时对总体质量矩阵与总体刚度矩阵进行组装，且阻尼矩阵与有效刚度矩阵的计算也在GPU上实现，因此极大地缩短了薄板动力响应分析过程中矩阵的组装计算时间；同时，本发明采用稀疏格式存储总体质量矩阵、总体刚度矩阵、阻尼矩阵与有效刚度矩阵，减小了存储需求，能够实现较大规模问题的计算。
(2)本发明采用Newmark法对薄板动力响应进行分析。在GPU上对薄板有效刚度矩阵仅进行一次三角分解，位移、速度与加速度响应值均在GPU上进行循环迭代计算，极大地减少了薄板动力响应分析的求解时间。
(3)本发明采用无网格Galerkin法，不需要网格，对薄板动力响应分析结果具有较高精度，而且本发明方法可广泛适用于二维任意几何形状的薄板动力响应分析问题。
附图说明
图1为本发明方法的主流程框图。
图2为本发明实施例的薄板示意图
图3为图2所示实施例的薄板节点示意图。
图4为本发明方法的薄板动力响应分析曲线图。
具体实施方式
下面结合附图和实施例对本发明作进一步详细的描述。
参见图1，本发明薄板无网格Galerkin结构动力响应分析的GPU加速方法，包括如下顺序的步骤：
(1)通过CPU读取薄板的节点坐标、背景积分网格、材料物理性能、约束条件和动态载荷数据，以及求解设置参数与动力响应的计算时间步长，并存储至主机内存中；在薄板内部和位移边界布置积分点，并求解薄板内部和位移边界上的积分点数据；通过对薄板的节点、积分点循环搜索，确定每个节点影响域内的积分点数据，以及每个积分点定义域内的节点数据；再对薄板的节点循环搜索，确定交叉节点对信息；然后，将薄板内部和位移边界上的积分点数据、交叉节点对信息存储至主机内存中。
(2)在CPU上将薄板的节点坐标、积分点数据、积分点定义域内的节点数据、交叉节点对信息由主机内存复制到GPU全局存储器中，设置GPU中的线程块与线程的数量，并行加速计算节点的形函数值、形函数一阶导数值及形函数二阶导数值；并行加速组装总体刚度矩阵K、总体惩罚刚度矩阵Kα和总体质量矩阵M，且以按行压缩存储CSR格式存储于GPU全局存储器中。
(3)在CPU中根据设置的时间步长与总时间，读取薄板的节点坐标、作用力以及力作用的时间，依次求得Newmark法每一时间步长下的总体外力向量，并对应存储于总体外力向量矩阵F中，最后由主机内存复制至GPU全局存储器中。
(4)读取GPU全局存储器中总体刚度矩阵K和总体惩罚刚度矩阵Kα数据，根据边界条件对总体刚度矩阵进行修正其中，为修正后的总体刚度矩阵。
(5)读取GPU全局存储器中总体质量矩阵M、修正后的总体刚度矩阵设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
(6)设置位移向量矩阵a、速度向量矩阵加速度向量矩阵以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、σ＝1/2、 c6＝Δt(1-δ)、c7＝δΔt，t为时间步长。
(7)读取GPU全局存储器中修正后的总体刚度矩阵总体质量矩阵M、阻尼矩阵C，在GPU上形成有效刚度矩阵即：并对有效刚度矩阵进行三角分解，即：其中，为有效刚度矩阵，L为三角分解后的下三角矩阵，LT为三角分解后的上三角矩阵。
(8)求解每一时间步长(t＝0,Δt,2Δt,…)的位移、速度、加速度；其中，t为计算时间，Δt为时间步长；
具体包括如下过程：
(a)在CPU中读取位移向量矩阵a、速度向量矩阵加速度向量矩阵计算质量矩阵系数向量Mx与阻尼矩阵系数向量Cx： 并复制到GPU的全局存储器中；其中，at为t时刻位移向量，为t时刻速度向量，为t时刻加速度向量；
(b)读取GPU全局存储器中总体力向量矩阵F、总体质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷向量即：其中，Ft+Δt为t+Δt时刻的总体力向量；
(c)在GPU上采用回代的方法求解线性方程组，计算t+Δt时刻的位移向量at+Δt，即其中，at+Δt为t+Δt时刻的位移向量，为t+Δt时刻的有效载荷向量，L为三角分解后的下三角矩阵，LT为三角分解后的上三角矩阵；
(d)在GPU中计算t+Δt时刻的加速度向量速度向量即
(e)读取t+Δt时刻的位移向量at+Δt、速度向量加速度向量并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵加速度向量矩阵所对应的位置；
(f)判断t+Δt时刻是否达到所设置的总时间，若达到结束步骤(8)，若没达到，转到步骤(a)。
(9)将主机内存中存储的所求节点每一时间步长的位移、速度、加速度的计算结果写入到输出文件中。
下面是将本发明方法应用于工程上的一个实例：
参见图2，是一块边长为2m、厚度为0.1m的四边固定约束的方形薄板，材料的弹性模量为2.0×108Pa，泊松比μ＝0.3，密度为ρ＝7800kg/m3，中间承受10kN恒力，力的作用时间是0.5s，时间步长0.001s。无网格Galerkin法动力响应分析的薄板包括9801个节点，其分布如图3所示。本发明方法针对本实施例的具体实施步骤如下：
第一步：将已经准备好的薄板数据(包括薄板的节点坐标、背景积分网格顶点的编号和坐标、边界节点对中两节点的编号、边界信息、惩罚系数)、材料物理参数(包括材料弹性模量、泊松比、材料密度)、节点影响域系数、作用力以及力的作用时间，输入到input_Dynamic.txt文件中，在程序运行时，系统自动将以上数据读入到主机内存中。
一个标准的input_Dynamic.txt文件格式如下表1：
表1
第二步：在CPU中开辟9604(本实施例中背景积分网格个数)*4(每个四边形背景积分网格内布置的积分点数)*4*sizeof(double)字节的存储空间用于存储薄板内部积分点信息的数组，然后根据本实施例四边形背景积分网格顶点的坐标与编号、位移边界节点，在薄板内部布置积分点，并求得积分点的坐标、权系数及雅可比系数。
第三步：在CPU中开辟396(本实施例中位移边界节点对个数)*4(每个单元节点对上布置的积分点数)*4*sizeof(double)字节的存储空间用于存储位移边界积分点数组，然后通过CPU读取主机内存中存储的薄板背景积分网格的顶点坐标与编号，边界节点对中两节点的编号，在边界上布置积分点，求得边界上积分点的坐标、权系数、雅克比系数。
第四步：在CPU中开辟9801(模型节点的个数)*2(x、y两个方向)*sizeof(double)字节的存储空间用于存储节点的动态影响域半径数组，然后对薄板的节点循环，确定其它节点和该节点的最小距离rmin，再利用节点影响域半径＝3.0(节点影响域系数)*rmin，求得该节点的影响域半径，同理求得积分点的定义域半径。
第五步：在CPU中读取薄板的节点坐标、积分点坐标、节点的影响域半径、积分点的定义域半径，根据节点影响域半径，确定组成交叉节点对信息、每个节点影响域内的积分点数据以及每个积分点定义域内的节点数据；并将薄板的节点坐标、积分点坐标、交叉节点对信息、每个节点影响域内的积分点数据以及每个积分点定义域内的节点数据由主机内存复制到GPU全局存储器中，并为它们分配和主机内存中同样大小的存储空间。
第六步：读取GPU全局存储器中存储的薄板的节点坐标、积分点数据、积分点定义域内的节点数据，建立线程与积分点一一映射模式，通过最小二乘逼近，并行求得每个积分点定义域内节点的形函数值、形函数对应于各个坐标轴的二阶导数值，在本实施例中，为形函数对x、y轴的二阶导数值。
第七步：在GPU中建立线程块和交叉节点对一一映射模式，并行加速组装总体质量矩阵M、总体刚度矩阵K和总体惩罚刚度矩阵Kα，并以按行压缩存储CSR格式存储于GPU全局存储器中。
第八步：在CPU中根据设置的时间步长0.001s与总时间0.5s，以及存储在主机中薄板的节点坐标、作用力以及力作用的时间，依次求得Newmark法每一步长下的总体外力向量，并对应存储于总体外力向量矩阵F中，最后复制到GPU全局存储器中。
第九步：在GPU中读取GPU全局存储器中存储的总体刚度矩阵K和总体惩罚刚度矩阵Kα数据，根据边界条件对总体刚度矩阵进行修正
第十步：在GPU中根据总体质量矩阵M、修正后的总体刚度矩阵设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
第十一步：在CPU中设置位移向量矩阵a、速度向量矩阵加速度向量矩阵以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、δ＝1/2、 c6＝Δt(1-δ)、c7＝δΔt，Δt为时间步长。
第十二步：在GPU中根据修正后的总体刚度矩阵总体质量矩阵M、阻尼矩阵C，在GPU上形成有效的刚度矩阵即：
第十三步：在GPU中对有效刚度矩阵进行三角分解，即：
第十四步：在CPU中读取位移向量矩阵a、速度向量矩阵加速度向量矩阵计算质量矩阵与阻尼矩阵的系数向量并赋值： 并复制到GPU的全局存储器中。
第十五步：在GPU中读取总体力向量矩阵F、质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷即
第十六步：在GPU中读取有效刚度矩阵三角分解数据L，采用回代的方法计算t+Δt时刻的位移at+Δt，即
第十七步：在GPU中计算t+Δt时刻的加速度速度即
第十八步：读取t+Δt时刻的位移向量at+Δt、速度向量加速度向量并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵加速度向量矩阵所对应的位置。
第十九步：判断t+Δt时刻是否达到所设置的总时间，若达到转到第二十步，若没达到，转到第十四步。
第二十步：将主机内存中存储的所求节点的每一时间步长的位移、速度、加速度的计算结果写入到输出文件中，本实施例薄板无网格Galerkin法结构动力响应分析的结果如图4表示(注：参考解为有限元在网格极细密情况下所求的解)。
采用上述方法进行薄板动力响应分析，整理CPU与GPU动力响应的计算时间(s)如表2所示：
表2
节点数 CPU GPU 加速比
9801 759.713 28.219 26.82
注：
虽然参考优先实施例对本发明进行描述，但以上所述并不限定本发明的保护范围，凡属于本发明思路下的技术方案均属于本发明的保护范围，任何在本发明的精神及原则内的修改、改进等，也应视为本发明的保护范围。

Claims (2)

1.一种薄板无网格Galerkin结构动力响应分析的GPU加速方法，其特征在于包括如下顺序的步骤：
(1)通过CPU读取薄板的节点坐标、背景积分网格、材料物理性能、约束条件和动态载荷数据，以及求解设置参数与动力响应的计算时间步长，并存储至主机内存中；在薄板内部和位移边界布置积分点，并求解薄板内部和位移边界上的积分点数据；通过对薄板的节点、积分点循环搜索，确定每个节点影响域内的积分点数据，以及每个积分点定义域内的节点数据；再对薄板的节点循环搜索，确定交叉节点对信息；然后，将薄板内部和位移边界上的积分点数据、交叉节点对信息存储至主机内存中；
(2)在CPU上将薄板的节点坐标、积分点数据、积分点定义域内的节点数据、交叉节点对信息由主机内存复制到GPU全局存储器中，设置GPU中的线程块与线程的数量，并行加速计算节点的形函数值、形函数一阶导数值及形函数二阶导数值；并行加速组装总体刚度矩阵K、总体惩罚刚度矩阵Kα和总体质量矩阵M，且以按行压缩存储CSR格式存储于GPU全局存储器中；
(3)在CPU中根据设置的时间步长与总时间，读取薄板的节点坐标、作用力以及力作用的时间，依次求得Newmark法每一时间步长下的总体外力向量，并对应存储于总体外力向量矩阵F中，最后由主机内存复制至GPU全局存储器中；
(4)读取GPU全局存储器中总体刚度矩阵K和总体惩罚刚度矩阵Kα数据，根据边界条件对总体刚度矩阵进行修正其中，为修正后的总体刚度矩阵；
(5)读取GPU全局存储器中总体质量矩阵M、修正后的总体刚度矩阵设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
(6)设置位移向量矩阵a、速度向量矩阵加速度向量矩阵以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、σ＝1/2、 c6＝Δt(1-δ)、c7＝δΔt，Δt为时间步长；
(7)读取GPU全局存储器中修正后的总体刚度矩阵总体质量矩阵M、阻尼矩阵C，在GPU上形成有效刚度矩阵即：并对有效刚度矩阵进行三角分解，即：其中，为有效刚度矩阵，L为三角分解后的下三角矩阵，LT为三角分解后的上三角矩阵；
(8)求解每一时间步长(t＝0,Δt,2Δt,…)的位移、速度、加速度；其中，t为计算时间，Δt为时间步长；
(9)将主机内存中存储的所求节点每一时间步长的位移、速度、加速度的计算结果写入到输出文件中。
2.根据权利要求1所述薄板无网格Galerkin结构动力响应分析的GPU加速方法，其特征在于：所述步骤(8)包括如下步骤：
(a)在CPU中读取位移向量矩阵a、速度向量矩阵加速度向量矩阵计算质量矩阵系数向量Mx与阻尼矩阵系数向量Cx： 并复制到GPU的全局存储器中；其中，at为t时刻位移向量，为t时刻速度向量，为t时刻加速度向量；
(b)读取GPU全局存储器中总体力向量矩阵F、总体质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷向量即：其中，Ft+Δt为t+Δt时刻的总体力向量；
(c)在GPU上采用回代的方法求解线性方程组，计算t+Δt时刻的位移向量at+Δt，即其中，at+Δt为t+Δt时刻的位移向量，为t+Δt时刻的有效载荷向量，L为三角分解后的下三角矩阵，LT为三角分解后的上三角矩阵；
(d)在GPU中计算t+Δt时刻的加速度向量速度向量即
(e)读取t+Δt时刻的位移向量at+Δt、速度向量加速度向量并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵加速度向量矩阵所对应的位置；
(f)判断t+Δt时刻是否达到所设置的总时间，若达到结束步骤(8)，若没达到，转到步骤(a)。
CN201910212828.6A 2019-03-20 2019-03-20 薄板无网格Galerkin结构动力响应分析的GPU加速方法 Active CN109960865B9 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910212828.6A CN109960865B9 (zh) 2019-03-20 2019-03-20 薄板无网格Galerkin结构动力响应分析的GPU加速方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910212828.6A CN109960865B9 (zh) 2019-03-20 2019-03-20 薄板无网格Galerkin结构动力响应分析的GPU加速方法

Publications (3)

Publication Number Publication Date
CN109960865A true CN109960865A (zh) 2019-07-02
CN109960865B CN109960865B (zh) 2023-03-21
CN109960865B9 CN109960865B9 (zh) 2023-05-05

Family

ID=67024651

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910212828.6A Active CN109960865B9 (zh) 2019-03-20 2019-03-20 薄板无网格Galerkin结构动力响应分析的GPU加速方法

Country Status (1)

Country Link
CN (1) CN109960865B9 (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112571150A (zh) * 2020-12-09 2021-03-30 中南大学 一种用于监测薄板齿轮的薄板加工状态的非线性方法
CN115329250A (zh) * 2022-10-13 2022-11-11 中国空气动力研究与发展中心计算空气动力研究所 基于dg处理数据的方法、装置、设备及可读存储介质
CN115995277A (zh) * 2023-03-22 2023-04-21 中国空气动力研究与发展中心计算空气动力研究所 一种材料动力学特性评估方法、装置、设备及介质

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20060103648A1 (en) * 2004-11-16 2006-05-18 Livermore Software Technology Corporation Method and system for mesh-free analysis of general three-dimensional shell structures
US20060139347A1 (en) * 2004-12-27 2006-06-29 Choi Min G Method and system of real-time graphical simulation of large rotational deformation and manipulation using modal warping
CN103970960A (zh) * 2014-05-23 2014-08-06 湘潭大学 基于gpu并行加速的无网格伽辽金法结构拓扑优化方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20060103648A1 (en) * 2004-11-16 2006-05-18 Livermore Software Technology Corporation Method and system for mesh-free analysis of general three-dimensional shell structures
US20060139347A1 (en) * 2004-12-27 2006-06-29 Choi Min G Method and system of real-time graphical simulation of large rotational deformation and manipulation using modal warping
CN103970960A (zh) * 2014-05-23 2014-08-06 湘潭大学 基于gpu并行加速的无网格伽辽金法结构拓扑优化方法

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112571150A (zh) * 2020-12-09 2021-03-30 中南大学 一种用于监测薄板齿轮的薄板加工状态的非线性方法
CN115329250A (zh) * 2022-10-13 2022-11-11 中国空气动力研究与发展中心计算空气动力研究所 基于dg处理数据的方法、装置、设备及可读存储介质
CN115329250B (zh) * 2022-10-13 2023-03-10 中国空气动力研究与发展中心计算空气动力研究所 基于dg处理数据的方法、装置、设备及可读存储介质
CN115995277A (zh) * 2023-03-22 2023-04-21 中国空气动力研究与发展中心计算空气动力研究所 一种材料动力学特性评估方法、装置、设备及介质

Also Published As

Publication number Publication date
CN109960865B (zh) 2023-03-21
CN109960865B9 (zh) 2023-05-05

Similar Documents

Publication Publication Date Title
CN109960865B (zh) 薄板无网格Galerkin结构动力响应分析的GPU加速方法
Farhat et al. Robust and provably second‐order explicit–explicit and implicit–explicit staggered time‐integrators for highly non‐linear compressible fluid–structure interaction problems
CN110457790A (zh) 用于结构变形分析的近场动力学非连续伽辽金有限元方法
JP5255714B2 (ja) 三次元の流体シミュレーション方法
CN105260581A (zh) 舰船机电控制设备虚拟振动和冲击试验方法
CN107480357B (zh) 一种多杆系结构有限元模型边界模拟方法
Yuan et al. An equivalent modeling method for honeycomb sandwich structure based on orthogonal anisotropic solid element
Goddeke et al. Co-processor acceleration of an unmodified parallel solid mechanics code with FEASTGPU
CN111400898B (zh) 基于主模态法和应变的阵列天线振动变形预测方法及设备
CN114792037B (zh) 一种超材料隔振器的序贯稳健性优化设计方法
Liu et al. A wavelet multiresolution interpolation Galerkin method with effective treatments for discontinuity for crack growth analyses
CN109948253B (zh) 薄板无网格Galerkin结构模态分析的GPU加速方法
Zhang et al. A universal quadrilateral shell element for the absolute nodal coordinate formulation
Rezaiee-Pajand et al. Analytical scheme for solid stress analysis
Ichimura et al. A fast scalable implicit solver with concentrated computation for nonlinear time-evolution problems on low-order unstructured finite elements
CN116227155B (zh) 一种推进剂细观力学性能研究方法
Zhang et al. Linear and geometrically nonlinear analysis of plates and shells by a new refined non-conforming triangular plate/shell element
CN112733415B (zh) 一种薄壁弹性体边界的无网格处理方法、装置、终端设备及计算介质
CN103745038A (zh) 摩天轮结构有限元模型的建模方法
Biswas et al. Global load balancing with parallel mesh adaption on distributed-memory systems
Zapata et al. A GPU parallel finite volume method for a 3D Poisson equation on arbitrary geometries
Bogdanov et al. Adaptation and optimization of basic operations for an unstructured mesh CFD algorithm for computation on massively parallel accelerators
Playne et al. Benchmarking multi-GPU communication using the shallow water equations
Okumoto et al. Finite element method
CN114580144B (zh) 近场动力学问题的gpu并行实施方法

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
GR01 Patent grant
GR01 Patent grant
CI03 Correction of invention patent
CI03 Correction of invention patent

Correction item: Description

Correct: Paragraphs 1-79 of the specification submitted on April 6, 2023

False: Paragraphs 1-79 of the specification for the original documents submitted on the application date

Number: 12-01

Page: ??

Volume: 39