《模拟维粘性液体非线性大幅晃动的ale有限元分步法》由会员分享,可在线阅读,更多相关《模拟维粘性液体非线性大幅晃动的ale有限元分步法(8页珍藏版)》请在金锄头文库上搜索。
1、衡以 三 b ibloikw t 占 性 液 体 非 线 性 大 幅 晃 动 的 A L E 有 限 元 分 步 法 岳宝增王照林 ( 清华大学工程力学系北京, 1 0 0 0 8 4 刘延柱 上海交通大学工程力学系,上海,2 0 0 0 3 0 摘 要针叶N a v i e r - S t o k e s 方程将人 L E ( A r b i tr a ry L a g a n g e - E i d e r ) 方法引 入 到有限元的速度压力 分步 格式中 , 利用G a le r k i n 加权余童推导了 求解三维液体非 线性大幅羌动的数值离散方程. 推导了三维液体自由液面上结点法向矢
2、t的教值 计算方法 通过叶国 筒形贮腔申 三维液体大 幅芜 动的教值 模拟、 揭示了 三维液体 火幅芜动的一些重 要非线 性特性. 通过对带有团 环形肋板圆筒 形贮腔中 三维液 体 大幅芜 动的 教值模拟,揭示了 圈环形隔板队 尼器 时免动的阻 尼机理和阻尼效应. 编制了 三维流 体非线性大幅芜动的教值模拟计算机程序 进行了算例 计算分析, 得出了一些重要结论,并与有关实验结果基本符合 关健词犬幅芜 动,N a v i e r - S to k “方 程, 环形 肋板, 有限 元分步 法、 数值 模扣 带有自由液面的非稳态不可压粘性液体晃动间题的求 解, 具有 广泛的 工程背 景, 特别 是
3、涉及 航 天 , 化 工, 贮 运 等 领 域, 已 有 大量 的 文 献 报 导 0 - 1 s l . 求 解 晃 动问 题 的 主 要困 难在 于 : 一方面荃本变量速度和压力的求解与自 由液面的位置 有关, 另 一方面自由 液面 的位置 又随着 时间的变化而变动, 而且自由 液面上的边界条件是由 非常 复杂 的非 线性方程 所描述的。 近年 来, 王 照 林 老 师 在 这 方面 进 行了 卓 有 成 效 的 研 究, 并 取 得了 一 系 列 重 要 的 成 果 r2 3 .1 5 . 1 6 1 . 与 传统的差分方法相比 较、 有限元 方法的最大突出优点是能比 较容易地处理各种复
4、杂的几何 形状和统一处理各种典型的边界条件且精度较高.利用有限元 方法对于三维问题的研究, 由 于 数 据 处 理 及 计 算 机 程 序 实 现 方 面 的 困 难 , 目 前 这 方 面 的 大 献 报 导 还 比 较 少 用 有 限 元 法求解自 由面粘性流动 时, 运动的流体与有限元网络之间的关系即运动学描述是极为重要 的.由于拉格朗日 方法和欧拉方法这两种经典的描述方法各有长短, 这就促使人们把这两 种方法结合起来使用.A L E ( A r b i tr a ry L a g r a n g e -E i d e r ) 方法正 是由 此而 逐渐发展完 善起来 并 应 用 于 有
5、限 元 分 析 中 “ ,III 。 一 般 来 说 , 利 用 有 限 元 方 法 求 解 非 稳 态N a v ie r- S to k e , 方 程 , 速 度 多 项 式 插 值函 数 必 须比 压 力 的 多 项 式 插 值函 数 高 一 阶 即 混 合 插 值 法 Is l u l , 这 不 仅 大 大 增 加了 公 式 的 复杂 程 度 , 也 使 计 算 机 程 序 实 现 及 数 据 输人 更 加 困 难 。 文 献 1 4 推 导 的 有限 元分步格式使得这些困 难 得到克 服, 本 文将A L E 描述引 人到这种 分步 格式中 从而 能更有效地 国 家自 然 科 学
6、 基 金( 重 点 项目 批 准 号 : 1 9 3 3 2 0 2 0 ) 和 高 等 学 校 博 士学 科 点专 项 科 研 基 金 资 助 项目 1 6 6 模拟大幅晃动问题. 分步法 ( F r a c t i o n a l S t e p M e t h o d) 最 早是1h A . J . C h o ri n f , 荃本 思想是 把一个时 间增量步分为 两 步或更多 步, 第 在有限差分法中提出的,其 众 后 求 解一 个 近 似 的 中 间 速 度 场 , 一 步在动 量方程中 略去 压力项 或压力增 量 它一 般不满足 连 续方程; 第 二步, 由中 间速 度场求出 相
7、应的 压力 场, 而 速度由 质 量守恒来修 正, 这 一步导出 压力 点是对时间增量分步, P o i s s n方 程. 这 种分步有 限元方法的特 然后建 立相 应的有限元方程.在提出适当 的压力梯度边界条件后 这种方法的速度和压力可用同阶线性插值函数. 尤 其是 采用集中质量 法后, 这就使得有限元方程在算法上结构简单, 可 以 得 到 求 解 速 度 的 显 式 格 式 . 分 步 法 的 另 一 特 点 是可 以 推 广 到 三维问 题的 研究中.本文将这种方法用于求解 A L E描述的N o v i e r - S t o k e s 三 维数值模拟时, 对压力和速度 采用相同的
8、八 结点 六面 体单元, 方程中,在进行 性插 值函 数, 编制了 计算机计 算程序, 插值函数采用同阶 一次线 进行了 算例 分析 并对结 果和有 关实验结 果进行了比 较. 基本方程 A L E描述下的随体导数可写为 D f f o f o f 司, = a r f、 十 a x ,-o il! c , = u 一 u , 称为A L E 描述 下 的 对 流 速 度 , 其中 u ,为 流 体 质 点 的 物 质 速 度 , u , 标 系 下 的 网 格 速 度 ,X为 物 质 坐 标 ,x 为 参 考 坐 标 下 的 参 考 坐 标. 推 导出A L E 不 可 压 枯 性 流 体
9、的 N av ie r -S to k e s 方 程 组 如 下 为参考坐 描述下的 连 续 方 程 : (9 u , 二 。 a x , 运 动 方 程 为 : D u,D t Ix 、 a u,a x; 二 aasa -XI 关 “ 方 “ “ : 口 , 二 一 p S 0 二 (会 + 会 ),9xj e x, 边 界 条 件 为: 七 月 , = 。 口 ,. n j二a , P=P 初始条件为:u ( x , 0 ) 二 u o ( x ) P ( x , 0 ) = P a ( x ) 其 中 S . , S r分 别 表 示湿 壁 面 和自 由 面, 数, s 为 体积力.
10、P , P和v 分别为流体的密 度、 压力和 运动 枯性系 . 一 , 2 A L E 分步有限元法数值离散方程 (l)(2) 首先对1 3 a v i e r - S t o k e s 方程 进行时间离散得 u r .尸二 0 u; 二 : 一 【一 “J 告 , ,“ 一 u ,0, + u j,,二 一 ”“ , 推导出A L E有限元分步格式如下 1 ) 求出 中 间 速 度衅+ 1 u “ +1 = u , 一 A t c ,“ u n , 一 v ( u ,“, +, 一 f i. w l 2 ) 压 力 R n + l 计 算公 式 ( 3 ) _ + I P - 1 P.r
11、r = 几 u rj 口 c ( 4 ) 需要指出的 是在求解此压力泊松 方程时,要考虑压力 梯度边界条件 p 犷 , 。 。 二 P f a“+ i n o n S w 3 ) 进 行 速 度 修 正 , ,4 出 速 度 u ,“ 1 计 算 公 式 n + 1 = n + Iu , = u r一 A t +1 p .r ( 5 ) 数据实验表明在引人中间速度时, 不是舍弃整个压力项, 而是舍弃压力 增量项,计 算效果 明显比 前者要好,由( 3 ) 一 ( 5 ) 式分别由G a l e r k i n 方法可推出改进后的数值离散方程 、 ! ” = M ho u ; 一 B ; 二 、
12、 十 去 C , ,p ; + D . , ,u j 一 F m “ 一 m A “+1 +1A Qa p e = 一 景 C 1 - + + A -# 十 Q a“+ 一 Q a n+1U n+1M u 0 二 、 1u +1 一 尝 lC v p ;+l 一 、 P ; ) 其中 各单元系 数矩阵的计算可参考文献 刀 . 在以上各 式中,a , 刀 表示单 元的 结点数目,: , .1 ( = 1 , 2 , 3 ) 表示空间维数. 3 三维液体自由液面上结点法矢量的数值计算 在 利 用A L E 描 述 方 法 跟 踪 三 维 液 体自 由 液面 时 , 网 格 更 新 速 度 和 流
13、体 质点 速 度 之间 在 自 由液面上应满足约束条件:舀 ” 二 护。 ,因此在每一时间步上,一 旦自由 液面上 每个 结 点的流体质点速度求得后,只 要求出 结点的法向矢量,就可以 求得网格 更新速度,对于 每 一自由液面单元的法向矢量的计算可参考文 献 3 ,当每一自由 液面 单元的法向矢鱼求出 后, 就可以计算结点的法向 矢量.本文用 加权系数平均法如图1 所示, 、= i 生 ,_i S , ( 6 ) 其中1 表示与结点k 相邻单元的数目, n 、表示 相 邻 的 第I 个 单 元 的 单 元 法向 矢 量, S ; 一 表 示该 单元 的面积.而单 元面积可根据下式由 G a u
14、 s s 数值积分求得. d S = I G Id d 。(7 ) 最后求 得结点k 的单位法向量为 结点 图 I 结点法向矢童的面积加权系数法 N, 刀 = 了 甲 了 - 片 Nk I ( g ) 4 三维液体非线性大幅晃动的数值模拟 算 例i s 圆 筒形 贮 腔 半 径为 r , = 0 .3 m , 充 液 深度( 沿二方向) 凡= 0 .3 m, 受 到 横 向 ( 沿o x 方 向 ) 激励.f * 二 ( A s i n ( 2 nf ) , 0 , 一 S ) , 其中 A = 0 .0 5 / s 2 , f = 1 2 H z , 而 一 阶 反对 称 模 态 。由 文
15、。中 的 公 式 给 出 : , 一 牛 。 U th (A . ) , ; 满 足 d J ,(D = 0 , J ,为 第 一 类 2 s Vr o r o ” ” r 一 d 贝塞尔函数 03 6 rr ; ? 而 右I*而 八 ;- 了 八 , h n八 叭 j一 刊四叫L J a 5,6 月叼,rjn乃 33仓, 00八以 日偏荀月8而七二58岛 T i me ( s ) 图2 1 - 腔左 、 右 壁面 处 的 波高 变 化 时间 历 程( v = 1 0 -6 1 气己旷 图 3 T i me ( s ) 液体 对 k -1 腔沿 0 Y 方向 的 羌 动 力 变 化 时间 沥程
16、 v = 1 0 0 R 户 . 。 , 电 . 甘 7 甸 叼。 已 呜 1 侣 。 , 启. 之 J 召 。 丁 , 0 心 月 臼 . , .口 之 当 ,公介四 家 户 心口二 , 、户 乏 巴二j 不同时刻三维自 由液面形状 ( v =1 0 -1 仲酗 一 已. 二 r , -e , , 、 0. Me , 价 5 带圆环形隔板圆筒形贮腔中三维大幅晃动的 数值模拟 算 例2: 在 算例1 的 基 础 上 加 上 环形 肋 板 , 其 中 各 参 数的 意义 如下 :。 为 肋 板的 厚 度 1 7 0 取为。 O l m,、为肋板的宽度,d 为 肋板距静液面的高度. 本文对不同的,和d 值进行了 数值模拟计算. 图5带圈环形隔板的固筒形贮腔 图6圈环形隔板上的网格划分示意图 0 .3 2 3 03 1