CN111209657B - 考虑液体表面张力的固体变形界面计算方法 - Google Patents

考虑液体表面张力的固体变形界面计算方法 Download PDF

Info

Publication number
CN111209657B
CN111209657B CN201911397801.5A CN201911397801A CN111209657B CN 111209657 B CN111209657 B CN 111209657B CN 201911397801 A CN201911397801 A CN 201911397801A CN 111209657 B CN111209657 B CN 111209657B
Authority
CN
China
Prior art keywords
equation
solid
load
surface tension
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
CN201911397801.5A
Other languages
English (en)
Other versions
CN111209657A (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.)
Nanjing University of Aeronautics and Astronautics
Original Assignee
Nanjing University of Aeronautics and Astronautics
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 Nanjing University of Aeronautics and Astronautics filed Critical Nanjing University of Aeronautics and Astronautics
Priority to CN201911397801.5A priority Critical patent/CN111209657B/zh
Publication of CN111209657A publication Critical patent/CN111209657A/zh
Application granted granted Critical
Publication of CN111209657B publication Critical patent/CN111209657B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

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

Abstract

本发明提供了一种考虑液体表面张力的固体变形界面计算方法,直接考虑位移‑面力解耦数值方法,利用三次样条曲线对有限元模型中的流‑固边界面的离散点进行拟合,计算出每一点的曲率,得到液‑固界面的整体曲率数据,解决具有表面张力的流‑固耦合载荷问题。本发明的边界曲率计算与有限元分析交替进行,能够充分利用现有有限元程序的前处理、后处理、非线性本构关系等功能。

Description

考虑液体表面张力的固体变形界面计算方法
技术领域
本发明涉及固体变形界面计算领域,具体是一种考虑液体表面张力的固体变形界面计算方法。
背景技术
土壤、混凝土、生物组织属于含液多孔固体介质结构,这种具有流-固耦合属性的结构广泛应用于人们的生产实践中。在流体与固体的分界面,液体表面由于液体分子受力不均匀,促使液体表面收缩,形成表面张力。对于含液多孔固体介质结构,孔尺寸与固体的弹性模量都较小时,液体表面张力的作用对固体变形产生较为明显的影响。另外,液体表面张力对固体产生的面力大小与固体曲率相关。当固体变形时,固体的曲率发生变化,使液体对固体的面力变化,导致固体进一步改变形状,这是一个典型的力-位移耦合问题。
液体对固体的面力与表面张力、液-固界面的曲率有关,而对于同一种液体,其表面张力系数近似为常数。因此,只需要确定液-固界面的曲率就能够确定出液体对固体面力。采用分子动力学理论模拟液滴夹杂问题,能够在微观尺度下进行细致的分析,这种方法受到时间和空间尺度的限制,不适合计算大规模问题,不能够满足原子数目非常多的情况。利用现有的有限元商业软件,只能计算出每一个节点的位移,不能计算网格节点边界的曲率,不能够给出整体边界的曲率,无法计算流体对固体的面力。因此,需要进一步发展有限元计算方法,计算液-固界面的曲率,得到液体对固体面力,从而分析含液多孔结构介质结构的运动状态与力学性能。
发明内容
本发明为了解决现有技术的问题,提供了一种考虑液体表面张力的固体变形界面计算方法,直接考虑位移-面力解耦数值方法,利用三次样条曲线对有限元模型中的流-固边界面的离散点进行拟合,计算出每一点的曲率,得到液-固界面的整体曲率数据,解决具有表面张力的流-固耦合载荷问题。
一种考虑液体表面张力的固体变形界面拟合方法,其特征在于步骤如下:
步骤1:建立流-固耦合几何模型。
步骤2:利用有限元软件划分网格。
步骤3:将外载荷F,分为N步,初始载荷为0,载荷增量△F=F/N。
步骤4:提取流-固边界面节点坐标,设共有n个节点,坐标为(x0,y0)、(x1,y1)、……、(xn,yn),根据公式(1),进行三次样条曲线拟合。
Figure GDA0003424091110000021
其中
Figure GDA0003424091110000022
步骤4.1:在区间[0,1]中,生成n+1个点,为0,1/n,2/n,……,n-1/n,1。
步骤4.2:令t0=0,t1=1/n,t2=2/n,……,tn-1=n-1/n,tn=1。
步骤4.3:将x0=x(t0),x1=x(t1),…,xn-1=x(tn-1),xn=x(tn)和t0,t1,……,tn-1,tn带入公式(2),得到n+1个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程。
Figure GDA0003424091110000023
其中
Figure GDA0003424091110000024
步骤4.4:将公式(2)求导,带入t0与tn,令x’(t0)=0,x’(tn)=0,得到2个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程。
步骤4.5:将步骤4.3与步骤4.4得到的n+3个方程进行求解,得到α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数,得到方程(2)的具体表达式。
步骤4.6:将y0,y1,…,yn-1,yn和t0,t1,……,tn-1,tn带入公式(3),得到n+1个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程。
Figure GDA0003424091110000025
其中
Figure GDA0003424091110000026
步骤4.7:将公式(3)求导,带入t0与tn,令y’(t0)=0,y’(tn)=0,得到2个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程。
步骤4.8:将步骤4.6与步骤4.7得到的n个方程进行求解,得到α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数,得到方程(3)的具体表达式。
步骤4.9:三次样条曲线的参数方程为:
Figure GDA0003424091110000031
步骤4.10:按照公式(5)求解曲线的一阶导数,按照公式(6)求解曲线的二阶导数,按照公式(7)求解曲线的曲率k。
Figure GDA0003424091110000032
Figure GDA0003424091110000033
Figure GDA0003424091110000034
步骤5:按照公式(8),求解各节点处的表面张力对固体产生的面力△P。
△P=kγ (8)
γ为液体表面张力系数。
步骤6:将面力△P与载荷增量的总载荷,施加在有限元模型上,进行分析计算,得到新的流-固边界面坐标。
步骤7:流-固耦合边界面的位置收敛,则外载荷F增加一个增量步;若不收敛,则返回步骤4,进行迭代。
步骤8:步骤7迭代后的载荷大于外载荷F,则有限元计算结束,得出模型最终平衡状态位置;步骤7迭代后的载荷小于或等于外载荷F,返回步骤4,进行迭代。
本发明的有益效果在于:
1.通过三次样条曲线对流-固边界面离散点进行拟合,能够计算出液-固边界面的曲率,确定出液体对表面张力作用对固体的面力。
2.边界曲率计算与有限元分析交替进行,能够充分利用现有有限元程序的前处理、后处理、非线性本构关系等功能。
附图说明
图1是考虑液体表面张力的固体变形界面拟合方法流程图;
图2是加载前有限元网格模型;
图3是表面张力引起的面力与外载荷共同作用后,变形后网格图。
具体实施方式
下面结合附图和具体实施方式对本发明作进一步说明。
实施例1:
步骤1:建立流-固耦合几何模型。
步骤2:利用有限元软件划分网格,如图2。
步骤3:将外载荷F,分为N步,初始载荷为0,载荷增量△F=F/N。
步骤4:提取流-固边界面节点坐标,设共有40个节点,坐标为(x0,y0)、(x1,y1)、……、(x40,y40),根据公式(1),进行三次样条曲线拟合。
Figure GDA0003424091110000041
其中
Figure GDA0003424091110000042
步骤4.1:在区间[0,1]中,生成40+1个点,为0,1/40,2/40,……,39/40,1。
步骤4.2:令t0=0,t1=1/40,t2=2/40,……,t39-1=39/40,t40=1。
步骤4.3:将x0=x(t0),x1=x(t1),…,x39=x(t39),x40=x(t40)和t0,t1,……,t39,t40带入公式(9),得到41个关于α1,α2,α3,α4,β1,β2,……,β38,β39参数的方程。
Figure GDA0003424091110000043
其中
Figure GDA0003424091110000044
步骤4.4:将公式(9)求导,带入t0与t40,令x’(t0)=0,x’(t40)=0,得到2个关于α1,α2,α3,α4,β1,β2,……,β38,β39参数的方程。
步骤4.5:将步骤4.3与步骤4.4得到的43个方程进行求解,得到α1,α2,α3,α4,β1,β2,……,β38,β39参数值,从而得到方程(9)的具体表达式。
步骤4.6:将y0,y1,…,y39,y40和t0,t1,……,t39,t40带入公式(10),得到41个关于α1,α2,α3,α4,β1,β2,……,β39,β40参数的方程。
Figure GDA0003424091110000045
Figure GDA0003424091110000046
步骤4.7:将公式(10)求导,带入t0与t40,令y’(t0)=0,y’(t40)=0,得到2个关于α1,α2,α3,α4,β1,β2,……,β38,β39参数的方程。
步骤4.8:将步骤4.6与步骤4.7得到的43个方程进行求解,得到α1,α2,α3,α4,β1,β2,……,β38,β39参数值,从而得到方程(10)的具体表达式。
步骤4.9:三次样条曲线的参数方程为:
Figure GDA0003424091110000051
步骤4.10:按照公式(5)求解曲线的一阶导数,按照公式(6)求解曲线的二阶导数,按照公式(7)求解曲线的曲率k。
Figure GDA0003424091110000052
Figure GDA0003424091110000053
Figure GDA0003424091110000054
步骤5:按照公式(8),求解各节点处的表面张力对固体产生的面力△P。
△P=kγ (8)
γ为液体表面张力系数。
步骤6:将面力△P与载荷增量的总载荷,施加在有限元模型上,进行分析计算,得到新的流-固边界面坐标。
步骤7:流-固耦合边界面的位置收敛,则外载荷F增加一个增量步;若不收敛,则返回步骤4,进行迭代。
步骤8:步骤7迭代后的载荷大于外载荷F,则有限元计算结束,得出模型最终平衡状态位置;步骤7迭代后的载荷小于或等于外载荷F,返回步骤4,进行迭代。流程图如图1。
本发明具体应用途径很多,以上所述仅是本发明的优选实施方式,应当指出,对于本技术领域的普通技术人员来说,在不脱离本发明原理的前提下,还可以作出若干改进,这些改进也应视为本发明的保护范围。

Claims (2)

1.一种考虑液体表面张力的固体变形界面计算方法,其特征在于:步骤如下:
步骤1:建立流-固耦合几何模型;
步骤2:利用有限元软件划分网格;
步骤3:将外载荷F,分为N步,初始载荷为0,载荷增量△F=F/N;
步骤4:提取流-固边界面节点坐标,设共有n个节点,坐标为(x0,y0)、(x1,y1)、……、(xn,yn),根据以下公式,进行三次样条曲线拟合,并计算节点处曲率k,
Figure FDA0003424091100000011
其中
Figure FDA0003424091100000012
步骤5:求解各节点处的表面张力对固体产生的面力△P,
△P=kγ; (8)
γ为液体表面张力系数;
步骤6:将面力△P与载荷增量的总载荷,施加在有限元模型上,进行分析计算,得到新的流-固边界面坐标;
步骤7:流-固耦合边界面的位置收敛,则外载荷F增加一个增量步;若不收敛,则返回步骤4,进行迭代;
步骤8:步骤7迭代后的载荷大于外载荷F,则有限元计算结束,得出模型最终平衡状态位置;步骤7迭代后的载荷小于或等于外载荷F,返回步骤4,进行迭代。
2.根据权利要求1所述的考虑液体表面张力的固体变形界面计算方法,其特征在于,步骤4所述的曲线拟合及节点处曲率k计算过程如下:
步骤4.1:在区间[0,1]中,生成n+1个点,为0,1/n,2/n,……,n-1/n,1;
步骤4.2:令t0=0,t1=1/n,t2=2/n,……,tn-1=n-1/n,tn=1;
步骤4.3:将x0=x(t0),x1=x(t1),…,xn-1=x(tn-1),xn=x(tn)和t0,t1,……,tn-1,tn带入以下公式,得到n+1个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程,
Figure FDA0003424091100000013
其中
Figure FDA0003424091100000014
步骤4.4:将公式(2)求导,带入t0与tn,令x’(t0)=0,x’(tn)=0,得到2个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程;
步骤4.5:将步骤4.3与步骤4.4得到的n+3个方程进行求解,得到α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数,得到方程(2)的具体表达式;
步骤4.6:将y0,y1,…,yn-1,yn和t0,t1,……,tn-1,tn带入以下公式,得到n+1个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程,
Figure FDA0003424091100000021
其中
Figure FDA0003424091100000022
步骤4.7:将步骤4.6得到的参数方程进行求导,带入t0与tn,令y’(t0)=0,y’(tn)=0,得到2个关于α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数的方程;
步骤4.8:将步骤4.6与步骤4.7得到的n个方程进行求解,得到α1,α2,α3,α4,β1,β2,……,βn-2,βn-1参数,得到4.6)参数方程的具体表达式;
步骤4.9:三次样条曲线的参数方程为:
Figure FDA0003424091100000023
步骤4.10:按照公式(5)求解曲线的一阶导数,按照公式(6)求解曲线的二阶导数,按照公式(7)求解曲线的曲率k,
Figure FDA0003424091100000024
Figure FDA0003424091100000025
Figure FDA0003424091100000026
CN201911397801.5A 2019-12-30 2019-12-30 考虑液体表面张力的固体变形界面计算方法 Active CN111209657B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911397801.5A CN111209657B (zh) 2019-12-30 2019-12-30 考虑液体表面张力的固体变形界面计算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911397801.5A CN111209657B (zh) 2019-12-30 2019-12-30 考虑液体表面张力的固体变形界面计算方法

Publications (2)

Publication Number Publication Date
CN111209657A CN111209657A (zh) 2020-05-29
CN111209657B true CN111209657B (zh) 2022-04-22

Family

ID=70786493

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911397801.5A Active CN111209657B (zh) 2019-12-30 2019-12-30 考虑液体表面张力的固体变形界面计算方法

Country Status (1)

Country Link
CN (1) CN111209657B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114046774B (zh) * 2022-01-05 2022-04-08 中国测绘科学研究院 综合cors网和多源数据的地面形变连续监测方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012122948A (ja) * 2010-12-10 2012-06-28 Jfe Steel Corp 張り剛性評価圧子モデル、その圧子モデルを使用した張り剛性解析装置及び解析方法
KR20120088382A (ko) * 2011-01-31 2012-08-08 한국생산기술연구원 표면장력을 구현하기 위한 유한요소 해석 방법
CN103760089A (zh) * 2014-01-29 2014-04-30 山东农业大学 非饱和土相对渗透系数的试验-数值分析联合测定法
CN104866695A (zh) * 2015-06-24 2015-08-26 武汉大学 一种经gpu加速的浸没边界-格子玻尔兹曼流固耦合模拟方法
CN106599364A (zh) * 2016-11-11 2017-04-26 长春理工大学 贯穿微型腔液滴的有限元流固耦合建模及其模态分析方法
CN107122571A (zh) * 2017-06-06 2017-09-01 大连理工大学 一种考虑水合物分解的沉积物多场耦合模型的建模方法
CN108052783A (zh) * 2018-01-29 2018-05-18 济南大学 一种基于自适应步长的非饱和土动力数值计算方法
CN108984933A (zh) * 2018-07-25 2018-12-11 太原科技大学 弹流润滑条件下计算滚动轴承载荷和压力的边界元法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2009107304A1 (ja) * 2008-02-29 2009-09-03 日本電気株式会社 モデル解析システム、モデル解析方法、及びモデル解析プログラム
US20150242545A1 (en) * 2014-02-21 2015-08-27 Junghyun Cho Method of Simulation of Moving Interfaces using Geometry-Aware Volume of Fluid Method

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012122948A (ja) * 2010-12-10 2012-06-28 Jfe Steel Corp 張り剛性評価圧子モデル、その圧子モデルを使用した張り剛性解析装置及び解析方法
KR20120088382A (ko) * 2011-01-31 2012-08-08 한국생산기술연구원 표면장력을 구현하기 위한 유한요소 해석 방법
CN103760089A (zh) * 2014-01-29 2014-04-30 山东农业大学 非饱和土相对渗透系数的试验-数值分析联合测定法
CN104866695A (zh) * 2015-06-24 2015-08-26 武汉大学 一种经gpu加速的浸没边界-格子玻尔兹曼流固耦合模拟方法
CN106599364A (zh) * 2016-11-11 2017-04-26 长春理工大学 贯穿微型腔液滴的有限元流固耦合建模及其模态分析方法
CN107122571A (zh) * 2017-06-06 2017-09-01 大连理工大学 一种考虑水合物分解的沉积物多场耦合模型的建模方法
CN108052783A (zh) * 2018-01-29 2018-05-18 济南大学 一种基于自适应步长的非饱和土动力数值计算方法
CN108984933A (zh) * 2018-07-25 2018-12-11 太原科技大学 弹流润滑条件下计算滚动轴承载荷和压力的边界元法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
Local limit load solutions for surface cracks in plates and cylinders using finite element analysis;Sattari-Far, I 等;《INTERNATIONAL JOURNAL OF PRESSURE VESSELS AND PIPING》;20040131;第81卷(第1期);第57-66页 *
Modeling for Solder Self-Assembly of MEMS Microstructure Powered by SurfaceTension;Kailin Pan 等;《2006 7th International Conference on Electronic Packaging Technology》;20061231;第1-4页 *
微注射成形中表面张力效应的数值模拟;石建军 等;《应用力学学报》;20101231;第27卷(第4期);第727-731+850页 *
流体黏性及表面张力对气泡运动特性的影响;艾旭鹏 等;《物理学报》;20171231;第66卷(第23期);第233-242页 *

Also Published As

Publication number Publication date
CN111209657A (zh) 2020-05-29

Similar Documents

Publication Publication Date Title
Benson et al. Isogeometric shell analysis: the Reissner–Mindlin shell
CN110298105B (zh) 饱和多孔介质大变形分析的ccpdi-impm方法
Gendre et al. Non-intrusive and exact global/local techniques for structural problems with local plasticity
CN105975645B (zh) 一种基于多步的含激波区域飞行器流场快速计算方法
CN103399992B (zh) 一种基于可靠寿命的结构耐久性优化设计方法
Ii et al. A full Eulerian fluid-membrane coupling method with a smoothed volume-of-fluid approach
Furukawa et al. Accurate cyclic plastic analysis using a neural network material model
CN112528411A (zh) 一种基于模态减缩的几何非线性结构噪声振动响应计算方法
Cui et al. Metal forming analysis using the edge-based smoothed finite element method
CN110705057B (zh) 各向同性固体材料的静态热弹性问题求解方法以及装置
CN111209657B (zh) 考虑液体表面张力的固体变形界面计算方法
JP2000275154A (ja) 応力−ひずみ関係シミュレート方法
CN107633106A (zh) 一种基于全局差分法的非均匀温度场下热模态灵敏度分析方法
KR101029468B1 (ko) 등가정하중을 이용한 동적 비선형 응답 구조 최적해산출방법
CN110188407B (zh) 多孔介质内液体流动参数的确定方法及装置
CN107153755B (zh) 一种页岩气井数值模拟的求解方法
CN113611377B (zh) 一种利用晶体塑性模型模拟混合控制蠕变疲劳变形的方法
CN109271655B (zh) 一种基于非对称有限元算法的材料尺度效应分析方法
CN104834795A (zh) 包带连接结构接触摩擦非线性特性模拟方法及系统
CN107844646A (zh) 一种细长体分布式载荷等效减缩方法
Biotteau et al. Three dimensional automatic refinement method for transient small strain elastoplastic finite element computations
CN116542082A (zh) 一种膜片热压成型时变形预测方法、装置、设备和介质
CN113343523B (zh) 一种基于空间网格的有限元数值模拟分析方法
JP2000312933A (ja) スプリングバック量予測方法
CN110333148B (zh) 一种基于振动衰减曲线精细化分析的土动剪切模量测试方法

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