CN116011301B - B样条等几何状态空间有限元方法 - Google Patents

B样条等几何状态空间有限元方法 Download PDF

Info

Publication number
CN116011301B
CN116011301B CN202310109283.2A CN202310109283A CN116011301B CN 116011301 B CN116011301 B CN 116011301B CN 202310109283 A CN202310109283 A CN 202310109283A CN 116011301 B CN116011301 B CN 116011301B
Authority
CN
China
Prior art keywords
equation
spline
stress
boundary
finite element
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
CN202310109283.2A
Other languages
English (en)
Other versions
CN116011301A (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.)
Donghua University
Original Assignee
Donghua 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 Donghua University filed Critical Donghua University
Priority to CN202310109283.2A priority Critical patent/CN116011301B/zh
Publication of CN116011301A publication Critical patent/CN116011301A/zh
Application granted granted Critical
Publication of CN116011301B publication Critical patent/CN116011301B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Complex Calculations (AREA)

Abstract

本发明涉及计算力学领域,公开了B样条等几何状态空间有限元方法,在Hellinger‑Reissner变分方程基础上,仅在复合材料叠层板各层的最上面划分单元,并使用B样条插值位移及应力分量,而沿厚度方向使用精细积分法,高效迭代得到转移矩阵,从而将基于变分方程的非齐次一阶微分方程转变为线形方程组,本发明保持了状态空间有限元法的降维特点,因而保证计算效率,同时等几何分析过程使得CAD建模与CAE分析过程一体化,提高前处理效率,本发明具有简单、高效、高精度的特点,具有较强通用性,易于在复合材料叠层结构失效分析中推广使用。

Description

B样条等几何状态空间有限元方法
技术领域
本发明属于计算力学领域,具体涉及B样条等几何状态空间有限元方法。
背景技术
等几何方法(Isogeometric Analysis)是一种新型的有限元增强手段。主要将计算机辅助设计(Computer Aided Design,CAD)中构造模型所采用的形函数,如B样条形函数或非均匀有理B样条(Non-Uniform Rational B-Splines,NURBS)形函数,直接引入至有限元的等参元插值中,从而消除了几何误差,并缩减了有限元划分单元的时间。同时,B样条或NURBS形函数具备的迭代升阶特性,使得描述结构几何外形更加方便、精准。
目前,强度高且性能好的复合叠层结构被广泛用于航空航天、汽车交通等工程领域,具有数以亿计的商用价值。需要注意的是,层间是叠层结构的薄弱环节,层间应力的准确计算是优化复合材料叠层板设计的必要过程。然而,商用有限元软件在计算层间应力时,应力场是基于位移及应变场而求解的,因而结果在层间两侧不连续,与实际情况相违背。同时,传统有限元法在处理复杂边界条件时存在模拟不精确的劣势,亟待提升模拟效果。为此,提升单元阶次是有效途径,而传统Lagrange单元在单元边界处仅有C^0连续,且无法自由升阶,导致复杂边界处物理场的模拟精度提升有限。
上述解析或数值方法已经充分展现了状态空间法在求解链式结构物理场的优势。为了更方便的结合前处理与分析过程,本发明引入了等几何分析方法,形成B样条状态空间有限元法。该方法结合了状态空间有限元法的离散方便的优势,同时兼具B样条自由升降阶的特点,可以准确处理更为复杂的边界条件,且B样条的引入使得结构设计及数值分析一体化,省略了划分单元的时间,提升计算效率
发明内容
为解决上述技术问题,本发明提供B样条等几何状态空间有限元方法,以解决现有技术中的问题,为实现上述发明目的,本发明所采用的技术方案是:
B样条等几何状态空间有限元方法,包括如下步骤:
步骤一,对于复合材料叠层板,对层间界面采用B样条建模,并提取该界面的控制点坐标,并计算对应的B样条形函数;
步骤二,确定复合材料叠层板的材料参数、边界条件以及叠层板所受的温度载荷;
步骤三,利用Hellinger-Reissner混合变分方程,得到位移、应力、形变之间的方程;
步骤四,引入几何方程、物理方程、平衡方程,将步骤三中的方程的形变变量变为应力变量,分别计算层间切应力及层间两侧的膜应力;
步骤五,利用步骤一所提取的B样条控制点坐标及B样条形函数,插值复合材料叠层结构层间界面上的位移及应力场;
步骤六,利用步骤二确定的叠层板边界条件,引入B样条插值,构建边界点与相同单元内其他节点间的插值关系,从而消去边界处的多余变量,得到满足特定边界条件的方程;
步骤七,对于步骤六中得到的方程,通过变分原理提取等式;
步骤八,将步骤七所得方程转变为非齐次一阶微分方程;
步骤九,根据步骤八中方程的通解形式,得到叠层材料中某一层材料上下表面分量之间的关系式,并通过连接条件,得到包含所有未知分量的线性方程组;
步骤十,求解线性方程组,得到叠层板层间界面的位移及应力场。
本发明具有以下有益效果:
1.保持了状态空间有限元法的降维特点的基础上,通过建模分析过程一体化,从而提高前处理分析效率,因而计算效率更高;
2.将B样条插值用于等几何分析,且保持了优秀的计算精度
3.通用性较好,在复合材料叠层结构失效分析中能够广泛使用。
附图说明
图1为本方法的流程图;
图2为本方法的叠层结构示意图;
图3为本方法温度载荷示意图;
图4为本方法使用3阶B样条插值时单元离散的示意图;
图5为本方法使用2阶B样条计算四边简支边界条件下的位移u云图与有限元结果的对比;
图6为本方法使用3阶B样条计算x边自由y边固支条件下的应力σ_yz的云图与有限元结果的对比。
具体实施方式
下面将结合本发明实施例中的图1-图6,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅是本发明一部分实施例,而不是全部的实施例,若未特别指明,实施例中所用的技术手段为本领域技术人员所熟知的常规手段。
本发明提供了基于B样条的新式计算热弹性力学的计算方法,该方法在水平方向使用传统有限元法并用等几何B样条单元离散,在竖直方向上使用状态空间法推导得到计算方程,该方法在计算精度保持与有限元软件一个水平的同时,计算效率大大优于有限元法。
本发明具体涉及一种B样条等几何状态空间有限元方法,包括如下步骤:
步骤一,对于复合材料叠层板,对层间界面采用B样条建模,并提取该界面的控制点坐标,并计算对应的B样条形函数;
步骤一中的复合材料叠层结构如图2所示,结构尺寸为长、宽均为2m,每层结构的高度h均相同,材料的上层和下层均为碳纤维树脂,材料的中间层为环氧树脂,材料均为正交各向异性的。
步骤二,确定复合材料叠层板的材料参数、边界条件以及叠层板所受的温度载荷;
步骤三,利用Hellinger-Reissner混合变分方程,得到位移、应力、形变之间的方程;
假设边界条件全部满足,得到两个Hellinger-Reissner混合变分方程,即式(1)和式(2):
其中,u是位移分量,ε及σ分别表示应变及应力分量,f是体力分量,δ是变分符号,为针对正交各向异性材料的微分算子,即微分算子/>和/>所组成的矩阵,dV表示对应的体积微元;
步骤四,引入几何方程、物理方程、平衡方程,将步骤三中的方程,即式(2)的形变变量变为应力变量,分别计算层间切应力及层间两侧的膜应力;
式(1)和式(2)可看作是平衡方程及几何方程,结合热弹性力学中的物理方程:
ε=Sσ+J#(3)
其中JT=[αxT,αyT,αzT,0,0,0]是温度载荷,S为柔度矩阵,也是C刚度矩阵的逆矩阵。
叠层结构所示的温度载荷如图3所示,将式(3)代入到式(2)得到几何方程式,代入后,式(1)和式(2)的分量均只有位移和应力。
分离两个方程,即式(1)和代入物理方程之后式(2)的剪切应力和弯曲应力,把剪切应力分量和位移分量组合在一起,并和平衡方程式(1)结合,得到:
其中A、B、C、D、E、G、H等均为常数矩阵,T是温度载荷。
步骤五,利用步骤一所提取的几何结构的B样条控制点坐标及形函数,插值复合材料叠层结构层间界面上的位移及应力场;
首先划分B样条等几何单元,如图4所示,当B样条阶数为2时,本方法中沿x轴和y轴使用的节点向量坐标均为:
[0,0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,1,1];
控制点向量坐标均为:
[0,0.05,0.15,0.25,0.35,0.45,0.55,0.65,0.75,0.85,0.95,1];
而当阶数为3时,沿x轴和y轴的节点坐标同样均为:
[0,0,0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1,1,1];
控制点向量的坐标均为:
[0,0.0333,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.96666];
B样条形函数的构成如式(4)所示。
其中,Ni,(ξ)为第i个p阶形函数在局部坐标ξ上的值。
本方法只需要在水平面上划分单元,不需要和有限元方法一样在竖直方向上划分单元,大大减少了单元数量和前处理繁琐程度。
划分单元之后,一个单元内部任意坐标的点的分量可以用作用于该单元的控制点乘上对应的基函数值求得。例如以第r个单元的位移u和应力σxz为例:
其中Qi(ξ,η)为平面单元下形函数的值,是x轴方向和y轴方向各自形函数值的乘积。
离散化之后,步骤五中就可以将式(1)、式(2)改写为:
其中:
pT=[uT,vT,wT]T,q=[σxz Tyz Tzz T]T,s=σxx Txy Tyy T]T
步骤六,利用步骤二确定的叠层板边界条件,引入B样条插值,构建边界点与相同单元内其他节点间的插值关系,从而消去边界处的多余变量,使得式(12)、式(13)能得到满足特定边界条件的方程;
本方法共研究了四种边界条件,分别是
x边与其对应边简支,y边与其对应边也为简支;
x边与其对应边自由,y边与其对应边简支;
x边与其对应边简支,y边与其对应边固支;
x边与其对应边自由,y边与其对应边固支。
简支、固支和自由边三种条件下,位移和应力的初始条件也不同。简支条件下不固定法向位移,固定其他方向位移为0;固支条件下固定所有位移均为0;自由边条件下不限制位移。
步骤七,对于步骤六中得到的方程,即已经满足边界条件的式(12)和式(13),通过变分原理提取等式;
代入边界条件,即分离已知的和未知的分量,再根据变分原理,步骤七中满足边界条件的式(12)和式(13)可以写为:
下标为f的矩阵和向量代表已经代入了边界条件。
步骤八,将步骤七所得方程,即式(14)和式(15)转变为非齐次一阶微分方程;
通过变量替换组合式(14)和式(15),得到
式中:
这是一个非齐次的一阶线性微分方程,对叠层材料的任意一层来说,都能列出形如此式的线性微分方法。在竖直方向上使用数值方法求解该微分方程。
步骤九,根据步骤八中方程,即式(16)的通解形式,得到叠层材料中某一层材料上下表面分量之间的关系式,并通过连接条件,得到包含所有未知分量的线性方程组;
对于步骤八中得到的微分方程,即式(16),其格式与状态方程一致,故可以使用状态方程的通解对其进行计算:
但通解法的问题是状态转移矩阵φi(z)需要用到泰勒展开式计算,由于本方法中未知数较多,计算状态转移矩阵所需的时间十分长,并且最后计算结果精确度也不高,故而本方法为了改进使用精细积分法。
对于叠层材料的第i层,如果该层厚度h足够小,那么ri(z)就可以近似表示为(ri(zi)+ri(zi-1))/2;而如果厚度不是足够小,那么可以人为地给该层分层,之后在子层中也可以这样表示ri(z)。我们将任意一层分为Ki=2k个子层,每一个子层的厚度就是Δi=h/2k,以叠层材料第i层的第一个子层为例,对微分方程的左右两边均取积分处理:
式中的O=ΔiMi/2,
考虑到层和层之间的连接条件,即上一个子层的下层和下一个子层的上层是一层,则它们的物理分量是一样的,于是对于一层内的任意子层都施加如上的积分处理,可以得到
本式的数值都可以只通过迭代k次高效求得。之后对叠层材料每一层都作如上处理并根据层与层之间的连接条件,可以得到最终需要计算的一个线性方程组,包含了所有未知的分量。精细积分法的特点在于不需要计算状态转移矩阵φi(z),即不需要一个需要展开很多项的复杂矩阵,从而大大减少了计算时间,并且计算精度更高。
需要注意的是本方法要求式(21)中的A是正定的,这需要未知分量中相应的位移和应力分量的个数一致,然而在自由边和固支边的边界条件下并不满足这一点,这就需要额外的操作。自由边下可以由式(8)进行消除多余的分量,而固支边可以由式(5)-(7)消除掉多余的分量。
步骤十,求解线性方程组,得到叠层板层间界面的位移及应力场。
求解线性方程组就能够求得叠层材料每一层上下表面的位移和应力分量的值,即叠层材料层间的受力情况。层间应力往往引起叠层材料边缘脱黏,形成层间裂纹,造成整个叠层材料刚度和强度下降,使结构过早失效,故准确地分析叠层材料层间的位移和应力在工程界具有重大意义。
图5和图6展示了使用本方法计算叠层材料力学特征以及使用传统有限元方法的结果云图,由图中可以看出,本方法的计算精度与有限元法十分接近。
本发明保持了状态空间有限元法的降维特点的基础上,通过建模分析过程一体化,从而提高前处理分析效率,因而计算效率更高;通过将B样条插值用于等几何分析,且保持了优秀的计算精度;通用性较好,在复合材料叠层结构失效分析中能够广泛使用。
本发明提供了一种分析复合材料叠层板力学场的B样条插值的等几何状态空间有限元方法。在Hellinger-Reissner变分方程基础上,仅在复合材料叠层板层间界面划分单元,并使用B样条插值位移及应力分量,而沿厚度方向使用精细积分法,高效迭代得到转移矩阵,从而将基于变分方程的非齐次一阶微分方程转变为线形方程组。本发明保持了状态空间有限元法的降维特点,因而保证计算效率。同时等几何分析过程使得CAD建模与CAE分析过程一体化,提高前处理效率。本发明具有简单、高效、高精度的特点,具有较强通用性,易于在复合材料叠层结构失效分析中推广使用。
进一步的,所述步骤五中,通过以下方法消除边界条件中的多余变量:
简直边的边界条件中无多余变量需要处理;
固支边的边界条件中切应力分量为多余变量,根据物理方程消除:
自由边的边界条件中仅应力分量σxx为多余变量,根据如下物理方程消除:
以上所述的实施例仅是对本发明的优选方式进行描述,并非对本发明的范围进行限定,在不脱离本发明设计精神的前提下,本领域普通技术人员对本发明的技术方案做出的各种变形、变型、修改、替换,均应落入本发明权利要求书确定的保护范围内。

Claims (2)

1.B样条等几何状态空间有限元方法,其特征在于,包括如下步骤:
步骤一,对于复合材料叠层板,对层间界面采用B样条建模,并提取该界面的控制点坐标,并计算对应的B样条形函数;
步骤二,确定复合材料叠层板的材料参数、边界条件以及叠层板所受的温度载荷;
步骤三,利用Hellinger-Reissner混合变分方程,得到位移、应力、形变之间的方程;
步骤四,引入几何方程、物理方程、平衡方程,将步骤三中的方程的形变变量变为应力变量,分别计算层间切应力及层间两侧的膜应力;
步骤五,利用步骤一所提取的该界面的B样条控制点坐标及B样条形函数,插值复合材料叠层结构层间界面上的位移及应力场;
步骤六,利用步骤二确定的叠层板边界条件,引入B样条插值,构建边界点与相同单元内其他节点间的插值关系,从而消去边界处的多余变量,得到满足边界条件的方程;
步骤七,对于步骤六中得到的方程,通过变分原理提取等式;
步骤八,将步骤七所得方程转变为非齐次一阶微分方程;
步骤九,根据步骤八中方程的通解形式,得到叠层材料中某一层材料上下表面分量之间的关系式,并通过连接条件,得到包含所有未知分量的线性方程组;
步骤十,求解线性方程组,得到叠层板层间界面的位移及应力场;
所述步骤三中的Hellinger-Reissner混合变分方程为两个,分别为:
其中,u是位移分量,ε及σ分别表示应变及应力分量,f是体力分量,δ是变分符号,为针对正交各向异性材料的微分算子,dV表示对应的体积微元;
结合以下热力学中的物理方程,
ε=Sσ+JT (3)
其中JT=[αxT,αyT,αzT,0,0,0]是温度载荷,αx表示在x方向上的热膨胀系数,αy表示在y方向上的热膨胀系数,αz表示在z方向上的热膨胀系数,S为柔度矩阵,也是C刚度矩阵的逆矩阵;
所述步骤一中,通过迭代构造B样条形函数:
其中,Ni,p(ξ)为第i个p阶形函数在局部坐标ξ上的值,ξi+p为第i+p个结点的局部坐标,ξ代表局部坐标,N代表B样条的形函数,ξi为第i个结点的局部坐标,Ni+1,p-1(ξ)为p-1阶B样条第i+1个形函数在局部坐标ξ处的值。
2.根据权利要求1所述的B样条等几何状态空间有限元方法,其特征在于,所述步骤五中,通过以下方法消除边界条件中的多余变量:
简直边的边界条件中无多余变量需要处理;
固支边的边界条件中切应力分量为多余变量,根据物理方程消除:
自由边的边界条件中仅应力分量σxx为多余变量,根据如下物理方程消除:
CN202310109283.2A 2023-02-13 2023-02-13 B样条等几何状态空间有限元方法 Active CN116011301B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310109283.2A CN116011301B (zh) 2023-02-13 2023-02-13 B样条等几何状态空间有限元方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310109283.2A CN116011301B (zh) 2023-02-13 2023-02-13 B样条等几何状态空间有限元方法

Publications (2)

Publication Number Publication Date
CN116011301A CN116011301A (zh) 2023-04-25
CN116011301B true CN116011301B (zh) 2024-01-23

Family

ID=86037371

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310109283.2A Active CN116011301B (zh) 2023-02-13 2023-02-13 B样条等几何状态空间有限元方法

Country Status (1)

Country Link
CN (1) CN116011301B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116757026B (zh) * 2023-06-08 2024-04-05 上海交通大学 一种基于等几何分析的板架结构极限强度分析实现方法
CN116629079B (zh) * 2023-07-21 2024-01-23 北京大学 混合有限元空间构造及求解线弹性力学问题的方法及装置

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112035980A (zh) * 2020-09-08 2020-12-04 南京航空航天大学 等几何混合Kirchhoff-Love壳单元的构造方法
CN114496123A (zh) * 2022-01-25 2022-05-13 北京航空航天大学 一种预测可折叠复合材料圆柱壳轴向压缩性能的方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107526898B (zh) * 2017-09-13 2019-12-27 大连理工大学 变刚度复合材料板壳结构建模分析与可靠度优化设计方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112035980A (zh) * 2020-09-08 2020-12-04 南京航空航天大学 等几何混合Kirchhoff-Love壳单元的构造方法
CN114496123A (zh) * 2022-01-25 2022-05-13 北京航空航天大学 一种预测可折叠复合材料圆柱壳轴向压缩性能的方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
基于状态空间法的叠层材料热力分析;韩志林;程长征;盛宏玉;;重庆大学学报(05);82-89 *
复合材料条形域问题混合状态Hamiltonian元的半解析解;邹贵平,唐立民,周哲玮;应用力学学报(03);1-6 *

Also Published As

Publication number Publication date
CN116011301A (zh) 2023-04-25

Similar Documents

Publication Publication Date Title
CN116011301B (zh) B样条等几何状态空间有限元方法
CN107526898B (zh) 变刚度复合材料板壳结构建模分析与可靠度优化设计方法
Sahoo et al. A new inverse hyperbolic zigzag theory for the static analysis of laminated composite and sandwich plates
JP2022521907A (ja) ハイブリッド繊維複合材料の板巻きシェル構造に対する高速協調最適化方法
CN109241559B (zh) 一种基于子结构的复合材料弹性参数识别方法
Belarbi et al. An efficient eight‐node quadrilateral element for free vibration analysis of multilayer sandwich plates
Zhuang et al. Modal and aeroelastic analysis of trapezoidal corrugated-core sandwich panels in supersonic flow
Sheng et al. A state space finite element for laminated composite plates
Trinh et al. A mixed inverse differential quadrature method for static analysis of constant-and variable-stiffness laminated beams based on Hellinger-Reissner mixed variational formulation
Ferreira et al. Radial basis functions collocation for the bending and free vibration analysis of laminated plates using the Reissner-Mixed Variational Theorem
CN109117504B (zh) 一种双向功能梯度曲壳振动分析方法
CN112487684B (zh) 一种力热耦合环境下层合板非概率可靠性拓扑优化方法
Jha et al. Free vibration of functionally graded plates with a higher-order shear and normal deformation theory
Icardi et al. Optimisation of sandwich panels with functionally graded core and faces
Ton That et al. Nonlinear bending analysis of functionally graded plates using SQ4T elements based on twice interpolation strategy
Dillinger et al. Static aeroelastic stiffness optimization of a forward swept composite wing with CFD-corrected aero loads
Shaterzadeh et al. Stability analysis of composite perforated annular sector plates under thermomechanical loading by finite element method
Shah et al. Through-the-thickness stress distributions near edges of composite laminates using stress recovery scheme and third order shear and normal deformable theory
CN110889253B (zh) 复合材料层合板等效方法
Sheng et al. A semi-analytical finite element for laminated composite plates
Ahmadi Three-dimensional stress analysis in torsion of laminated composite bar with general layer stacking
Helfen et al. Numerical multiscale modelling of sandwich plates
CN117292777A (zh) 一种基于锯齿效应的复合材料层合板层间应力预测方法、装置、介质及设备
CN117352108A (zh) 一种基于d-d铺层的复合材料层合板设计方法
Afshin et al. Free-edge effects in a cylindrical sandwich panel with a flexible core and laminated composite face sheets

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