CN116341419A - 一种流固耦合的数值确定方法、系统及电子设备 - Google Patents

一种流固耦合的数值确定方法、系统及电子设备 Download PDF

Info

Publication number
CN116341419A
CN116341419A CN202310552265.1A CN202310552265A CN116341419A CN 116341419 A CN116341419 A CN 116341419A CN 202310552265 A CN202310552265 A CN 202310552265A CN 116341419 A CN116341419 A CN 116341419A
Authority
CN
China
Prior art keywords
fluid
matrix
solid
unit
grid
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
CN202310552265.1A
Other languages
English (en)
Other versions
CN116341419B (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.)
Institute of Mountain Hazards and Environment IMHE of CAS
Original Assignee
Institute of Mountain Hazards and Environment IMHE of CAS
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 Institute of Mountain Hazards and Environment IMHE of CAS filed Critical Institute of Mountain Hazards and Environment IMHE of CAS
Priority to CN202310552265.1A priority Critical patent/CN116341419B/zh
Publication of CN116341419A publication Critical patent/CN116341419A/zh
Application granted granted Critical
Publication of CN116341419B publication Critical patent/CN116341419B/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
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • 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
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • 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
    • G06F30/25—Design optimisation, verification or simulation using particle-based methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/08—Fluids
    • 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/14—Force analysis or force optimisation, e.g. static or dynamic forces
    • 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)
  • General Engineering & Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Geometry (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • Evolutionary Computation (AREA)
  • Data Mining & Analysis (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Computational Mathematics (AREA)
  • Fluid Mechanics (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开一种流固耦合的数值确定方法、系统及电子设备，涉及流固耦合数值计算领域。本发明对流体域和固体域进行网格划分和参数赋值，能够将流体域与固体域集成为统一区域，采用相同的变量进行描述，同时进行求解，无需传递流‑固边界处的物质信息，提高了流固耦合的计算效率与计算精度，并且克服了“附加质量”的问题。基于网格质量确定是否网格重划分，有效避免了流体大变形所产生的网格畸变问题。通过组合构建系统矩阵，形成了统一构架下的显式算法框架，既满足了显式算法在大型流固系统方程求解时的解耦优势，又避免了流固信息相互传递的数据交换耗时与精度丢失。

Description

一种流固耦合的数值确定方法、系统及电子设备
技术领域
本发明涉及流固耦合数值计算领域，特别是涉及一种流固耦合的数值确定方法、系统及电子设备。
背景技术
流固耦合问题广泛存在于生活中，例如灾害流体与防护结构之间的动态相互作用问题、高速列车运行中的气动弹性问题、管道输送流体以及人体血液输送等。流固耦合问题计算的准确性，直接关乎着流体运动范围判断与结构动态响应的精准性，影响工业设备与工程结构设计的可靠性。流固耦合问题计算的难点在于其伴随的流体大变形、流体自由边界以及流固非定常边界等问题。
流固耦合的数值求解方法按运动描述可分为拉格朗日方法、欧拉方法以及任意欧拉拉格朗日方法三个大类。其中，如何捕捉流体自由边界与流固边界分别是欧拉方法与任意欧拉拉格朗日方法的难点，故而拉格朗日方法对于自由流体边界与流固边界的流固耦合问题更为适合。目前，国内外发展了较多的流固耦合数值求解方法，包括不同的数值算法、不同的积分模型、不同的时间离散格式，但是当前流固耦合数值算法主要集中于发展隐式迭代求解方法的研究，包括不同时间离散模型以及不同的积分算法。同时，也有针对多种不同的数值算法进行耦合求解的研究。某些数值求解方法虽简化了求解过程，但只适用于特定的流固耦合边界条件，总体而言，存在以下问题：
（1）耦合不同求解算法时，在每个载荷步中需要传递两者界面中的状态量，进行交错求解，易引发收敛性问题，同时由于方法间的数据映射、数据交换，降低了流固耦合计算的精度。
（2）隐式求解算法在处理大型的工程问题时需求解大规模的系统方程，严重制约了计算效率；隐式迭代过程中可能存在单元反转而导致负的雅可比矩阵；隐式算法难以发展并行求解策略，制约了大型工程问题流固耦合求解的发展前景。
（3）显式求解方法中，多采用弱耦合的求解策略，将流体域与固体域分开求解，需传递流-固边界处的物质信息，同时可能存在“附加质量”的问题。
基于此，目前极其缺乏流固耦合的统一构架显式计算方法。
发明内容
为解决现有技术存在的上述问题，本发明提供了一种流固耦合的数值确定方法、系统及电子设备，以实现流固耦合数值的统一构架显式计算。
为实现上述目的，本发明提供了如下方案：
一种流固耦合的数值确定方法，包括：
对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元；
确定流体网格单元的单元矩阵和固体网格单元的单元矩阵；所述流体网格单元的单元矩阵包括：单元质量矩阵、单元粘度矩阵、单元压力质量矩阵、流体压力场的单元梯度算子和单元外力向量；所述固体网格单元的单元矩阵包括：单元质量矩阵、单元阻尼矩阵、单元刚度矩阵和单元外力向量；
基于固体域的节点编号和流体域的节点编号确定自由度，并基于所述固体域的节点编号、所述流体域的节点编号和所述自由度组装流体网格单元的单元矩阵和固体网格单元的单元矩阵得到系统矩阵；所述系统矩阵包括：单元质量矩阵、单元粘度/阻尼矩阵、单元梯度算子、单元刚度矩阵、单元压力质量矩阵和单元外力向量；所述固体域的节点编号和流体域的节点编号均在对流体域和固体域进行网格划分的过程中生成；
基于所述系统矩阵确定流固耦合的数值；所述流固耦合的数值包括：中间时刻的速度场、下一时刻的位移场、下一时刻的压力场、下一时刻的加速度场和下一时刻的速度场；
基于所述流固耦合的数值更新所述流体域的节点位置得到更新后的流体网格单元；
确定更新后的流体网格单元的网格质量；
当所述网格质量不满足设定要求时，则对流体域重新进行网格划分和参数赋值得到新的流体网格单元，并返回确定流体网格单元的单元矩阵和固体网格单元的单元矩阵的步骤；
当所述网格质量满足设定要求时，则将上一时刻的流体网格单元替换为更新后的流体网格单元，并返回确定流体网格单元的单元矩阵和固体网格单元的单元矩阵的步骤。
可选地，对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元，具体包括：
采用Delaunay三角网格划分固体得到固体域的三角形单元；
根据固体的几何边界，提取固体域的边界节点，将固体域的边界节点与流体域节点共同进行三角形网格划分，再进行边界识别，删除识别结果中的虚假网格，得到流体域的三角形单元及流体与固体中间的接触单元；
对流体域的三角形单元进行流体参数赋值、对固体域的三角形单元进行固体材料参数赋值，并对流体与固体中间的接触单元赋予流体单元材料性质；所述流体参数包括流体密度、粘度系数和体积模量；固体材料参数包括密度、弹性模量和泊松比。
可选地，采用alpha-shape方法识别所述固体的几何边界。
可选地，基于所述系统矩阵确定流固耦合的数值，具体包括：
基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定中间时刻的速度场，并根据当前时刻的位移场、预设时步和所述中间时刻的速度场确定下一时刻的位移场；
基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定下一时刻的加速度场，并根据所述下一时刻的加速度场和所述中间时刻的速度场确定下一时刻的速度场。
可选地，基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定中间时刻的速度场，具体包括：
确定临界时步，并基于所述临界时步确定预设时步；
基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定初始加速度场；
基于初始加速度和所述预设时步确定上一时刻的加速度场；
采用中心差分的时间离散格式，根据上一时刻的速度场和当前时刻的加速度场确定中间时刻的速度场。
可选地，对于不可压缩流体，基于中间时刻的速度场以及所述系统矩阵中的单元压力质量矩阵、单元梯度算子确定下一时刻的压力场；
对于可压缩流体，根据流体的密度变化确定下一时刻的压力场。
可选地，所述设定要求为流体域的三角形单元中三角形的角度小于设定值。
根据本发明提供的具体实施例，本发明公开了以下技术效果：
本发明通过对流体域和固体域进行网格划分和参数赋值，能够将流体域与固体域集成为统一区域，采用相同的变量进行描述、进行同时求解，无需传递流-固边界处的物质信息，提高了流固耦合的计算效率与计算精度，并且克服了“附加质量”的问题。本发明基于网格质量确定是否网格重划分，有效地避免了流体大变形所产生的网格畸变问题。本发明通过组合构建系统矩阵，形成了统一构架下的显式算法框架，既满足了显式算法在大型流固系统方程求解时的解耦优势，又避免了流固信息相互传递的数据交换耗时与精度丢失。此外，本发明在确定流固耦合的数值时，压力场、位移场、速度场以及加速度场都可以单独求解，极易扩展并行算法，为将来工程中大型流固耦合问题的快速求解提供了良好的算法构架基础。
本发明还提供了以下两种实施结构：
一种流固耦合的数值确定系统，应用于上述提供的流固耦合的数值确定方法；所述系统包括：
网格单元划分模块，用于对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元；
单元矩阵确定模块，用于确定流体网格单元的单元矩阵和固体网格单元的单元矩阵；所述流体网格单元的单元矩阵包括：单元质量矩阵、单元粘度矩阵、单元压力质量矩阵、流体压力场的单元梯度算子和单元外力向量；所述固体网格单元的单元矩阵包括：单元质量矩阵、单元阻尼矩阵、单元刚度矩阵和单元外力向量；
系统矩阵确定模块，用于基于固体域的节点编号和流体域的节点编号确定自由度，并基于所述固体域的节点编号、所述流体域的节点编号和所述自由度组装流体网格单元的单元矩阵和固体网格单元的单元矩阵得到系统矩阵；所述系统矩阵包括：单元质量矩阵、单元粘度/阻尼矩阵、单元梯度算子、单元刚度矩阵、单元压力质量矩阵和单元外力向量；所述固体域的节点编号和流体域的节点编号均在对流体域和固体域进行网格划分的过程中生成；
耦合数值确定模块，用于基于所述系统矩阵确定流固耦合的数值；所述流固耦合的数值包括：中间时刻的速度场、下一时刻的位移场、下一时刻的压力场、下一时刻的加速度场和下一时刻的速度场；
网格单元更新模块，用于基于所述流固耦合的数值更新所述流体域的节点位置得到更新后的流体网格单元；
网格质量确定模块，用于确定更新后的流体网格单元的网格质量；
第一循环执行模块，用于当所述网格质量不满足设定要求时，则对流体域重新进行网格划分和参数赋值得到新的流体网格单元，并返回所述单元矩阵确定模块中执行的步骤；
第二循环执行模块，用于当所述网格质量满足设定要求时，则将上一时刻的流体网格单元替换为更新后的流体网格单元，并返回所述单元矩阵确定模块中执行的步骤。
一种电子设备，包括：
存储器，用于存储计算机程序；
处理器，与所述存储器连接，用于调取并执行所述计算机程序，以实施上述提供的流固耦合的数值确定方法。
可选地，所述存储器为计算机可读存储介质。
因本发明上述提供的两种实施结构实现的技术效果与本发明提供的流固耦合的数值确定方法实现的技术效果相同，故在此不再进行赘述。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明提供的流固耦合的数值确定方法的实施流程图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明的目的是提供一种流固耦合的数值确定方法、系统及电子设备，以实现流固耦合数值的统一构架显式计算，进而能够提高流固耦合的计算效率与计算精度，并且克服“附加质量”的问题，能够有效避免流体大变形所产生的网格畸变问题，同时，既能够满足显式算法在大型流固系统方程求解时的解耦优势，又能够避免出现流固信息相互传递的数据交换耗时与精度丢失的问题。
为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。
如图1所示，本发明提供的流固耦合的数值确定方法，包括：
步骤1：流体域及固体域网格划分及参数赋值。在该步骤中，将流体域及固体域划分为三角形单元，流体域与固体域中间的网格通过边界识别剔除虚假的单元，形成初始网格。获取流体域及固体域的材料参数及边界条件，分别对流体域及固体域的单元赋值，并施加相应的力与位移边界条件。具体的：
步骤1-1：固体域内部无需进行网格重划分，故只需进行一次Delaunay三角网格划分，得到固体域的三角形单元。
步骤1-2：根据固体几何边界，提取固体域边界节点，将固体边界节点与流体域节点共同进行三角形网格划分，再进行边界识别，删除其中的虚假网格，得到流体域的三角形单元及流体与固体中间的接触单元。
其中，边界识别采用alpha-shape方法，对三角形单元的外接圆半径R进行判断：
Figure SMS_1
。
式中，h为三角形单元的平均外接圆半径长度，
Figure SMS_2
为设定的参数，取值范围一般为1</>
Figure SMS_3
<2。
获取参数：流体参数包括流体密度
Figure SMS_4
、粘度系数/>
Figure SMS_5
和体积模量/>
Figure SMS_6
，固体材料参数包括密度/>
Figure SMS_7
、弹性模量/>
Figure SMS_8
和泊松比/>
Figure SMS_9
，固体域与流体域中间的接触单元赋予流体单元材料性质。
基于上述步骤得到流体网格单元与固体网格单元。
步骤2：计算每个流体网格单元与固体网格单元的单元矩阵并进行组装。该步骤中，主要是计算流体域的单元质量矩阵
Figure SMS_10
、单元粘度矩阵/>
Figure SMS_13
、单元压力质量矩阵/>
Figure SMS_16
、流体压力场的单元梯度算子/>
Figure SMS_12
和单元外力向量/>
Figure SMS_14
。计算固体域的单元质量矩阵/>
Figure SMS_17
、单元阻尼矩阵/>
Figure SMS_18
、单元刚度矩阵/>
Figure SMS_11
和单元外力向量/>
Figure SMS_15
。
为形成统一构架，根据固体域与流体域的节点编号及自由度将流体域与固体域的单元矩阵及向量进行共同组装，形成总体的系统矩阵。系统矩阵包括单元质量矩阵
Figure SMS_19
、单元粘度/阻尼矩阵/>
Figure SMS_20
、单元梯度算子/>
Figure SMS_21
、单元刚度矩阵/>
Figure SMS_22
、单元外力向量/>
Figure SMS_23
和单元压力质量矩阵/>
Figure SMS_24
。
步骤2-1：流体与固体的单元质量矩阵表达式为：
Figure SMS_25
。
式中，
Figure SMS_26
、/>
Figure SMS_27
、/>
Figure SMS_28
、/>
Figure SMS_29
和/>
Figure SMS_30
分别为密度、形函数、面积、厚度与速度。需要注意的是，对于可压缩的流体单元，其密度是变化的，故而求解单元质量矩阵时，需采用实时更新的节点密度值，流体网格单元内的密度/>
Figure SMS_31
也采用插值的形式，为：
Figure SMS_32
。
式中，
Figure SMS_33
为密度插值所采用的形函数，/>
Figure SMS_34
与/>
Figure SMS_35
分别为节点i的密度及对应节点i的形函数，/>
Figure SMS_36
为单元密度向量，/>
Figure SMS_37
和/>
Figure SMS_38
分别为节点i和节点k的形函数。
步骤2-2：结合流体域的单元质量矩阵及流体域的单元内密度插值，将可压缩流体的单元质量矩阵
Figure SMS_39
变为：
Figure SMS_40
步骤2-3：流体单元粘度矩阵
Figure SMS_41
为：
Figure SMS_42
。
式中，
Figure SMS_43
与/>
Figure SMS_44
分别为应变矩阵与流体本构矩阵，/>
Figure SMS_45
表达式如下：
Figure SMS_46
。
式中，
Figure SMS_47
为流体粘度系数。
步骤2-4：流体域的单元压力质量矩阵
Figure SMS_48
通过对角化处理后变为：
Figure SMS_49
。
式中，
Figure SMS_50
为体积模量。
步骤2-5：流体域的压力场的单元梯度算子
Figure SMS_51
为：
Figure SMS_52
。
步骤2-6：固体单元刚度矩阵
Figure SMS_53
为：
Figure SMS_54
。
固体单元根据固体材料性质选择相应的阻尼矩阵，常用的瑞利阻尼可采用质量矩阵与刚度矩阵的线性组合。固体本构矩阵
Figure SMS_55
为：
Figure SMS_56
。
式中，拉梅常数
Figure SMS_57
与/>
Figure SMS_58
为：
Figure SMS_59
。
步骤2-7：再根据节点编号及自由度确定单元矩阵中的元素在总体矩阵中的位置，对每个单元进行循环，进行叠加即可组装为总体系统矩阵。
步骤3：计算中间时刻的速度场与下一时刻的位移场。
步骤3-1：确定临界时步
显式算法是条件稳定的，所选取的时步
Figure SMS_60
需小于临界值，否则显式算法不稳定，临界时步/>
Figure SMS_61
的确定公式为：
Figure SMS_62
。
式中，
Figure SMS_63
为安全系数，/>
Figure SMS_64
为单元特征尺寸，/>
Figure SMS_65
为声音在单元介质中的传播速度。
步骤3-2：计算初始加速度
Figure SMS_66
：
第一步计算时，初始加速度
Figure SMS_67
可看作系统的初始内力与外力引起系统不平衡的加速度趋势，计算公式如下：
Figure SMS_68
。
其中，
Figure SMS_69
、/>
Figure SMS_70
和/>
Figure SMS_71
分别表示初始速度、初始位移和初始压力。
步骤3-3：采用中心差分的时间离散格式，中间时刻的速度场通过上一时刻的速度场和加速度场进行计算。中间时刻的速度场
Figure SMS_72
为：
Figure SMS_73
。式中，/>
Figure SMS_74
和/>
Figure SMS_75
分别为上一时刻的速度场和加速度场。
步骤3-4：下一时刻的位移场通过上一时刻的位移场与中间时刻的速度场进行计算。下一时刻的位移场
Figure SMS_76
为：
Figure SMS_77
。/>
Figure SMS_78
为上一时刻的位移场。
步骤4：计算下一时刻的压力场，下一时刻的压力场
Figure SMS_79
。
下一时刻的压力场需根据可压缩流体与不可压缩流体两种进行计算，具体如下：
a)对于不可压缩流体，根据中间时刻的速度场
Figure SMS_80
进行计算。中间时刻的质量守恒方程为：
Figure SMS_81
。
式中，
Figure SMS_82
为中间时刻的压力场，根据中心差分离散格式，有：
Figure SMS_83
。
式中，
Figure SMS_84
为下一时刻的压力场，/>
Figure SMS_85
为上一时刻的压力场。
结合上述两式即可得下一时刻的压力场
Figure SMS_86
为：
Figure SMS_87
。
b)对于可压缩流体，根据流体的密度变化进行计算。质量守恒方程可简化为：
Figure SMS_88
。
式中，
Figure SMS_89
为当前时步节点密度插值矩阵的集中形式，/>
Figure SMS_90
为当前时步节点密度向量，/>
Figure SMS_91
为节点等效密度的向量形式，/>
Figure SMS_92
为上一时步节点密度插值矩阵的集中形式，/>
Figure SMS_93
为上一时步节点密度向量。
故而，下一时刻的节点压力场
Figure SMS_94
可根据节点的密度变化，通过状态方程求出，如下式所示：
Figure SMS_95
。
式中，
Figure SMS_96
、/>
Figure SMS_97
与/>
Figure SMS_98
分别为上一时刻的压力值、体积模量与密度，系数/>
Figure SMS_99
取值一般为7。
步骤5：计算下一时刻的加速度场，并更新下一时刻的速度场。
此处与步骤3都用到了流固耦合系统质量矩阵的逆矩阵，故需采用集中质量矩阵，不然大型矩阵求逆将极为困难。根据流固耦合系统的平衡方程
Figure SMS_100
，可推出下一时刻加速度场/>
Figure SMS_101
的计算公式为：
Figure SMS_102
。
式中，采用的速度项为中间时刻的速度向量，故计算完下一时刻的速度场
Figure SMS_103
后，可采用下式对加速度场进行修正：
Figure SMS_104
。
下一时刻的速度场可通过中间时刻的速度场与下一时刻的加速度场计算，为：
Figure SMS_105
。
步骤6：判断所有流体单元的网格质量，判断是否需要进行网格重划分。
若发生单元反转或单元畸变达到临界值，则删除所有流体网格，并对流体域的节点进行网格重划分，再进入计算步骤2，进行下一时刻计算。
具体的，前面步骤计算完成后，更新所有网格单元的节点位置，对流体网格单元进行循环，判断每个流体网格单元的畸形度，可根据每个流体三角形单元的角度进行判断，当存在三角形角度小于设定值时（此时即为畸形度超出预期），则进入步骤1的网格重新划分，网格划分完成后则进入步骤2进行下一时刻的计算。当所有单元都满足要求时，则无需重划分网格，直接进入步骤2进行下一时刻的计算。
基于上述描述，相对于现有技术，本发明具有以下优点：
（1）本发明采用拉格朗日框架描述，可有效识别流体的自由界面与非定常流固边界，且将流体域与固体域集成为统一区域，采用相同的变量进行描述，进行同时求解，无需传递流-固边界处的物质信息，提高了流固耦合的计算效率与计算精度，并且克服了“附加质量”的问题。
（2）本发明采用基于拉格朗日方法的粒子有限元方法，结合了有限元方法与网格重划分技术，有效地避免了流体大变形所产生的网格畸变问题，且所有信息可均存储在节点位置上。
（3）本发明采用了统一构架下的显式算法框架，既满足了显式算法在大型流固系统方程求解时的解耦优势，又避免了流固信息相互传递的数据交换耗时与精度丢失。
（4）本发明的压力场、位移场、速度场以及加速度场都可以单个方程进行求解，在本发明的基础上，极易扩展并行算法，为将来工程中大型流固耦合问题的快速求解提供了良好的算法构架基础。
本发明还提供了以下两种实施结构：
一种流固耦合的数值确定系统，应用于上述提供的流固耦合的数值确定方法。系统包括：
网格单元划分模块，用于对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元。
单元矩阵确定模块，用于确定流体网格单元的单元矩阵和固体网格单元的单元矩阵。流体网格单元的单元矩阵包括：单元质量矩阵、单元粘度矩阵、单元压力质量矩阵、流体压力场的单元梯度算子和单元外力向量。固体网格单元的单元矩阵均包括：单元质量矩阵、单元阻尼矩阵、单元刚度矩阵和单元外力向量。
系统矩阵确定模块，用于基于固体域的节点编号和流体域的节点编号确定自由度，并基于固体域的节点编号、流体域的节点编号和自由度组装流体网格单元的单元矩阵和固体网格单元的单元矩阵得到系统矩阵。系统矩阵包括：单元质量矩阵、单元粘度/阻尼矩阵、单元梯度算子、单元刚度矩阵、单元压力质量矩阵和单元外力向量。固体域的节点编号和流体域的节点编号均在对流体域和固体域进行网格划分的过程中生成。
耦合数值确定模块，用于基于系统矩阵确定流固耦合的数值。流固耦合的数值包括：中间时刻的速度场、下一时刻的位移场、下一时刻的压力场、下一时刻的加速度场和下一时刻的速度场。
网格单元更新模块，用于基于流固耦合的数值更新流体域的节点位置得到更新后的流体网格单元。
网格质量确定模块，用于确定更新后的流体网格单元的网格质量。
第一循环执行模块，用于当网格质量不满足设定要求时，则对流体域重新进行网格划分和参数赋值得到新的流体网格单元，并返回单元矩阵确定模块中执行的步骤。
第二循环执行模块，用于当网格质量满足设定要求时，则将上一时刻的流体网格单元替换为更新后的流体网格单元，并返回单元矩阵确定模块中执行的步骤。
一种电子设备，包括：
存储器，用于存储计算机程序。
处理器，与存储器连接，用于调取并执行计算机程序，以实施上述提供的流固耦合的数值确定方法。
此外，上述的存储器中的计算机程序通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备（可以是个人计算机、服务器或者网络设备等）执行本发明各个实施例方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器、随机存取存储器、磁碟或者光盘等各种可以存储程序代码的介质。
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。对于实施例公开的系统和电子设备而言，由于其与实施例公开的方法相对应，所以描述的比较简单，相关之处参见方法部分说明即可。
本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处。综上所述，本说明书内容不应理解为对本发明的限制。

Claims (10)

1.一种流固耦合的数值确定方法，其特征在于，包括：
对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元；
确定流体网格单元的单元矩阵和固体网格单元的单元矩阵；所述流体网格单元的单元矩阵包括：单元质量矩阵、单元粘度矩阵、单元压力质量矩阵、流体压力场的单元梯度算子和单元外力向量；所述固体网格单元的单元矩阵包括：单元质量矩阵、单元阻尼矩阵、单元刚度矩阵和单元外力向量；
基于固体域的节点编号和流体域的节点编号确定自由度，并基于所述固体域的节点编号、所述流体域的节点编号和所述自由度组装流体网格单元的单元矩阵和固体网格单元的单元矩阵得到系统矩阵；所述系统矩阵包括：单元质量矩阵、单元粘度/阻尼矩阵、单元梯度算子、单元刚度矩阵、单元压力质量矩阵和单元外力向量；所述固体域的节点编号和流体域的节点编号均在对流体域和固体域进行网格划分的过程中生成；
基于所述系统矩阵确定流固耦合的数值；所述流固耦合的数值包括：中间时刻的速度场、下一时刻的位移场、下一时刻的压力场、下一时刻的加速度场和下一时刻的速度场；
基于所述流固耦合的数值更新所述流体域的节点位置得到更新后的流体网格单元；
确定更新后的流体网格单元的网格质量；
当所述网格质量不满足设定要求时，则对流体域重新进行网格划分和参数赋值得到新的流体网格单元，并返回确定流体网格单元的单元矩阵和固体网格单元的单元矩阵的步骤；
当所述网格质量满足设定要求时，则将上一时刻的流体网格单元替换为更新后的流体网格单元，并返回确定流体网格单元的单元矩阵和固体网格单元的单元矩阵的步骤。
2.根据权利要求1所述的流固耦合的数值确定方法，其特征在于，对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元，具体包括：
采用Delaunay三角网格划分固体得到固体域的三角形单元；
根据固体的几何边界，提取固体域的边界节点，将固体域的边界节点与流体域节点共同进行三角形网格划分，再进行边界识别，删除识别结果中的虚假网格，得到流体域的三角形单元及流体与固体中间的接触单元；
对流体域的三角形单元进行流体参数赋值、对固体域的三角形单元进行固体材料参数赋值，并对流体与固体中间的接触单元赋予流体单元材料性质；所述流体参数包括流体密度、粘度系数和体积模量；固体材料参数包括密度、弹性模量和泊松比。
3.根据权利要求2所述的流固耦合的数值确定方法，其特征在于，采用alpha-shape方法识别所述固体的几何边界。
4.根据权利要求1所述的流固耦合的数值确定方法，其特征在于，基于所述系统矩阵确定流固耦合的数值，具体包括：
基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定中间时刻的速度场，并根据当前时刻的位移场、预设时步和所述中间时刻的速度场确定下一时刻的位移场；
基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定下一时刻的加速度场，并根据所述下一时刻的加速度场和所述中间时刻的速度场确定下一时刻的速度场。
5.根据权利要求4所述的流固耦合的数值确定方法，其特征在于，基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定中间时刻的速度场，具体包括：
确定临界时步，并基于所述临界时步确定预设时步；
基于所述系统矩阵中的单元质量矩阵、单元外力向量、单元粘度/阻尼矩阵、单元刚度矩阵和单元梯度算子确定初始加速度场；
基于初始加速度和所述预设时步确定上一时刻的加速度场；
采用中心差分的时间离散格式，根据上一时刻的速度场和当前时刻的加速度场确定中间时刻的速度场。
6.根据权利要求4所述的流固耦合的数值确定方法，其特征在于，对于不可压缩流体，基于中间时刻的速度场以及所述系统矩阵中的单元压力质量矩阵、单元梯度算子确定下一时刻的压力场；
对于可压缩流体，根据流体的密度变化确定下一时刻的压力场。
7.根据权利要求2所述的流固耦合的数值确定方法，其特征在于，所述设定要求为流体域的三角形单元中三角形的角度小于设定值。
8.一种流固耦合的数值确定系统，其特征在于，应用于如权利要求1-7任意一项所述的流固耦合的数值确定方法；所述系统包括：
网格单元划分模块，用于对流体域和固体域进行网格划分和参数赋值得到流体网格单元和固体网格单元；
单元矩阵确定模块，用于确定流体网格单元的单元矩阵和固体网格单元的单元矩阵；所述流体网格单元的单元矩阵包括：单元质量矩阵、单元粘度矩阵、单元压力质量矩阵、流体压力场的单元梯度算子和单元外力向量；所述固体网格单元的单元矩阵包括：单元质量矩阵、单元阻尼矩阵、单元刚度矩阵和单元外力向量；
系统矩阵确定模块，用于基于固体域的节点编号和流体域的节点编号确定自由度，并基于所述固体域的节点编号、所述流体域的节点编号和所述自由度组装流体网格单元的单元矩阵和固体网格单元的单元矩阵得到系统矩阵；所述系统矩阵包括：单元质量矩阵、单元粘度/阻尼矩阵、单元梯度算子、单元刚度矩阵、单元压力质量矩阵和单元外力向量；所述固体域的节点编号和流体域的节点编号均在对流体域和固体域进行网格划分的过程中生成；
耦合数值确定模块，用于基于所述系统矩阵确定流固耦合的数值；所述流固耦合的数值包括：中间时刻的速度场、下一时刻的位移场、下一时刻的压力场、下一时刻的加速度场和下一时刻的速度场；
网格单元更新模块，用于基于所述流固耦合的数值更新所述流体域的节点位置得到更新后的流体网格单元；
网格质量确定模块，用于确定更新后的流体网格单元的网格质量；
第一循环执行模块，用于当所述网格质量不满足设定要求时，则对流体域重新进行网格划分和参数赋值得到新的流体网格单元，并返回所述单元矩阵确定模块中执行的步骤；
第二循环执行模块，用于当所述网格质量满足设定要求时，则将上一时刻的流体网格单元替换为更新后的流体网格单元，并返回所述单元矩阵确定模块中执行的步骤。
9.一种电子设备，其特征在于，包括：
存储器，用于存储计算机程序；
处理器，与所述存储器连接，用于调取并执行所述计算机程序，以实施如权利要求1-7任意一项所述的流固耦合的数值确定方法。
10.根据权利要求9所述的电子设备，其特征在于，所述存储器为计算机可读存储介质。
CN202310552265.1A 2023-05-17 2023-05-17 一种流固耦合的数值确定方法、系统及电子设备 Active CN116341419B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202310552265.1A CN116341419B (zh) 2023-05-17 2023-05-17 一种流固耦合的数值确定方法、系统及电子设备

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202310552265.1A CN116341419B (zh) 2023-05-17 2023-05-17 一种流固耦合的数值确定方法、系统及电子设备

Publications (2)

Publication Number Publication Date
CN116341419A true CN116341419A (zh) 2023-06-27
CN116341419B CN116341419B (zh) 2023-08-01

Family

ID=86895185

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202310552265.1A Active CN116341419B (zh) 2023-05-17 2023-05-17 一种流固耦合的数值确定方法、系统及电子设备

Country Status (1)

Country Link
CN (1) CN116341419B (zh)

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20140257765A1 (en) * 2013-03-05 2014-09-11 Livermore Software Technology Corporation Numerical Simulation of FSI Using The Space-Time CE/SE Solver With A Moving Mesh For The Fluid Domain
CN104850689A (zh) * 2015-04-30 2015-08-19 昆明理工大学 一种基于固定网格技术的流固耦合计算方法
CN106383930A (zh) * 2016-08-31 2017-02-08 中国石油大学(华东) 一种尾轴承‑转子系统多重流固耦合计算方法
CN107506562A (zh) * 2017-09-29 2017-12-22 西安科技大学 一种水润滑橡胶轴承双向热流固耦合计算方法
CN113378440A (zh) * 2021-06-23 2021-09-10 四川大学 一种流固耦合数值模拟计算方法、装置及设备
CN114912333A (zh) * 2022-07-19 2022-08-16 上海索辰信息科技股份有限公司 一种轻流体下和结构的声振耦合的模拟方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20140257765A1 (en) * 2013-03-05 2014-09-11 Livermore Software Technology Corporation Numerical Simulation of FSI Using The Space-Time CE/SE Solver With A Moving Mesh For The Fluid Domain
CN104850689A (zh) * 2015-04-30 2015-08-19 昆明理工大学 一种基于固定网格技术的流固耦合计算方法
CN106383930A (zh) * 2016-08-31 2017-02-08 中国石油大学(华东) 一种尾轴承‑转子系统多重流固耦合计算方法
CN107506562A (zh) * 2017-09-29 2017-12-22 西安科技大学 一种水润滑橡胶轴承双向热流固耦合计算方法
CN113378440A (zh) * 2021-06-23 2021-09-10 四川大学 一种流固耦合数值模拟计算方法、装置及设备
CN114912333A (zh) * 2022-07-19 2022-08-16 上海索辰信息科技股份有限公司 一种轻流体下和结构的声振耦合的模拟方法

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
WOLFGANG A. WALL 等: "A strong coupling partitioned approach for fluid–structure interaction with free surfaces", COMPUTERS & FLUIDS, vol. 36, pages 169, XP005720663, DOI: 10.1016/j.compfluid.2005.08.007 *
崔涛 等: "高速列车流固耦合振动的研究方法及其应用", 铁道学报, vol. 35, no. 04, pages 16 - 22 *
张强 等: "圆筒内旋转细长管与流体界面耦合的数值方法", 计算力学学报, vol. 28, no. 04, pages 601 - 606 *
张胜涛 等: "高超声速进气道前缘流场-热-结构耦合分析", 空气动力学学报, vol. 35, no. 03, pages 436 - 443 *
李娜 等: "基于CBS有限元的流动-传热-变形耦合计算方法", 南京航空航天大学学报, no. 05, pages 622 - 626 *
苏丹 等: "基于CFD技术的高效叶栅耦合颤振分析方法", 航空学报, vol. 35, no. 12, pages 3232 - 3243 *

Also Published As

Publication number Publication date
CN116341419B (zh) 2023-08-01

Similar Documents

Publication Publication Date Title
Long et al. A meshless local Petrov-Galerkin method for solving the bending problem of a thin plate
Dobrev et al. The target-matrix optimization paradigm for high-order meshes
CN111859766B (zh) 可变计算域的拉格朗日积分点有限元数值仿真系统及方法
CN116245049B (zh) 节点式非结构网格的边界修正方法、装置、设备及介质
Li et al. An adaptive SVD–Krylov reduced order model for surrogate based structural shape optimization through isogeometric boundary element method
Pappalardo et al. Development of ANCF tetrahedral finite elements for the nonlinear dynamics of flexible structures
Wang et al. Semi-implicit formulation of the immersed finite element method
CN109740182A (zh) 一种基于再生核粒子的无网格物理变形仿真方法
CN109977475A (zh) 一种用于三维复杂管路流固耦合计算的动网格更新方法
Shontz et al. A robust solution procedure for hyperelastic solids with large boundary deformation
Liu et al. A topology optimization method with constant volume fraction during iterations for design of compliant mechanisms
CN114781610A (zh) 一种数据的处理方法、神经网络的训练方法以及相关设备
CN111125963A (zh) 基于拉格朗日积分点有限元的数值仿真系统及方法
Zheng et al. An improved local remeshing algorithm for moving boundary problems
CN106354954A (zh) 一种基于叠层基函数的三维力学模态仿真模拟方法
JP5589198B2 (ja) 有限要素法において８ノード六面体エレメントの剪断ロッキングを低減する方法
CN102819454A (zh) 基于gpu的有限元显式并行求解仿真方法
CN116341419B (zh) 一种流固耦合的数值确定方法、系统及电子设备
CN111881629B (zh) 一种气动热-结构热传导耦合非线性降阶模型方法
CN113378440A (zh) 一种流固耦合数值模拟计算方法、装置及设备
JP6043146B2 (ja) 骨に沿った筋肉の運動および関節のまわりの筋肉の運動を数値的にシミュレートする方法およびシステム
CN110147571B (zh) 一种组件结构的拓扑优化方法及装置
McDaniel et al. Efficient mesh deformation for computational stability and control analyses on unstructured viscous meshes
Sheng et al. Multiblock approach for calculating incompressible fluid flows on unstructured grids
CN112396154A (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