Microsoft Word - 9-N doc

Similar documents
第 48 卷第 9 期中南大学学报 ( 自然科学版 ) Vol.48 No 年 9 月 Journal of Central South University (Science and Technology) Sep DOI: /j.issn

24 26,,,,,,,,, Nsho [7] Nakadokoro [8],,,, 2 (Tradtonal estmaton of mage Jacoban matrx), f(t 1 ) p(t 2 ) : f(t 1 ) = [f 1 (t 1 ), f 2 (t 1 ),, f m (t

第 29 卷第 9 期 Vol. 29 NO. 9 重庆工商大学学报 ( 自然科学版 ) J Chongqing Technol Business Univ. Nat Sci Ed Sept X * ABAQUS 1 2

34 7 S R θ Z θ Z R A B C D PTP θ t 0 = θ 0 θ t 0 = 0 θ t 0 = 0 θ t = θ θ t = 0 θ t = 0 θ t V max θ t a max 3 θ t A θ t t 0 t / V max a max A = 3 4 S S

33 5 Vol.33,No JournalofHebeiUniversityofScienceandTechnology Oct.2012 : (2012) /,, ( 河北科技大学机械工程学院, 河北石家庄 ) : 利用计算流体

本 土 天 蝗 傳 奇 - 台 灣 大 蝗 生 活 史 及 生 態 習 性 的 研 究 摘 要 台 灣 大 蝗 在 交 配 時 警 覺 性 降 低, 蝗 會 背 著 蝗 跳 到 遠 處, 但 不 會 飛, 肚 子 餓 時 會 進 食, 但 蝗 不 會 交 配 後 蝗 會 選 擇 土 質 堅 實 植

5 551 [3-].. [5]. [6]. [7].. API API. 1 [8-9]. [1]. W = W 1) y). x [11-12] D 2 2πR = 2z E + 2R arcsin D δ R z E = πr 1 + πr ) 2 arcsin



/MPa / kg m - 3 /MPa /MPa 2. 1E ~ 56 ANSYS 6 Hz (a) 一阶垂向弯曲 (b) 一阶侧向弯曲 (c) 一阶扭转 (d) 二阶侧向弯曲 (e) 二阶垂向弯曲 (f) 弯扭组合 2 6 Hz

BISQ理论模型与声波测井响应研究

一 大聖 起源與簡述: 1 本宮溯自清光貳年(西元一八二二年) 有位唐先賢黃迎 奉請 北極玄上帝金尊東渡來台 由烏石港上陸至四圍保埤口 後遷到 外澳石空中路嶺上 茅蓋神堂供人膜拜 化蒼民墾土農耕 後澤 被農信仰日深 為祈求石空 風調雨順 保佑黎民平安 先賢數 人等誠心合建永固廟堂 因高水長 所稱為(

Vol. 36 ( 2016 ) No. 6 J. of Math. (PRC) HS, (, ) :. HS,. HS. : ; HS ; ; Nesterov MR(2010) : 90C05; 65K05 : O221.1 : A : (2016)

于达等 :-, : ; τ = ( ) / ρ ; ; ; ρ ; θ ; ; τ ; ;,,,,,, : = ( +.. ) ρ = μ = π = ρ πμ ; ; ; μ ( ), : = -.. ( ), ( ),, :. = [ -. ] ( ) ( ) 试验条件,,. ~. /,....

2 Abstract 厦门大学博硕士论文摘要库 1

T e = K 1 Φ m I 2 cosθ K 1 Φ m I cosθ 2 1 T 12 e Φ / 13 m I 4 2 Φ m Φ m 14 I 2 Φ m I 2 15 dq0 T e = K 2 ΦI a 2 16

第 卷第 期 年 月 杭州师范大学学报 自然科学版!"# "$%% &'" $((& ) 31-,- 方程的多辛 " &&!! 格式 王俊杰 王连堂 杨宽德 普洱学院数学系 云南普洱 西北大学数学系 陕西西安 摘 要?4 3!3 方程作为一类重要的非线性方程有着许多广泛的应用前景 基于 :#"$7&

[1] Nielsen [2]. Richardson [3] Baldock [4] 0.22 mm 0.32 mm Richardson Zaki. [5-6] mm [7] 1 mm. [8] [9] 5 mm 50 mm [10] [11] [12] -- 40% 50%


Fig1 Theforceappliedtothetrainwhenrunning :w = w j +w q (3) :w = w = w 0 +w j (4) w i 121 基本阻力 w r = 600 R ( N/kN) (8) :R : [2] w s [3] w s =0

07-3.indd

cm /s c d 1 /40 1 /4 1 / / / /m /Hz /kn / kn m ~

福建中联房地产开发集团有限公司

自然科学版 预处理 视盘粗定位 视盘垂直坐标的粗定位 视盘水平坐标的粗定位

Ashdgsahgdh


~ 10 2 P Y i t = my i t W Y i t 1000 PY i t Y t i W Y i t t i m Y i t t i 15 ~ 49 1 Y Y Y 15 ~ j j t j t = j P i t i = 15 P n i t n Y

Specification for structura design of buried prestressed concrete pipeine of water suppy and sewerage engineering CECS 40:

[9] R Ã : (1) x 0 R A(x 0 ) = 1; (2) α [0 1] Ã α = {x A(x) α} = [A α A α ]. A(x) Ã. R R. Ã 1 m x m α x m α > 0; α A(x) = 1 x m m x m +

#4 ~ #5 12 m m m 1. 5 m # m mm m Z4 Z5

黑嘴子公示版本.doc

CHIPS Oaxaca - Blinder % Sicular et al CASS Becker & Chiswick ~ 2000 Becker & Chiswick 196

(H ~z 2,.3 V T A(HlH -H+BH - (A(z 2lz 2 -z 2 +Bz 2 (5,H ADCP,z 2 R ε,(, (0~z V B A(lz -+B (6 R (V 2 i - 珚 2 V (H i -H - 2 i i u m u * z 0,A

1 引 言 大 陆 与 台 湾 两 地 之 间 的 交 流 与 日 剧 增, 大 量 与 台 湾 有 关 的 信 息 进 入 了 大 陆 居 民 的 生 活 随 着 交 流 的 不 断 深 入, 我 们 发 现 台 湾 国 语 和 我 们 所 使 用 的 普 通 话 存 在 一 定 的 差 别 台

南部县定水镇汽车加气站项目环评表送审--公式本.doc

* CUSUM EWMA PCA TS79 A DOI /j. issn X Incipient Fault Detection in Papermaking Wa

Revit Revit Revit BIM BIM 7-9 3D 1 BIM BIM 6 Revit 0 4D 1 2 Revit Revit 2. 1 Revit Revit Revit Revit 2 2 Autodesk Revit Aut

99710b43ZW.PDF

, (, ) ABSTRACT Ths paper ntroduces the western theores of dscrmnaton n labor market, ncludng man representatve's theores, ant-dscrmnaton measur

京都大学防災研究所年報 第56号(平成24年度)

國家圖書館典藏電子全文

PCA+LDA 14 1 PEN mL mL mL 16 DJX-AB DJ X AB DJ2 -YS % PEN

Y m G C I IMC C II IMC R Y = + C I IMC F f G - G^ + - F f G^ C I IMC D + C I IMC F f G - G^ 9 D() R() C IMC() Ⅱ CIMC() U() Ⅰ G() Y() Y m() - G

Mnq 1 1 m ANSYS BEAM44 E0 E18 E0' Y Z E18' X Y Z ANSYS C64K C70C70H C /t /t /t /mm /mm /mm C64K

<4D F736F F D20B0AAB6AFA5ABBEF4C059B0EAA4A BEC7A67EABD7B2C4A440BEC7B4C1AAC0B77CBBE2B0ECB2C4A454A6B8B77CC4B3B04FBFFD2E646F63>

(,00);,, (,,00);,,,, (,00) (,, 00;,00),, (00) IPO, IPO,,,, ( ),,,, (Loughran,Rtter,00;Rtter,003), IPO,IPO, (Rtter,003;Jenknson et al.,006),, IPO,, 5%(

Microsoft Word - 系统建设1.doc

f 2 f 2 f q 1 q 1 q 1 q 2 q 1 q n 2 f 2 f 2 f H = q 2 q 1 q 2 q 2 q 2 q n f 2 f 2 f q n q 1 q n q 2 q n q n H R n n n Hessian

Microsoft Word - 11-秦华伟.doc

2009 年第 54 卷第 12 期 : 1779 ~ csb.scichina.com SCIENCE IN CHINA PRESS,,, ;, , 200

1.0 % 0.25 % 85μm % U416 Sulfate expansion deformation law and mechanism of cement stabilized macadam base of saline areas in Xinjiang Song

Microsoft Word - A _ doc

1.2 编 制 依 据 法 律 法 规 (1) 中 华 人 民 共 和 国 环 境 保 护 法 ( ); (2) 中 华 人 民 共 和 国 大 气 污 染 防 治 法 ( ); (3) 中 华 人 民 共 和 国 水 污 染 防 治 法 ( )

85% NCEP CFS 10 CFS CFS BP BP BP ~ 15 d CFS BP r - 1 r CFS 2. 1 CFS 10% 50% 3 d CFS Cli

Microsoft Word - 0 目录及封面1.doc

% GIS / / Fig. 1 Characteristics of flood disaster variation in suburbs of Shang

2 : 237.,. [6 7] (Markov chan Monte Carlo, MCMC). MCMC, [8 9].,,, [0 ].,, : ),,,.,, ; 2),,.,.,. : ),.,,. ; 2),.,,. ; 3), EM, EM,.,, EM, EM. K M,.,. A

1 引言

2 期 张 慧 : 海 南 东 北 部 地 区 地 壳 地 震 各 向 异 性 特 征 初 步 研 究 构 造 背 景 海 南 岛 位 于 中 国 大 陆 东 南 沿 海 地 震 带 的 雷 琼 地 震 亚 区, 包 括 海 南 岛 和 雷 州 半 岛 及 附 近 海 域, 大 致 以

Microsoft Word 張嘉玲-_76-83_

IPCC CO (IPCC2006) 1 : = ( 1) 1 (kj/kg) (kgc/gj) (tc/t)

mm 5 1 Tab 1 Chemical composition of PSB830 finishing rolled rebars % C Si Mn P S V 0 38 ~ 1 50 ~ 0 80 ~ ~

Microsoft PowerPoint - aspdac_presentation_yizhu

kg/m3 1.55g/cm3 m3 580kg 37.42% %

untitled


资源 环境 生态 土壤 气象

SVM OA 1 SVM MLP Tab 1 1 Drug feature data quantization table

试卷

Applied Mathematics and Mechanics Vol. 34 No. 9 Sep ISSN GPU Boltzmann *? GPU Boltz

untitled

目 录 中 文 摘 要 1 英 文 摘 要 第 一 章 综 述 3 第 二 章 风 险 分 解 方 法 介 绍 6.1 风 险 分 解 预 备 知 识 6. 本 文 采 用 的 风 险 分 解 方 法 6 第 三 章 风 险 序 列 估 计 方 法 及 所 用 数 据 描 述 10 第 四 章 数

12-1b T Q235B ML15 Ca OH Table 1 Chemical composition of specimens % C Si Mn S P Cr Ni Fe

FZ1.s92

1 VLBI VLBI 2 32 MHz 2 Gbps X J VLBI [3] CDAS IVS [4,5] CDAS MHz, 16 MHz, 8 MHz, 4 MHz, 2 MHz [6] CDAS VLBI CDAS 2 CDAS CDAS 5 2

建设项目基本情况

Microsoft Word - 报告报批全文公示版0119


é ê


Microsoft PowerPoint - 概率统计Ch02.ppt [Compatibility Mode]

作为市场化的人口流动


N E LaCoste ±7 µgal ±10 µgal (1 Gal = 0.01 m/s 2 ) Kringe [8] ξ η = 1 2π

Thesis for the Master degree in Engineering Research on Negative Pressure Wave Simulation and Signal Processing of Fluid-Conveying Pipeline Leak Candi

2011第1期第二部分

successful and it testified the validity of the designing and construction of the excavation engineering in soft soil. Key words subway tunnel

目 录 1 前 言 任 务 由 来 环 境 影 响 评 价 工 作 过 程 关 注 的 主 要 环 境 问 题 环 评 主 要 结 论 7 2 总 则 评 价 目 的 评 价 原 则 编 制 依 据 环

11 25 stable state. These conclusions were basically consistent with the analysis results of the multi - stage landslide in loess area with the Monte


2.PFA 사출 Valve 지윤주

untitled

84 Ab M Ⅰ B K Ⅱ Ba M Ⅲ Ⅰ


2015-6

2 决 (1) (2) (3) 408AD-537AD AD-509AD 2.66 ( 1990) 衆 ( ) 6.5 1) 1)

标题

《自动化学报》作者加工稿件须知

untitled

中华人民共和国国家标准 气体灭火系统施工及验收规范GB50263-97 Code for installation and acceptance of gas fire-extinguishing systems 主编部门:中华人民共和国公安部 批准部

Transcription:

017 年第 6 卷第 13 期 :139 ~ 1401 中国科学 杂志社 SCIENCE CHINA PRESS 基于多尺度混合有限元的离散裂缝两相渗流数值模拟 张庆福 黄朝琴 * 姚军 * 王月英 李阳 中国石油大学 ( 华东 ) 石油工程学院 青岛 66580 * 联系人 E-mal: em.group.up@gmal.om; RCOGFR_UPC@16.om 016-05-09 收稿 016-07-5 修回 016-07-5 接受 016-11-04 网络版发表国家自然科学基金 (514049 5134007 51490654) 山东省自然科学基金(ZR014EEQ010) 中央高校基本科研业务费专项(15CX05037A 14CX0507A) 和青岛市博士后项目 (01548) 资助 摘要裂缝性介质通常具有多尺度特征 离散裂缝模型虽具有计算精度高 拟真性好的优点 但传统数值方法在解决此类多尺度流动问题时 难以突破计算量大的瓶颈 不利于实际应用. 对此 本文将离散裂缝模型和多尺度混合有限元相结合 仅需进行宏观大尺度计算 通过多尺度基函数来刻画小尺度裂缝精细流动特征 在保证计算精度的同时降低了计算量. 在小尺度上 采用模拟有限差分法构建离散裂缝模型的多尺度基函数 该方法不仅具有良好的局部守恒性 而且适用于任何复杂离散裂缝网格. 文章详细阐述了离散裂缝模型多尺度混合有限元两相流动数值计算格式的建立 重点介绍了如何使用模拟有限差分法构建离散裂缝模型的多尺度基函数 并采用超样本技术进一步提高计算准确性. 数值结果表明 本文计算方法不仅能够准确捕捉离散裂缝性介质中的精细流动特征 而且具有很高的计算效率. 关键词 裂缝性介质 多尺度混合有限元法 离散裂缝模型 模拟有限差分 两相渗流 现代很多工程和科学问题都有多尺度解 如复合材料 非均质多孔介质问题等 是典型的多尺度问题. 裂缝性油藏中 裂缝空间尺度跨越大 [1~5] 流体在裂缝性油藏中的流动呈现显著多尺度特征. 目前 主要采用双重介质 等效连续介质和离散裂缝模型作为裂缝性介质的流动模型. [67] 双重介质模型认为裂缝性介质中存在两个平行的渗流系统 : 裂缝系统和基岩系统. 由于双重介质模型中窜流函数难以确定 很难应用于两相和多相流问题 [89]. 等效连续介质模型将介质视为一个假想的连续体 其操作简单且计算效率高 [10]. 然而 对于两相或多相流问题 尚未有成熟的理论和方法获取相应的等效参数 [11]. 离散裂缝模型则对介质中的裂缝予以显示表征 由于其计算精度高 拟真性好 该模型越来越受到人们的重视 本文选取离散裂缝模 型作为流动模型. 裂缝性油藏中 裂缝具有强烈多尺度性 为了得到较精细的解 需要建立高精度的地质模型 可达数百万甚至数亿个网格单元. 此时 如果基于传统的数值计算方法求解离散裂缝模型 巨大的计算量将超出当今计算机的计算能力. 尺度升级 [113] 方法虽然可以降低计算量 但浪费了大量的精细模型的数据信息 无法反映小尺度流动特征. 因为小尺度信息对整体流动模拟有重大影响 所以尺度升级法模拟精度不高. 因此寻找既能降低计算量又能保证计算精度的新型数值方法势在必行. 近年来 国内外学者提出一些多尺度数值模拟方法 [14~17] 如多尺度有限元法 [14] 变分多尺度方法 [16] [17] 广义多尺度有限元法等. 多尺度方法是在粗网格上求解全局流动方程 并且利用多尺度基函数来捕捉小尺度特征. 通过使用多尺度基函数离散全 引用格式 : 张庆福 黄朝琴 姚军 等. 基于多尺度混合有限元的离散裂缝两相渗流数值模拟. 科学通报 017 6: 139 1401 Zhang Q F Huang Z Q Yao J et al. Two-phase numeral smulaton o dsrete rature model based on multsale mxed nte element method (n Chnese). Chn S Bull 017 6: 139 1401 do: 10.1360/N97016-00584 016 中国科学 杂志社 www.shna.om sb.shna.om

局流动方程 可以得到一个规模更小的方程组从而减少计算量 同时可以反映介质内的非均质性. 多尺度有限元法虽然可以得到反映小尺度特征的解 但 [18] 是这些解却并非局部守恒的. Chen 等人进一步提出了多尺度混合有限元法 该方法基于混合有限元来构造多尺度基函数 除了具有多尺度有限元的优点外 同时在粗网格上也满足局部守恒 适于处理多相复杂流动过程. 多尺度方法近年来也开始应用于裂缝性油藏数 [19] 值模拟中. Gulbransen 等人使用多尺度混合有限元法求解缝洞型油藏单相流问题. 其基函数采用 Stokes-Brnkman 方程计算 导致其难以应用于两相 [0] 流流动模拟. Zhang 等人将多尺度有限元法应用于离散裂缝模型的研究中 并利用全局信息提高计算精度. 最近 广义多尺度有限元法被应用于裂缝性介质页岩气流动模拟中 [1]. 本文将多尺度混合有限元推广到裂缝性油藏两相流流动模拟中 在减少计算量的同时保证计算精度. 文章阐述了多尺度混合有限元的基本原理 建立了离散裂缝模型多尺度基函数的模拟有限差分计算格式 并采用 IMPES 方法对其两相流问题进行了求解 最后通过对比多尺度解和参考解验证了方法的正确性和程序的鲁棒性. 1 两相渗流数学模型 仅考虑不可压缩流体的等温渗流过程 给出了湿相 ( 下标 w) 以及非湿相 ( 下标 n) 的控制方程. 如下 : S q nw (1) t kr K p k gz nw () S n Sw 1 (3) p Sw pn pw (4) 式中 为孔隙度 p 为流体压力 为相 的速度 ρ 是流体密度 K 为渗透率张量 kr 为相对渗透率 q 为源汇项 S 为饱和度 是流体黏度 g 是重力加速度 p 为毛管力. 式 (1)~(4) 可写为 KpK G (5) w w n n q (6) t 其中 = n w 为总速度 qt qn qw 为总源汇项 k r / 是相 的流度 t n w为总流度. 图 1 离散裂缝模型 Fgure 1 Grd shemats o dsrete rature model 本文在粗网格上使用 Dary 方程对速度和压力进行近似 在细网格上基于模拟有限差分法求解离散裂缝模型获得多尺度基函数. 在离散裂缝模型中 基于流量等效原理 沿裂缝开度方向流体相关参数不变 因此可对裂缝进行降维处理. 对于二维问题 裂缝简化为一维线单元 如图 1 所示. 多尺度混合有限元法.1 混合变分形式及离散 在粗网格上 使用混合有限元对 Dary 方程进行 d 离散. 令 R d 3 是有界闭区域 L 为 上的平方可积空间. 定义如下函数空间 : L0 p pl pd 0 d d H L L 1d 令 U 和 V 是 L 和 H 0 (6) 相应的离散形式为求 p UV 积分方程组 : 1 Sw d 的有限维子空间 式 (5) 和 u K pud w w w n n n 使下述等效 u S S G d (7) 对所有的 ( u l) U V都成立. 设 和速度的单元近似表达式如下 : l d lqd (8) t 和 分别是 U 和 V 的一组基函数. 压力 ψ p p (9) k k 为保证 n 在相邻单元交界面处连续 引入单元表面压力. 把上式代入式 (7) 和 (8) 可以得到线性方 1393

017 年 5 月第 6 卷第 13 期 程组 : 式中 B C D Q g T C 0 0 p q T D 0 0 λ 0 Q 是外向流量向量 为表面压力向量. B Sw p 为单元压力向量 1 (10) λ d K ψ j C d j ψ D ψ n j d. 当局部流动问题包含裂缝时 采用模拟有限差分法对离散裂缝模型的控制方程进行离散 []. 基于模拟有限差分法 可以准确计算基岩 - 裂缝与裂缝 - 裂缝之间的流速. 由 Dary 定律可知 对于任一网格单元 其边界上的法向速度 可以写为 T e p (11) 其中 1 m T 为传导矩阵 T 界面数 e 1 1 T m 为网格单元 p 为单元压力 为边界面压力. 将基岩系统 ( 下标为 m) 和裂缝系统 ( 下标为 ) 耦合在一起可以得到离散裂缝模型的模拟有限差分数值计算格式 Bm Cm Dm 0 0 m gm T m 0 0 0 0 C pm qm T T Dm 0 0 C 0 π m q (1) 0 0 C B D g T 0 0 0 D 0 π 0 上述方程的系数矩阵具体如下 : 1 T 1 e 1 B C 1 T N e e N e I 1 D I N e 其中 N e 为网格单元总数 I E m. 相比于混合有限元法 模拟有限差分法不仅拥有良好的局部守恒性 而且对复杂网格有很高的适应性 适用于任何复杂 网格系统.. 多尺度混合有限元基函数 令 = 为区域 的粗网格部分 h k K E k 的细网格部分 只要满足 E 0 就有 Ek. 定义粗网格交界面 j 是 每一个 对应一个速度基函数 每一个粗网格单元 对应一个压力基函数. 如果区域 j 内仅包含基岩单元 速度基函数 满足下面局部流 动问题 : t K (13) ω x x ωj x xj (14) 0 x 满足 n 0 x. ( x) 相当于 上的加权函数 将 d( ) 分布到细网格上 满足 ω ( x)ds =1. ( x) 通常取为 1/ qdx 0 ωx (15) qx/ q d qdx 0 但这种取值会导致低渗透区域计算的流速偏高 本 文选取 ω x x d qdx 0 qx q d qdx 0 (16) 其中 x trae( K ( x)) d trae( A ) 表示矩阵 A 的特 征值的和. 这种选取方法可以有效避免均一化产生 的误差 提高基函数计算精度. 然后通过上面介绍的 混合有限元法对式 (13) 和 (14) 进行求解. 如果粗网格包含裂缝单元时 使用离散裂缝模 型表征裂缝 多尺度基函数满足局部流动问题 K (17) a t K (18) n ω x x a ωj x xj 0 x (19) 使用前面介绍的模拟有限差分法对上述方程组 进行数值求解 基函数满足 n 0 x..3 超样本技术 为了避免封闭边界影响基函数计算的精确性 采用超样本技术计算多尺度基函数. 超样本技术就 是在计算多尺度基函数时 使目标区域包含更多的 小尺度网格 即在更大的区域上计算多尺度基函数 如图 所示. 超样本技术可以明显提高强非均质性油 藏数值模拟的精确性. 1394

图 ( 网络版彩色 ) 超样本单元示意图 Fgure (Color onlne) Shemat desrpton o oersampled regon 多尺度方法计算时 裂缝性介质内的流体流动特征主要通过多尺度基函数体现出来. 当裂缝穿过粗网格边界或者裂缝端点靠近粗网格顶点时 由于封闭边界的影响 基函数对裂缝附近流体流动的表征效果较差. 此时在超样本单元 K E 上计算基函数 计算得到的基函数可以充分考虑裂缝导流作用 提高多尺度数值模拟的精确性. 但是 使用超样本技术时 随着超样本单元体积的增大 精度的提升会变缓 同时计算量会越来越大. 因此在处理裂缝性介质时 应综合考虑计算精度与计算量 结合裂缝分布 选择合适的超样本单元..4 全局大尺度方程 由于多尺度基函数是在相邻粗网格上计算 因此为了建立大尺度方程 首先将局部流动问题求解得到的多尺度基函数分为两部分 其中 H H (0) j ( E) E / H j ( E) 0 E ( E) E H j j ( E) 0 E j H (1) 令 Ψ 是以所有基函数 作为列向量的矩阵 I 是粗 网格到细网格的变换矩阵 ( 若第 j 个粗网格包含第 个细网格 I 1 否则 I 0 ). 根据多尺度方法原 理 小尺度速度和压力场可由相应的基函数近似表 示出来 那么 p 可近似表示为 粗网格表面压力 ψq p Ip () λ 可表示为 λ ψ n d s (3) 令 J 是粗网格表面到细网格表面的变换矩阵 ( 若第 j 个粗网格表面包含第 个细网格表面 J 1 否则 J 0 ). 根据式 () 对小尺度方程进行组装得到大尺度方程组 T ψ 0 0 B C D 0 T T T 0 I 0 C 0 0 p I q (4) T T 0 0 J D 0 0 λ 0 综合式 () 和 (3) 大尺度方程组可写为 B C D Q 0 T C 0 0 p q (5) T D 0 0 λ 0 T T T 其中 B Ψ B Ψ C Ψ C I q I q D T ψ C J. 小尺度速度可由相应的基函数表示出来 即 ψq. 3 饱和度方程求解 3.1 有限体积法计算格式 由式 (1) 和 () 可得湿相饱和度方程 Sw w qw t (6) K p K G w w n n w n (7) 本文采用有限体积法对饱和度方程进行离散. 首先在 上对方程进行积分 可得 S d w n p n w n t K K G n d = q d. (8) w 方便起见 这里省略掉了 S w 下标 w. 应用 准则 可得到下述离散格式 : k1 k 1 k1 k S S F ( S ) (1 ) F ( S ) t γ ( k q ) w S (9) 其中 γ w n F ( S) S n K p n K ( ) Gn d n w o (30) 1395

017 年 5 月第 6 卷第 13 期 式中 上标 k 表示第 k 个时间步. 再边界面 A k 上 采用上游迎风格式计算 如下 : w ( S) n 0 w ( S) (31) w ( Sj) n 0 对饱和度方程显示求解 即 0. 为达到计算稳定 性 时间步长采用 CFL 条件 如下 : 式中 t n max{ ( S)} w 0 S 1 max( q 0) mn( 0) n γ S 1 * w w w * * S S S 1 Sw Sor S (3) (33) 式中 S * 为归一化后的水相饱和度 ; S w 为束缚水饱和 度 ; S or 为残余油饱和度. 3. 裂缝交叉处饱和度计算 当裂缝相交时 关键问题是裂缝交叉处饱和度 的计算. 主要有两种处理方法 : 一是基于 Delta-Star 传导率计算的上游迎风格式 [3] 另一种则是上游迎 风加权计算格式 [4] 本文采用精度较高的加权格式. 假设有 N I 个裂缝单元 e 相交于 I. 每个裂缝单元相应的分流函数为 w e 在交叉 处的渗流速度为 e ; 定义 I 处的流入和流出如下 : e 0 0 N 流出 e 0 N N 流入. I 由质量守恒定律可知 N N e e N1 1 (34) I. (35) 进一步 由上游迎风计算格式定义 可得 N N N I (36) w e w I e w I e N1 1 1 因此 裂缝交叉处 I 的上游迎风加权分流函数为 4 数值算例 w I NI w e N1 N 1 e. (37) 本节首先通过一个数值算例与实验结果的对比 验证了程序的正确性 ; 然后通过复杂离散裂缝模型数 值算例进一步验证了方法的正确性和程序的鲁棒性. 4.1 简单离散裂缝模型算例 考虑图 3 和 4 所示的物理模型 模型大小为 1 m 1 m 0.05 m 内部流动可视为平面流动. 图 3 中包含 一条垂直裂缝 图 4 中包含两条交叉裂缝. 物理模型 均采用玻璃沙 (160~180 目 ) 结合环氧树脂胶结压实而 成 然后由透明有机玻璃板封装. 基岩可视为是均质 的 渗透率 K m =10 m 孔隙度 0.4. 裂缝开度均 为 1 mm 渗透率 K =a /1=8.33 10 4 m. 水的黏度 w =1 mpa s 油的黏度 0 =5 mpa s 水的密度 w =1000 kg/m 3 油的密度 0 =800 kg/m 3. 模型初始时刻饱和油 油水相的流度为 w S w w 0 1 Sw 0. 针对上述两个物理模型采用本文方法进行数值模拟 结果如图 3 和 4 所示. 通过对 比数值模拟结果和实验结果 可以发现数值计算结 果与实验结果基本一致 验证了本方法的正确性. 通 过对比多尺度数值结果与离散裂缝模型的数值结果可以看出 多尺度结果与完全的离散裂缝模型的结果也基本一致. 图 5 给出了多尺度方法求解的速度与传统数值方法求解离散裂缝模型得到的速度的对比曲线. 图 5(a) 表示单裂缝模型中沿 x=0.5 m 的速度分布 图 5(b) 表示交叉裂缝模型中沿 y=0.5 m 的速度分布. 图 5 说明速度误差较小 多尺度方法求解的结果较为准确. 同时 多尺度方法计算时间相比于传统数值方法缩短了 40%. 对比结果表明 多尺度方法在保证计算精度的同时可以大幅缩减计算时间. 4. 复杂离散裂缝模型算例 考虑一个复杂裂缝介质模型 研究区域为 180 m 180 m 模型包含 10 条裂缝 如图 6 所示. 基岩为均质各向同性 孔隙度 =0. 渗透率 K m =1 m. 所有裂缝开度均为 a=1 mm 渗透率 K =a /1=8.33 10 7 (10 3 m ). 油水的物性参数与 4.1 算例一致. 油藏初始含水饱和度为零 基岩和裂缝的油相流度 o 1 Sw o 水相流度 w Sw w. 采用 10 10 粗网格划分 细网格采用三角网格. 采用相对 L 范数表示速度相对误差 其中 e ms 表示参考解得到的速度 (38) ms 表示多尺度方 1396

图 3 ( 网络版彩色 ) 单裂缝模型及其含水饱和度结果对比. (a) 物理实验 [] ; (b) 多尺度解 ; () 参考解 Fgure 3 (Color onlne) Comparson o water saturaton proles between expermental results and multsale results. (a) Physal expermental results [] ; (b) MsMFEM results; () DFN results 法得到的速度. 图 7 给出了 0.5 PV 时的含水饱和度分布对比 其数值计算结果表明 多尺度混合有限元能够准确描述离散裂缝油藏整体流动特征. 图 8 给出了参考解和多尺度解得含水率曲线对比 两者基本一致 进一步验证了多尺度方法的正确性与程序的鲁棒性. 针对图 6 所示的裂缝介质模型 图 9 给出了多尺度混合有限元方法在不同粗网格剖分下得到的速度相对误差. 从误差分布图可以看出 粗网格划分对多尺度模拟结果有一定影响 当 x 方向和 y 方向上的粗网格数量相近时 速度误差较小 模拟效果较好. 5 结论 (1) 离散裂缝模型对介质中的每条裂缝予以显 示表征 具有计算精度高 拟真性好的优点. 但随着现代地质建模技术的发展 地质模型可能包含数百万甚至数亿个网格 采用传统的数值方法对其进行求解 计算量巨大 将超出当今计算机的计算能力. 本文提出离散裂缝模型的多尺度混合有限元计算格式 在保证计算精度的同时大幅减小计算量 非常适合对裂缝性油藏进行精细流动模拟. 本文以 D 裂缝介质为例 为了进行 3D 离散裂缝模型的多尺度模拟 可以使用自适应的网格剖分技术. 这种网格剖分技术使用非结构化的粗网格 具有较高的灵活性 可以适应裂缝复杂的几何形态 进而进行 3D 的多尺度离散裂缝模拟. () 结合离散裂缝模型和模拟有限差分法构建速度基函数 可以准确捕捉基岩和裂缝之间的相互 1397

017 年 5 月第 6 卷第 13 期 图 4 ( 网络版彩色 ) 交叉裂缝模型及其含水饱和度结果对比. (a) 物理实验 [] ; (b) 多尺度解 ; () 参考解 Fgure 4 (Color onlne) Comparson o water saturaton proles between expermental results and multsale results. (a) Physal expermental results [] ; (b) MsMFEM results; () DFN results 图 5 ( 网络版彩色 ) 多尺度解与参考解速度对比. (a) 单裂缝模型中沿 x=0.5 m 的速度分布 ; (b) 交叉裂缝模型中沿 y=0.5 m 的速度分布 Fgure 5 (Color onlne) Veloty omparson or ne-sale and MsMFEM solutons. (a) Veloty dstrbutons along y-dreton on x=0.5 m o sngle rature model; (b) eloty dstrbutons along x-dreton on y=0.5 m o rossed rature model 1398

图 6 ( 网络版彩色 ) 离散裂缝模型及多尺度网格剖分 Fgure 6 (Color onlne) Dsrete rature expermental model and the multsale grd used n numeral smulaton 图 7 ( 网络版彩色 ) 注入体积为 0.5 PV 时饱和度分布图. (a) 参考解 ; (b) 多尺度解 Fgure 7 (Color onlne) Comparson o njeted lud saturaton proles between ne-sale numeral results and MsMFEM results at 0.5 PV. (a) Reerene soluton; (b) MsMFEM soluton 图 8 ( 网络版彩色 ) 含水率曲线对比 Fgure 8 (Color onlne) Water-ut omparson or MsMFEM 图 9 ( 网络版彩色 ) MsMFEM 在不同粗网格剖分下得到的速度误差 Fgure 9 (Color onlne) MsMFEM eloty error or derent oarse meshes 1399

017 年 5 月第 6 卷第 13 期 作用 从而可以对离散裂缝油藏进行精细流动模拟. 而且模拟有限差分法有很高的灵活性 原则上适用于任意复杂网格 具有明显优势. 使用超样本技术计算长裂缝和复杂裂缝基函数 可以提高基函数计算 的精度. (3) 多尺度混合有限元的基函数可以采用并行计算得到 进一步减少计算量. 因此 文章提出的方法对于裂缝性油藏数值模拟有很高的潜在价值. 参考文献 1 Huang Z Q Yao J L Y J et al. Numeral alulaton o equalent permeablty tensor or ratured uggy porous meda based on homogenzaton theory. Commun Comput Phys 011 9: 180 04 Huang Z Q Yao J Wang Y Y et al. Numeral smulaton on water loodng deelopment o ratured reserors n a dsrete-rature model. Chn J Comput Phys 011 8: 41 49 3 Yao J Huang Z Q Wang Z S et al. Mathematal model o lud low n ratured uggy reserors based on dsrete rature-ug network (n Chnese). Ata Petrol Sn 010 31: 815 819 [ 姚军 黄朝琴 王子胜 等. 缝洞型油藏的离散缝洞网络流动数学模型. 石油学报 010 31: 815 819] 4 Huang Z Q Yao J Wang Y Y et al. Numeral study on two-phase low through ratured porous meda. S Chna Teh S 011 54: 41 40 5 Huang Z Q Yao J L Y J et al. Permeablty analyss o ratured uggy porous meda based on homogenzaton theory. S Chna Teh S 010 53: 839 847 6 Barenblatt G I Zhelto I P Kohna I N. Bas onepts n the theory o seepage o homogeneous lquds n ssured roks [strata]. J Appl Math Meh 1960 4: 186 1303 7 Warren J E Root P J. The behaor o naturally ratured reserors. So Pet Eng J 1963 3: 45 55 8 Karm-Fard M Gong B Durlosky L J. Generaton o oarse-sale ontnuum low models rom detaled rature haraterzatons. Water Resour Res 006 4: 1 9 9 Bourbaux B. Fratured reseror smulaton: A hallengng and rewardng ssue. Ol Gas S Teh 010 65: 7 38 10 Zhou Z F. Theory on Dynams o Fluds n Fratured Medum (n Chnese). Beng: Hgher Eduaton Press 007 [ 周志芳. 裂隙介质水动力学原理. 北京 : 高等教育出版 007] 11 Wang Y Y Yao J Huang Z Q. Reew on low models o ratured meda (n Chnese). J Daqng Pet Inst 011 35: 4 48 [ 王月英 姚军 黄朝琴. 裂隙岩体流动模型综述. 大庆石油学院学报 011 35: 4 48] 1 Chrste M A. Upsalng or reseror smulaton. J Pet Teh 1996 48: 1004 1010 13 Durlosky L J. Numeral alulaton o equalent grd blok permeablty tensors or heterogeneous porous meda. Water Resour Res 1991 7: 699 708 14 Hou T Y Wu X H. A multsale nte element method or ellpt problems n omposte materals and porous meda. J Comput Phys 1997 134: 169 189 15 Wenan E Engqust B. The heterogeneous mult-sale method. Leture Notes Comput S Eng 00 347: 89 110 16 Juanes R. A aratonal multsale nte element method or multphase low n porous meda. Fnte Elem Anal Des 005 41: 763 777 17 Eende Y Gals J Hou T Y. Generalzed multsale nte element methods (GMsFEM). J Comput Phys 013 51: 116 135 18 Chen Z Hou T Y. A mxed multsale nte element method or ellpt problems wth osllatng oeents. Math Comput 003 7: 541 576 19 Gulbransen A Hauge V L Le K A. A multsale mxed nte element method or uggy and naturally ratured reserors. SPE J 009 15: 395 403 0 Zhang N Yao J Huang Z Q et al. Aurate multsale nte element method or numeral smulaton o two-phase low n ratured meda usng dsrete-rature model. J Comput Phys 013 4: 40 438 1 Akkutlu I Y Eende Y Vaslyea M. Multsale model reduton or shale gas transport n ratured meda. Comput Geos 015 0: 1 1 Huang Z Q Gao B Wang Y Y. Two-phase low smulaton o dsrete rature model usng a noel mmet nte derene method (n Chnese). J Chna Un Pet 014 38: 97 105 [ 黄朝琴 高博 王月英. 基于模拟有限差分法的离散裂缝模型两相流动模拟. 中国石油大学学报 : 自然科学版 014 38: 97 105] 3 Karm-Fard M Durlosky L J Azz K. An eent dsrete-rature model applable or general-purpose reseror smulators. SPE J 004 9: 7 36 4 Monteagudo J E P Froozabad A. Control-olume method or numeral smulaton o two-phase mmsble low n two-and three-dmensonal dsrete-ratured meda. Water Resour Res 004 40: 196 1 1400

Summary or 基于多尺度混合有限元的离散裂缝两相渗流数值模拟 Two-phase numeral smulaton o dsrete rature model based on multsale mxed nte element method ZHANG QngFu HUANG ZhaoQn * YAO Jun * WANG YueYng & LI Yang Shool o Petroleum Engneerng Chna Unersty o Petroleum (Huadong) Qngdao 66580 Chna * Correspondng authors E-mal: em.group.up@gmal.om; RCOGFR_UPC@16.om Fratured reserors are haraterzed by omplex ratures on multple sales and are qute dult to model. Numeral smulaton o ratured meda s usually done based on dual-porosty model and equalent ontnuum model. Dual-porosty model treats matrx system and rature system as two parallel ontnuous systems oupled by rosslow unton. Ths model s only ald or reserors wth hghly deeloped ratures. Equalent ontnuum model treats ratured meda as a ontnuum meda and equalent absolute permeablty tensors or eah grd blok are alulated to desrbe the heterogenety o the reseror. Ths model s eent only when there exst representate element olume and the equalent permeablty are dult to dede. Both models treat the ratured meda as a smpled model and annot desrbe the multsale lows exatly beause they annot presely onsder the derson eet o the ratures. Although the dsrete rature network (DFN) model an prode a detaled representaton o low haraterst tradtonal numeral method does not sutable or DFN. The major dulty s the sze o the omputaton. A tremendous amount o omputer memory and CPU tme are requred and ths an exeed lmt o today s omputer resoures. Upsalng methods are generally used to redue the omputatonal ost. Howeer t s not possble to hae a pror estmates o the errors that are present when omplex low proesses are nestgated usng oarse models onstruted a smpled settngs. In ths paper multsale mxed nte element method (MsMFEM) s proposed to smulate water/ol two phase low n dsrete rature meda. By ombnng MsMFEM wth the dsrete rature model we am towards a numeral sheme that altates ratured reseror smulaton wthout upsalng. The MsMFEM uses a standard Dary model to approxmate pressure and luxes on a oarse grd. The multsale bass untons are onstruted numerally by solng loal derental equatons on the ne-sale grd. The adantage o MsMFEM s that the bass untons apable o reletng normaton about ratures wthn elements. Thereore ths method an apture the ne-sale eets on the oarse grd that s multsale method an redue the omputatonal ost and keep hgh alulaton auray at the same tme. Tradtonal numeral methods generally dult to deal wth omplex grd element n ths paper mmet nte derene (MFD) method s used to onstrut the multsale bass untons due to ts loal onserateness and applablty o omplex grds. Compared wth tradtonal multsale mxed nte element methods ths method s sutable or arbtrary omplex grd system. Ths paper ntrodued undamental prnples o the multsale mxed nte element method and desrbed the numeral sheme o dsrete rature model based on MsMFEM n detal. Then we dedued dsrete rature model omputng ormulate or the multsale bass unton by usng mmet nte derene method. Oersamplng tehnque s appled to get more aurate small-sale detals. IMPES sheme s used n the two-phase low smulaton. Physal experment s used to proe the aldty the multsale method. The numeral results show that ompared wth tradtonal numeral method the MsMFEM an represent the ne-sale low n rature networks exatly and meanwhle has a hgher omputatonal eeny. ratured meda multsale mxed nte element method dsrete rature model mmet nte derene method two-phase low do: 10.1360/N97016-00584 1401