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

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

Info

Publication number
CN109960865B9
CN109960865B9 CN201910212828.6A CN201910212828A CN109960865B9 CN 109960865 B9 CN109960865 B9 CN 109960865B9 CN 201910212828 A CN201910212828 A CN 201910212828A CN 109960865 B9 CN109960865 B9 CN 109960865B9
Authority
CN
China
Prior art keywords
matrix
gpu
thin plate
time
node
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.)
Active
Application number
CN201910212828.6A
Other languages
English (en)
Other versions
CN109960865B (zh
CN109960865A (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

Images

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

Abstract

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

Description

薄板无网格Galerkin结构动力响应分析的GPU加速方法
技术领域
本发明属于计算机辅助工程中薄板的结构动力仿真分析技术领域，具体涉及一种薄板无网格Galerkin(伽辽金)结构动力响应分析的GPU(GraphicProcessing 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α数据，根据边界条件对总体刚度矩阵进行修正
Figure GDA0004164508540000031
其中，
Figure GDA0004164508540000032
为修正后的总体刚度矩阵；
(5)读取GPU全局存储器中总体质量矩阵M、修正后的总体刚度矩阵
Figure GDA0004164508540000033
设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
Figure GDA0004164508540000034
(6)设置位移向量矩阵a、速度向量矩阵
Figure GDA0004164508540000035
加速度向量矩阵
Figure GDA0004164508540000036
以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、σ＝1/2、
Figure GDA0004164508540000037
Figure GDA0004164508540000038
c6＝Δt(1-δ)、c7＝δΔt，Δt为时间步长；
(7)读取GPU全局存储器中修正后的总体刚度矩阵
Figure GDA0004164508540000039
总体质量矩阵M、阻尼矩阵C，在GPU上形成有效刚度矩阵
Figure GDA0004164508540000041
即：
Figure GDA0004164508540000042
并对有效刚度矩阵进行三角分解，即：
Figure GDA0004164508540000043
其中，
Figure GDA0004164508540000044
为有效刚度矩阵，L为
Figure GDA0004164508540000045
三角分解后的下三角矩阵，LT为
Figure GDA0004164508540000046
三角分解后的上三角矩阵；
(8)求解每一时间步长(t＝0,Δt,2Δt,…)的位移、速度、加速度；其中，t为计算时间，Δt为时间步长；
(9)将主机内存中存储的所求节点每一时间步长的位移、速度、加速度的计算结果写入到输出文件中。
具体的，所述步骤(8)包括如下步骤：
(a)在CPU中读取位移向量矩阵a、速度向量矩阵
Figure GDA0004164508540000047
加速度向量矩阵
Figure GDA0004164508540000048
计算质量矩阵系数向量Mx与阻尼矩阵系数向量Cx：
Figure GDA0004164508540000049
Figure GDA00041645085400000410
并复制到GPU的全局存储器中；其中，at为t时刻位移向量，
Figure GDA00041645085400000411
为t时刻速度向量，
Figure GDA00041645085400000412
为t时刻加速度向量；
(b)读取GPU全局存储器中总体力向量矩阵F、总体质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷向量
Figure GDA00041645085400000413
即：
Figure GDA00041645085400000414
其中，Ft+Δt为t+Δt时刻的总体力向量；
(c)在GPU上采用回代的方法求解线性方程组，计算t+Δt时刻的位移向量at+Δt，即
Figure GDA00041645085400000415
其中，at+Δt为t+Δt时刻的位移向量，
Figure GDA00041645085400000416
为t+Δt时刻的有效载荷向量，L为
Figure GDA00041645085400000417
三角分解后的下三角矩阵，LT为
Figure GDA00041645085400000418
三角分解后的上三角矩阵；
(d)在GPU中计算t+Δt时刻的加速度向量
Figure GDA00041645085400000419
速度向量
Figure GDA00041645085400000420
即
Figure GDA00041645085400000421
(e)读取t+Δt时刻的位移向量at+Δt、速度向量
Figure GDA00041645085400000422
加速度向量
Figure GDA00041645085400000423
并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵
Figure GDA00041645085400000424
加速度向量矩阵
Figure GDA0004164508540000051
所对应的位置；
(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α数据，根据边界条件对总体刚度矩阵进行修正
Figure GDA0004164508540000061
其中，
Figure GDA0004164508540000062
为修正后的总体刚度矩阵。
(5)读取GPU全局存储器中总体质量矩阵M、修正后的总体刚度矩阵
Figure GDA0004164508540000063
设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
Figure GDA0004164508540000064
(6)设置位移向量矩阵a、速度向量矩阵
Figure GDA0004164508540000071
加速度向量矩阵
Figure GDA0004164508540000072
以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、σ＝1/2、
Figure GDA0004164508540000073
Figure GDA0004164508540000074
c6＝Δt(1-δ)、c7＝δΔt，Δt为时间步长。
(7)读取GPU全局存储器中修正后的总体刚度矩阵
Figure GDA0004164508540000075
总体质量矩阵M、阻尼矩阵C，在GPU上形成有效刚度矩阵
Figure GDA0004164508540000076
即：
Figure GDA0004164508540000077
并对有效刚度矩阵进行三角分解，即：
Figure GDA0004164508540000078
其中，
Figure GDA0004164508540000079
为有效刚度矩阵，L为
Figure GDA00041645085400000710
三角分解后的下三角矩阵，LT为
Figure GDA00041645085400000711
三角分解后的上三角矩阵。
(8)求解每一时间步长(t＝0,Δt,2Δt,…)的位移、速度、加速度；其中，t为计算时间，Δt为时间步长；
具体包括如下过程：
(a)在CPU中读取位移向量矩阵a、速度向量矩阵
Figure GDA00041645085400000712
加速度向量矩阵
Figure GDA00041645085400000713
计算质量矩阵系数向量Mx与阻尼矩阵系数向量Cx：
Figure GDA00041645085400000714
Figure GDA00041645085400000715
并复制到GPU的全局存储器中；其中，at为t时刻位移向量，
Figure GDA00041645085400000716
为t时刻速度向量，
Figure GDA00041645085400000717
为t时刻加速度向量；
(b)读取GPU全局存储器中总体力向量矩阵F、总体质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷向量
Figure GDA00041645085400000718
即：
Figure GDA00041645085400000719
其中，Ft+Δt为t+Δt时刻的总体力向量；
(c)在GPU上采用回代的方法求解线性方程组，计算t+Δt时刻的位移向量at+Δt，即
Figure GDA00041645085400000720
其中，at+Δt为t+Δt时刻的位移向量，
Figure GDA00041645085400000721
为t+Δt时刻的有效载荷向量，L为
Figure GDA00041645085400000722
三角分解后的下三角矩阵，LT为
Figure GDA00041645085400000723
三角分解后的上三角矩阵；
(d)在GPU中计算t+Δt时刻的加速度向量
Figure GDA0004164508540000081
速度向量
Figure GDA0004164508540000082
即
Figure GDA0004164508540000083
(e)读取t+Δt时刻的位移向量at+Δt、速度向量
Figure GDA0004164508540000084
加速度向量
Figure GDA0004164508540000085
并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵
Figure GDA0004164508540000086
加速度向量矩阵
Figure GDA0004164508540000087
所对应的位置；
(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
Figure GDA0004164508540000091
第二步：在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α数据，根据边界条件对总体刚度矩阵进行修正
Figure GDA0004164508540000111
第十步：在GPU中根据总体质量矩阵M、修正后的总体刚度矩阵
Figure GDA0004164508540000112
设置质量比例系数α、刚度比例系数β，在GPU上求解阻尼矩阵C，即
Figure GDA0004164508540000113
第十一步：在CPU中设置位移向量矩阵a、速度向量矩阵
Figure GDA0004164508540000114
加速度向量矩阵
Figure GDA0004164508540000115
以及参数α0、δ，并根据参数α0、δ计算积分常数c0、c1、c2、c3、c4、c5、c6、c7；其中，α0＝1/4、δ＝1/2、
Figure GDA0004164508540000116
Figure GDA0004164508540000117
c6＝Δt(1-δ)、c7＝δΔt，Δt为时间步长。
第十二步：在GPU中根据修正后的总体刚度矩阵
Figure GDA0004164508540000118
总体质量矩阵M、阻尼矩阵C，在GPU上形成有效的刚度矩阵
Figure GDA0004164508540000119
即：
Figure GDA00041645085400001110
第十三步：在GPU中对有效刚度矩阵进行三角分解，即：
Figure GDA00041645085400001111
第十四步：在CPU中读取位移向量矩阵a、速度向量矩阵
Figure GDA0004164508540000121
加速度向量矩阵
Figure GDA0004164508540000122
计算质量矩阵与阻尼矩阵的系数向量并赋值：
Figure GDA0004164508540000123
Figure GDA0004164508540000124
并复制到GPU的全局存储器中。
第十五步：在GPU中读取总体力向量矩阵F、质量矩阵M、阻尼矩阵C，计算t+Δt时刻的有效载荷
Figure GDA0004164508540000125
即
Figure GDA0004164508540000126
第十六步：在GPU中读取有效刚度矩阵三角分解数据L，采用回代的方法计算t+Δt时刻的位移at+Δt，即
Figure GDA0004164508540000127
第十七步：在GPU中计算t+Δt时刻的加速度
Figure GDA0004164508540000128
速度
Figure GDA0004164508540000129
即
Figure GDA00041645085400001210
第十八步：读取t+Δt时刻的位移向量at+Δt、速度向量
Figure GDA00041645085400001211
加速度向量
Figure GDA00041645085400001212
并由GPU的全局存储中复制到主机内存中位移向量矩阵a、速度向量矩阵
Figure GDA00041645085400001213
加速度向量矩阵
Figure GDA00041645085400001214
所对应的位置。
第十九步：判断t+Δt时刻是否达到所设置的总时间，若达到转到第二十步，若没达到，转到第十四步。
第二十步：将主机内存中存储的所求节点的每一时间步长的位移、速度、加速度的计算结果写入到输出文件中，本实施例薄板无网格Galerkin法结构动力响应分析的结果如图4表示(注：参考解为有限元在网格极细密情况下所求的解)。
采用上述方法进行薄板动力响应分析，整理CPU与GPU动力响应的计算时间(s)如表2所示：
表2
节点数 CPU GPU 加速比
9801 759.713 28.219 26.82
注：
Figure GDA0004164508540000131
虽然参考优先实施例对本发明进行描述，但以上所述并不限定本发明的保护范围，凡属于本发明思路下的技术方案均属于本发明的保护范围，任何在本发明的精神及原则内的修改、改进等，也应视为本发明的保护范围。

Claims (2)

1.一种薄板无网格Galerkin结构动力响应分析的GPU加速方法，其特征
(1)通过CPU读取薄板的节点坐标、背景积分网格、材料物理性能、约
(2)在CPU上将薄板的节点坐标、积分点数据、积分点定义域内的节点
(3)在CPU中根据设置的时间步长与总时间，读取薄板的节点坐标、作
(4)读取GPU全局存储器中总体刚度矩阵K和总体惩罚刚度矩阵Kα数
(5)读取GPU全局存储器中总体质量矩阵M、修正后的总体刚度矩阵
(6)设置位移向量矩阵a、速度向量矩阵
(7)读取GPU全局存储器中修正后的总体刚度矩阵
(8)求解每一时间步长(t＝0,Δt,2Δt,…)的位移、速度、加速度；其
(9)将主机内存中存储的所求节点每一时间步长的位移、速度、加速度
2.根据权利要求1所述薄板无网格Galerkin结构动力响应分析的GPU加
(a)在CPU中读取位移向量矩阵a、速度向量矩阵
(b)读取GPU全局存储器中总体力向量矩阵F、总体质量矩阵M、阻
(c)在GPU上采用回代的方法求解线性方程组，计算t+Δt时刻的位移
(d)在GPU中计算t+Δt时刻的加速度向量
(e)读取t+Δt时刻的位移向量at+Δt、速度向量
(f)判断t+Δt时刻是否达到所设置的总时间，若达到结束步骤(8)，若
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 CN109960865A (zh) 2019-07-02
CN109960865B CN109960865B (zh) 2023-03-21
CN109960865B9 true 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)

Families Citing this family (3)

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

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103970960A (zh) * 2014-05-23 2014-08-06 湘潭大学 基于gpu并行加速的无网格伽辽金法结构拓扑优化方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7499050B2 (en) * 2004-11-16 2009-03-03 Livermore Software Technology Corporation Method and system for mesh-free analysis of general three-dimensional shell structures
US7493243B2 (en) * 2004-12-27 2009-02-17 Seoul National University Industry Foundation Method and system of real-time graphical simulation of large rotational deformation and manipulation using modal warping

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103970960A (zh) * 2014-05-23 2014-08-06 湘潭大学 基于gpu并行加速的无网格伽辽金法结构拓扑优化方法

Also Published As

Publication number Publication date
CN109960865B (zh) 2023-03-21
CN109960865A (zh) 2019-07-02

Similar Documents

Publication Publication Date Title
CN109960865B9 (zh) 薄板无网格Galerkin结构动力响应分析的GPU加速方法
CN110457790A (zh) 用于结构变形分析的近场动力学非连续伽辽金有限元方法
JP5255714B2 (ja) 三次元の流体シミュレーション方法
CN105260581A (zh) 舰船机电控制设备虚拟振动和冲击试验方法
Jain et al. Hyper-reduction over nonlinear manifolds for large nonlinear mechanical systems
Goddeke et al. Co-processor acceleration of an unmodified parallel solid mechanics code with FEASTGPU
Eftekhari et al. Mixed finite element and differential quadrature method for free and forced vibration and buckling analysis of rectangular plates
US8548776B2 (en) Parallel physics solver
Sanfui et al. A two-kernel based strategy for performing assembly in FEA on the graphics processing unit
CN107480323B (zh) 基于LabVIEW和FPGA的向量式有限元并行求解仿真方法
CN109948253B (zh) 薄板无网格Galerkin结构模态分析的GPU加速方法
CN110147571B (zh) 一种组件结构的拓扑优化方法及装置
Oberhuber et al. TNL: Numerical library for modern parallel architectures
CN114792037A (zh) 一种超材料隔振器的序贯稳健性优化设计方法
CN102819454A (zh) 基于gpu的有限元显式并行求解仿真方法
CN103902764B (zh) 基于Householder变换的无约束结构静力分析方法
Ichimura et al. A fast scalable implicit solver with concentrated computation for nonlinear time-evolution problems on low-order unstructured finite elements
CN113467945B (zh) 基于无网格法拓扑优化的灵敏度并行及其gpu加速方法
Huang et al. Parallel Performance and Optimization of the Lattice Boltzmann Method Software Palabos Using CUDA
Wan et al. An efficient communication strategy for massively parallel computation in CFD
CN110781621A (zh) 一种含几何形状缺陷薄壁结构承载响应的快速重分析方法
Kuo et al. Development of a parallel explicit finite-volume Euler equation solver using the immersed boundary method with hybrid MPI-CUDA paradigm
CN114580144B (zh) 近场动力学问题的gpu并行实施方法
Wei et al. Acceleration of free-vibrations analysis with the Dual Reciprocity BEM based on ℋ-matrices and CUDA
CN114998091A (zh) 无网格法无矩阵对角预处理pcg求解及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