CN110457790B - 用于结构变形分析的近场动力学非连续伽辽金有限元方法 - Google Patents

Abstract

本发明属于近场动力学数据处理技术领域，公开了一种用于结构变形分析的近场动力学非连续伽辽金有限元方法，建立实体模型，确定各材料区域并赋予相应材料属性；采用标准有限单元法的网格生成算法划分单元网格；结点复制、生成匹配新方法的非连续单元；初始计算与实时更新结点与高斯点坐标；引入权函数构造键型近场动力学的弱形式，通过高斯点计算单元结点等效载荷、等效外载荷和质量矩阵，确定单元的受力变形计算方法；设置初始条件，施加应力和位移边界条件；提交静动力计算；应力、应变分析。本发明能够在传统有限元框架内实现运用键型近场动力学理论求解结构静动力变形破坏问题。

Description

用于结构变形分析的近场动力学非连续伽辽金有限元方法

技术领域

本发明属于近场动力学数据处理技术领域，尤其涉及一种用于结构变形分析的近场动力学非连续伽辽金有限元方法。

背景技术

固体结构断裂破坏模拟是计算力学界长期关注的问题，对该类不连续问题的常规建模分析方法是基于偏微分方程的传统连续介质力学理论及其有限元法，但该方法由于缺少描述断裂的长度尺度参数、需要预设裂纹路径与裂纹扩展准则，故面临不连续处奇异性问题、网格依赖性问题、断裂描述不准确与计算精度瓶颈；后又发展出扩展有限元法和各类无网格法，但扩展有限元法在处理复杂三维断裂时面临裂纹面追踪困难，无网格方法在结构大变形模拟分析中具有优势，但仍存在计算精度和计算效率不高、以及不连续处的奇异性问题。

新兴的非局部近场动力学(Peridynamics,简记为PD)理论采用空间积分方程描述物质的力学行为，引入了描述结构断裂的长度尺度参数，避免了传统数值计算方法在求解不连续问题时的奇异性，在分析裂纹扩展与断裂破坏问题时表现出先天的优势。

然而，现有近场动力学的数值求解多采用无网格粒子法，该求解方法要求将结构离散为高密度的均匀离散点阵，离散粒子或结点数量众多以及PD的非局部特性导致现有方法的计算量巨大，对于计算机主频、内存和存储空间提出更高要求；无网格粒子离散结构的方法，相较于网格离散方法存在边界处的离散误差问题；且相比于传统局部有限元，也存在边界效应问题与应力边界条件施加困难的不足。连续伽辽金有限单元法是实现传统连续介质力学理论与开展数值模拟的主要方法，但传统伽辽金有限元方法在分析不连续力学问题时存在不足。因此，有必要发展基于非局部近场动力学理论的非连续伽辽金有限元方法，该方法保持近场动力学分析不连续问题的优势、也可以利用成熟的有限元计算模式并提高计算效率，具有很大的优势，但是近场动力学的非连续伽辽金有限元解法研究不足，缺乏对于一般三维问题的理论建模和数值实施的算法细节，在Chen X,Gunzburger M.Continuousand discontinuous finite element methods for a peridynamics model ofmechanics.Computer Methods in Applied Mechanics and Engineering.2011；200(9-12):1237-1250首次提出一维算法后，Ren B,Wu CT,Askari E.A 3D discontinuousGalerkin finite element method with the bond-based peridynamics model fordynamic brittle failure analysis.International Journal of ImpactEngineering.2017；99:14-25对该方法做了进一步探索，但缺乏具体的算法和数值实施细节、不能有效指导工程实践，并且LS-DYNA商用软件对用户存在极大限制。

综上所述，现有技术存在的问题是：

(1)传统连续介质力学理论与有限元等数值方法在处理不连续问题时面临精度和效率的瓶颈。

(2)现有近场动力学的无网格粒子求解方法计算量巨大、对计算机性能提出更高要求；无网格粒子方法离散结构存在复杂边界的离散误差问题；现有方法还存在边界效应显著、应力边界条件准确施加的困难。

(3)近场动力学的非连续伽辽金有限元解法研究不足，缺乏对于一般三维问题的理论建模和数值实施的算法细节，实例应用与分析也不足。

解决上述技术问题的难点在于：针对一般性三维结构变形问题，需要构建近场动力学的非连续伽辽金有限元的理论模型与算法公式，给出完善的数值求解体系，包括矩阵构造计算方法、边界条件施加方法、各变量和应力应变的计算方法等。

解决上述技术问题的意义：发展近场动力学的有限元解法能够完善近场动力学的数值求解体系、增强与传统有限元方法的对比联系，更为深入开展近场动力学方法与传统有限元法的耦合建模奠定基础，有效减小计算工作量、提高计算精度与计算效率，同时方法的适用性广泛，能满足各种复杂模型的设计计算需求。

发明内容

针对现有技术存在的问题，本发明提供了一种用于结构变形分析的近场动力学非连续伽辽金有限元方法。

本发明是这样实现的，一种用于结构变形分析的近场动力学非连续伽辽金有限元方法。所述用于结构变形分析的近场动力学非连续伽辽金有限元方法，包括如下步骤：

步骤一，建立实体模型，确定各材料区域并赋予相应材料属性；

步骤二，采用标准有限单元法的网格生成算法划分单元网格，要求单元尺寸满足计算精度要求；

步骤三，采用结点复制算法、生成匹配近场动力学的非连续伽辽金有限元这种新方法的非连续单元；

步骤四，初始时刻计算、并在后续计算步中实时更新单元结点与高斯点坐标；

步骤五，引入权函数构造键型近场动力学的弱形式方程，引入等参单元插值形函数以便用结点变量值表征任意点变量，根据高斯点受其近场域范围内的其它高斯点的非局部作用、以及高斯点间相对位移计算单元结点等效载荷，计算应力边界和体力的结点等效外载荷，计算结点质量矩阵，给出整个构型的受力变形计算方法；

步骤六，设置初始条件，采用拉格朗日乘子法施加体力、应力和位移边界条件；

步骤七，提交静动力计算；

步骤八，由变形结果开展应力、应变分析。

进一步，上述步骤二和步骤三中，可按如下方式生成单元网格：

1)对于整个模型进行标准有限元网格划分、共有N个单元，其中，二维模型采用平面四结点等参元，三维模型采用六面体八节点等参元，要求单元尺寸大致相当，且单元尺寸平均为d＝L/50～L/400，L为所述实体模型最大边长；

2)提取标准单元和结点信息，采用结点复制算法生成匹配新方法的非连续单元，即在同一位置生成属于各单元的结点，则共有4N或8N个结点。

进一步，上述步骤四中，根据结点初始坐标计算高斯点初始坐标，根据结点位移实时计算高斯点的新坐标和位移，所述高斯点坐标及位移计算方法如下：

其中为高斯点坐标和位移，i表示高斯点编号，X_i为高斯点i所在单元的所有结点坐标，d_i为高斯点i所在单元的结点位移，为第i个高斯点所在单元的形函数矩阵；

其中，高斯点所在单元形函数矩阵和单元结点位移列向量如下：

进一步，上述步骤五中，引入权函数构造键型近场动力学的基于单元离散化的伽辽金弱形式控制方程为：

其中，ng为单元内高斯点的总数，ng'是x^g的紧支域内的高斯点数，ξ(x^g)是高斯点x^g与在其近场范围内的相邻高斯点x^g'之间的键，ΔV^g是高斯点x^g体积的增量；考虑某高斯点近场范围内的所有高斯点相互作用完成第一重求和，对高斯点所在单元其它高斯点完成第二重求和，则所述的结构静动力变形的矩阵形式控制方程为：

[K]_3n×3n[d]_3n×1+[F]_3n×1＝0，

其中，K为结构整体劲度矩阵，d为所有结点的位移列向量，F为结点等效外荷载列阵，M为结点质量矩阵。

质量矩阵M、整体劲度矩阵K和载荷向量F的计算方法如下：

进一步，整体劲度矩阵由单元劲度矩阵按结点编号组装得到，单元劲度矩阵可由矩阵形式的单元平衡方程对应获得，根据高斯点受其近场域范围内的其它高斯点的非局部作用计算矩阵形式的键力表达和单元劲度矩阵；所述方法的单元平衡方程如下：

[K^e]_12×24[d]_24×1+[b]_12×1＝0，

其中，

其中，i遍历计算域内所有的高斯点，j遍历高斯点i的近场范围内的高斯点；为高斯点所在单元的形函数矩阵，c为键型近场动力学的微观模量，T为局部坐标到整体坐标的转置矩阵。

进一步，两个高斯点之间存在非局部相互作用或本构力，在整体坐标下矩阵形式表达具体如下：

其中，f_ij为整体坐标下高斯点i和j之间的键力；为两高斯点构成的键的劲度矩阵，为局部坐标到整体坐标的转置矩阵，且有 和分别为整体坐标下高斯点j和i的位移向量；

且

进一步，上述步骤六中，对于静力问题，采用拉格朗日乘子法施加体力、应力和位移边界条件GU+U^*＝0，其中G为已知约束方程系数矩阵，U^*为已知位移约束值；步骤七采用大型稀疏方程组求解器PARDISO等求解代数方程组即可获得未知量U和拉格朗日乘子λ；对于动力变形问题，给定初始时刻位移和速度，步骤七通过显式中心差分算法迭代，获得t+Δt时间的位移解答，即

进一步，在由步骤七获得静动力计算结果后，在步骤八中采用近场动力学微分算子定义非局部变形梯度和应力，开展应变、应力结果分析；所述的非局部变形梯度和应力如下：

其中，g为近场动力学函数，λ和μ为拉梅常数。

综上所述，本发明的优点及积极效果为：

(1)本发明首次给出固体变形的键型近场动力学的非连续伽辽金静动力有限元解法，拓展完善了近场动力学的数值求解体系，为近场动力学方法与传统有限元的混合建模奠定基础；

(2)本发明给出了近场动力学弱形式控制方程的有限元解法，相比于传统的近场动力学强形式的无网格粒子解法，有限元法可以直接施加应力和位移边界条件、避免边界条件施加困难，形函数插值方法可以获得全场结果；

(3)本发明采用有限元网格离散结构，避免了无网格粒子方法在边界处的离散误差；

(4)本发明能够高效准确得求解结构静力变形、动力变形、弹性波传播与开裂破坏问题，计算精度和效率高、适用性广泛，满足各种复杂模型的设计计算需求。

附图说明

图1是本发明实施例提供的用于结构变形分析的近场动力学非连续伽辽金有限元的方法流程图。

图2是本发明实施例提供的用于结构变形分析的近场动力学非连续伽辽金有限元方法的原理图。

图3是本发明实施例提供的实体模型网格单元划分示意图。

图4(a)是本发明实施例提供的实体模型的传统有限元网格划分示意图；图4(b)是本发明实施例提供的非连续网格示意图。

图5是本发明实施例提供的计算结点等效载荷的关键积分示意图；表示高斯点受其近场域内其他高斯点的作用。

图6(a)为本发明实施例提供的计算结果的水平位移分布示意图；图6(b)为本发明实施例提供的计算结果的竖向位移分布示意图。

图7为本发明实施例提供的计算结果的水平正应变分布示意图。

图8为本发明实施例提供的计算结果的水平正应力分布示意图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

针对现有技术存在的问题，本发明提供了一种用于结构变形分析的近场动力学非连续伽辽金有限元方法，下面结合附图对本发明作详细的描述。

本发明的目的在于解决传统连续介质力学体系处理不连续问题面临的瓶颈，以及现有近场动力学数值求解体系不完善、及其非连续伽辽金有限元解法研究缺乏的问题，本发明提供一种用于固体结构变形分析的键型近场动力学的非连续伽辽金有限元建模与求解方法，该方法结合了近场动力学方法在求解不连续问题和有限元方法计算效率的优势，提高了计算精度与效率，并为近场动力学与传统有限元的混合建模奠定基础。

下面结合附图对本发明的技术方案作详细描述。

如图1所示，本发明实施例提供的用于结构变形分析的近场动力学非连续伽辽金有限元方法，包括如下步骤：

S101：建立实体模型，确定各材料区域并赋予相应材料属性；

S102：采用标准有限单元法的网格生成算法划分单元网格，要求单元尺寸满足计算精度要求；

S103：采用结点复制算法、生成匹配近场动力学的非连续伽辽金有限元这种新方法的非连续单元；

S104：初始时刻计算、并在后续计算步中实时更新单元结点与高斯点坐标；

S105：引入权函数构造键型近场动力学的弱形式方程，引入等参单元插值形函数以便用结点变量值表征任意点变量，根据高斯点受其近场域范围内的其它高斯点的非局部作用、以及高斯点间相对位移计算单元结点等效载荷，计算应力边界和体力的结点等效外载荷，计算结点质量矩阵，给出整个构型的受力变形计算方法；

S106：设置初始条件，采用拉格朗日乘子法施加体力、应力和位移边界条件；

S107：提交静动力计算；

S108：由变形结果开展应力、应变分析。

进一步，上述步骤(2)和(3)中，可按如下方式生成单元网格：

1)对于整个模型进行标准有限元网格划分、共有N个单元，其中，二维模型采用平面四结点等参元，三维模型采用六面体八节点等参元，要求单元尺寸大致相当，且单元尺寸平均为d＝L/50～L/400，L为所述实体模型最大边长；

2)提取标准单元和结点信息，采用结点复制算法生成匹配新方法的非连续单元，即在同一位置生成属于各单元的结点，则共有4N或8N个结点。

进一步，上述步骤(4)中，根据结点初始坐标计算高斯点初始坐标，根据结点位移实时计算高斯点的新坐标和位移，所述高斯点坐标及位移计算方法如下：

其中为高斯点坐标和位移，i表示高斯点编号，X_i为高斯点i所在单元的所有结点坐标，d_i为高斯点i所在单元的结点位移，为第i个高斯点所在单元的形函数矩阵；

其中，高斯点所在单元形函数矩阵和单元结点位移列向量如下：

进一步，上述步骤(5)中，引入权函数构造键型近场动力学的基于单元离散化的伽辽金弱形式控制方程为：

其中，ng为单元内高斯点的总数，ng'是x^g的紧支域内的高斯点数，ξ(x^g)是高斯点x^g与在其近场范围内的相邻高斯点x^g'之间的键，ΔV^g是高斯点x^g体积的增量；考虑某高斯点近场范围内的所有高斯点相互作用完成第一重求和，对高斯点所在单元其它高斯点完成第二重求和，则所述的结构静动力变形的矩阵形式控制方程为：

[K]_3n×3n[d]_3n×1+[F]_3n×1＝0，

其中，K为结构整体劲度矩阵，d为所有结点的位移列向量，F为结点等效外荷载列阵，M为结点质量矩阵。

质量矩阵M、整体劲度矩阵K和载荷向量F的计算方法如下：

进一步，整体劲度矩阵由单元劲度矩阵按结点编号组装得到，单元劲度矩阵可由矩阵形式的单元平衡方程对应获得，根据高斯点受其近场域范围内的其它高斯点的非局部作用计算矩阵形式的键力表达和单元劲度矩阵；所述方法的单元平衡方程如下：

[K^e]_12×24[d]_24×1+[b]_12×1＝0，

其中，

其中，i遍历计算域内所有的高斯点，j遍历高斯点i的近场范围内的高斯点；为高斯点所在单元的形函数矩阵，c为键型近场动力学的微观模量，T为局部坐标到整体坐标的转置矩阵。

进一步，两个高斯点之间存在非局部相互作用或本构力，在整体坐标下矩阵形式表达具体如下：

其中，f_ij为整体坐标下高斯点i和j之间的键力；为两高斯点构成的键的劲度矩阵，为局部坐标到整体坐标的转置矩阵，且有 和分别为整体坐标下高斯点j和i的位移向量；

且

进一步，上述步骤(6)中，对于静力问题，采用拉格朗日乘子法施加体力、应力和位移边界条件GU+U^*＝0，其中G为已知约束方程系数矩阵，U^*为已知位移约束值；步骤(7)采用大型稀疏方程组求解器PARDISO等求解代数方程组即可获得未知量U和拉格朗日乘子λ；对于动力变形问题，给定初始时刻位移和速度，步骤(7)通过显式中心差分算法迭代，获得t+Δt时间的位移解答，即

进一步，在由步骤(7)获得静动力计算结果后，在步骤(8)中采用近场动力学微分算子定义非局部变形梯度和应力，开展应变、应力结果分析；所述的非局部变形梯度和应力如下：

其中，g为近场动力学函数，λ和μ为拉梅常数。

下面结合具体实施例对本发明的技术方案作进一步说明。

本发明的用于结构变形分析的近场动力学非连续伽辽金有限元法，建模与求解方法流程如图2，该方法用于结构变形分析具有很高的可靠性。

实施例：

本实施例以研究尺寸为1000mm×500mm的二维含孔矩形板在水平均匀单轴拉伸载荷作用下的静力变形为例，利用本发明的方法进行近场动力学建模与分析，圆孔直径为100mm；杨氏模量为E＝200GPa，泊松比为υ＝1/3，质量密度为ρ＝7850kg/m³；载荷形式为沿长边方向的拉伸位移0.5mm，上下水平边界自由无约束，如图3。

建模方法主要步骤如下：

(1)建立一个实体模型，模型的外轮廓尺寸为1000mm×500mm、小孔直径为100mm，赋予材料相同的材料参数E、v等；

(2)对于整个模型进行标准有限元网格划分，如图4(a)采用平面四结点等参元、取网格平均尺寸d＝L/50＝20mm，L为所述实体模型最大边长，圆孔周围单元较细小，共离散得到N＝1312个标准单元、1400个结点；

(3)提取标准单元和结点信息，采用结点复制算法生成匹配新方法的非连续单元，如图4(b)，共生成1312个非连续单元、5248个结点，最终划分的网格模型如图3；

(4)单元内采用2×2高斯点，初始时刻计算、并在后续计算步中实时更新单元结点与高斯点坐标，具体为

其中，高斯点所在单元形函数矩阵和单元结点位移列向量如下：

(5)设置初始条件：所有结点初始位移u＝0，采用拉格朗日乘子法施加左边界固定约束与右边界位移边界条件0.5mm；

引入权函数构造键型近场动力学的非连续伽辽金有限元弱形式的控制方程为

其中，ng为单元内高斯点的总数，ng'是x^g的紧支域内的高斯点数，ξ(x^g)是高斯点x^g与在其近场范围内的相邻高斯点x^g'之间的键，ΔV^g是高斯点x^g体积的增量。

结构静动力变形的矩阵形式控制方程为：

其中，d为结构所有结点的位移列向量，K为整体劲度矩阵，F为结点等效外荷载列阵，M为质量矩阵：

其中，整体劲度矩阵由单元劲度矩阵按结点编号组装得到，单元劲度矩阵可由矩阵形式的单元平衡方程对应获得；

在整体坐标下，单元的平衡方程为：[K^e]_12×24[d]_24×1+[b]_12×1＝0，其中：

其中，i遍历计算域内所有的高斯点，j遍历高斯点i的近场范围内的高斯点；c为键型近场动力学的微观模量，为高斯点所在单元的形函数矩阵；

根据高斯点受其近场域范围内的其它高斯点的非局部作用，如图5，计算高斯点i和j之间的矩阵形式键力表达

其中，两高斯点构成的键的劲度矩阵局部到整体的转换矩阵

且

(7)提交静动力计算；对于本实施例的静力变形问题，采用拉格朗日乘子施加位移边界条件GU+U^*＝0，其中G为已知约束方程系数矩阵，U^*为已知位移约束值；采用大型稀疏方程组求解器求解代数方程组即可获得未知量U和拉格朗日乘子λ；水平位移和竖向位移分布如图6。

对于动力变形问题，给定初始时刻位移和速度，通过显式中心差分算法迭代，重复上述步骤，即可获得t+Δt时间的位移解答，即

由静动力计算的位移结果，计算非局部变形梯度F、应变ε和应力σ，进行应变、应力结果分析，如图7和图8，应变应力结果精度高。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，所述用于结构变形分析的近场动力学非连续伽辽金有限元方法包括如下步骤：

步骤一，建立实体模型，确定各材料区域并赋予相应材料属性；

步骤二，采用标准有限单元法的网格生成算法划分单元网格；

步骤三，采用结点复制算法、生成匹配近场动力学的非连续伽辽金有限元这种新方法的非连续单元；

步骤四，初始时刻计算、并在后续计算步中实时更新单元结点与高斯点坐标；

步骤五，引入权函数构造键型近场动力学的弱形式方程，引入等参单元插值形函数以便用结点变量值表征任意点变量，根据高斯点受其近场域范围内的其它高斯点的非局部作用、以及高斯点间相对位移计算单元结点等效载荷，计算应力边界和体力的结点等效外载荷，计算结点质量矩阵，给出整个构型的受力变形计算方法；

步骤六，设置初始条件，采用拉格朗日乘子法施加体力、应力和位移边界条件；

步骤七，提交静动力计算；

步骤八，由变形结果开展应力、应变分析。

2.如权利要求1所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，在步骤二和步骤三中，按如下方式生成单元网格：

1)对于整个模型进行标准有限元网格划分、共有N个单元，其中，二维模型采用平面四结点等参元，三维模型采用空间六面体八节点等参元，要求单元尺寸大致相当，且单元尺寸平均为d＝L/50～L/400，L为所述实体模型最大边长；

2)提取标准单元和结点信息，采用结点复制算法、在同一位置生成属于各单元的结点，则共有4N或8N个结点，即生成匹配新方法的非连续单元。

3.如权利要求1所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，步骤四中，根据结点初始坐标计算高斯点初始坐标，根据结点位移实时计算高斯点的新坐标和位移，所述高斯点坐标及位移计算方法如下：

其中分别为高斯点坐标和位移，i表示高斯点编号，X_i为高斯点i所在单元的所有结点坐标，d_i为高斯点i所在单元的结点位移，为第i个高斯点所在单元的形函数矩阵。

4.如权利要求3所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，高斯点所在单元形函数矩阵和单元结点位移列向量如下：

5.如权利要求1所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，步骤五中，引入权函数构造键型近场动力学的基于单元离散化的伽辽金弱形式控制方程为：

其中，ng为一个单元内高斯点的总数，ng'是x^g的紧支域内的高斯点数，ξ(x^g)是高斯点x^g与在其近场范围内的相邻高斯点x^g'之间的键，ΔV^g是高斯点x^g体积的增量；考虑某高斯点近场范围内的所有高斯点相互作用完成第一重求和，对高斯点所在单元其它高斯点完成第二重求和，结构静动力变形的矩阵形式控制方程为：

[K]_3n×3n[d]_3n×1+[F]_3n×1＝0，

其中，K为结构整体劲度矩阵，d为所有结点的位移列向量，F为结点等效外荷载列阵，M为结点质量矩阵。

6.如权利要求5所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，质量矩阵M、整体劲度矩阵K和载荷向量F的计算方法如下：

7.如权利要求6所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，整体劲度矩阵由单元劲度矩阵按结点编号组装得到，单元劲度矩阵可由矩阵形式的单元平衡方程对应获得，根据高斯点受其近场域范围内的其它高斯点的非局部作用计算矩阵形式的键力表达和单元劲度矩阵，所述方法的单元平衡方程如下：

[K^e]_12×24[d]_24×1+[b]_12×1＝0

其中，

其中，i遍历计算域内所有的高斯点，j遍历高斯点i的近场范围内的高斯点；为高斯点所在单元的形函数矩阵，c为键型近场动力学的微观模量，T为局部坐标到整体坐标的转置矩阵。

8.如权利要求7所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，两个高斯点之间存在非局部相互作用或本构力，在整体坐标下矩阵形式表达具体如下：

其中，f_ij为整体坐标下高斯点i和j之间的键力；为两高斯点构成的键的劲度矩阵，为局部坐标到整体坐标的转置矩阵，且有 和分别为整体坐标下高斯点j和i的位移向量，且有

9.如权利要求1所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，步骤六中，对于静力问题，采用拉格朗日乘子法施加体力、应力和位移边界条件GU+U^*＝0，其中G为已知约束方程系数矩阵，U^*为已知位移约束值；步骤七采用大型稀疏方程组求解器PARDISO等求解代数方程组即可获得未知量U和拉格朗日乘子λ；对于动力变形问题，给定初始时刻位移和速度，步骤七通过显式中心差分算法迭代，获得t+Δt时间的位移解答，即

10.如权利要求1所述的用于结构变形分析的近场动力学非连续伽辽金有限元方法，其特征在于，在由步骤七获得静动力计算结果后，在步骤八中采用近场动力学微分算子定义非局部变形梯度和应力，开展应变、应力结果分析；所述的非局部变形梯度和应力如下：

其中，g为近场动力学函数，λ和μ为拉梅常数。

Images