CC BY
44
ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
УДК 536.2:519.6.56 В. И. Васильев, В. В. Попов
ПАРАЛЛЕЛЬНАЯ РЕАЛИЗАЦИЯ ГИДРОДИНАМИЧЕСКОЙ МОДЕЛИ НЕФТЯНЫХ МЕСТОРОЖДЕНИЙ
В работе численно реализована математическая модель процесса добычи нефти. Использован итерационный метод сопряженных градиентов с применением попеременно-треугольного метода, и разработан алгоритм распараллеливания.
Результаты вычислительного эксперимента, проведенные на вычислительном кластере университета, показали эффективность распараллеливания.
Ключевые слова: модель, алгоритм, метод, попеременно-треугольный, масса, неразрывность, проницаемость, вычисление, кластер, параллельный.
V. I. Vasilyev, Ш V. Popov
The parallel realization of oil fields hydrodynamic models
In this article it is numerically realized the mathematical model of oil production. This includes iterative method of conjugate gradients in connection with alternate-triangular method; it is developed parallelization algorithm.
Results of computational experiments were carried out on the University supercomputer cluster and presented the effectiveness of parallelization.
Key words: model, algorithm, method, alternate-triangular, mass, continuity, permeability, scaling, cluster, parallel.
Математическое моделирование процесса разработки нефтяного месторождения с помощью вторичного метода добычи нефти посредством вытеснения из пористого коллектора одной жидкости другой жидкостью является достаточно наукоемкой проблемой. Наиболее адекватные модели можно получить, базируясь на фундаментальных законах механики гетерогенных сплошных сред [1]. Математическая модель процесса добычи нефти численно реализована с помощью итерационного метода сопряженных градиентов с
предобуславливателем Якоби [2]. В настоящей работе численная реализация получающихся на каждом временном слое систем разностных уравнений для давления проводится с помощью итерационного метода сопряженных градиентов с использованием попеременно-треугольного метода, и разработан алгоритм, позволяющий эффективное
распараллеливание.
ВАСИЛЬЕВ Василий Иванович - д.ф.-м.н., профессор, первый проректор по учебной и научной работе СВФУ.
E-mail: vasilMttnail.ru
ПОПОВ Василий Васильевич - к.ф.-м.н.. директор СКЦ ИМИ СВФУ.
E-mail: pwya3c@rami5ler.ru
Пусть нефтяная залежь представляет собой область произвольной геометрии Q,. и через нагнетательные скважины с координатами (xj,yi) закачивается вода с
дебитами qi, / = 1, Ж, ас эксплуатационных скважин с
координатами (х;, yi ) и дебитами qj , i =Nx+t, N2
добывается нефть. Приведем формулировку задачи из статьи [2].
Уравнение неразрывности для воды с учетом метода фиктивных областей имеет вид [3, 4]:
Э sp Nl
т —-1 = ciiv(p. v:) + У q,8{x - х,, у - у,),
Vt 2=1
(х,у)е Ql5r>0, (1)
Э(1-5)/7т ^ чо,
т—---------^ = div{p2v2)+ > ^<)(.v-.v.T v-;>). (2)
Уравнение движения представим в виде обобщенного закона Дарси
v, = -к ^ grad р, i = 1,2 . (3)
Mi
При этом учтем как слабую сжимаемость флюидов, так и упругую деформируемость коллектора [5]
Pi = Р-.>Л\+{р- р,.)! К,).
(4)
Суммируя уравнения (1) и (2) относительно функции давления, получаем невырождающееся уравнение параболического типа
т1 =т0(\ + (р-р0)/Кс), / = 1,2. (5)
Поскольку константы , Кс являются величинами
порядка 10' Па> то пренебрегая квадратичными слагаемыми, имеем
тр, = »> Р , (1 + (Р~Рд! ^),
где
^ =
Г
т,
s 1-5 +
VK, К2 j
др Э ( др
Э t дх\ Эх
Э
+ —
(6)
Здесь
f1 + 1 ' -1
V Ki + KC J , /=1,2, (7) 4
р = т(_
s 1-5
ЧК7 + ”КГ у
(8)
Подставив полученное выражение в уравнение неразрывности потока воды (1), имеем
ds
т---------------------1- т,
К, Э t Эх \ Эх
+ ■
Э t
_э ?
ЭДГ1 Эуу ,=1
(хл'Н Ол >().
+
\kf(s)/julAx,y)eQ, £2,(x,v)e пх \ Q.
где
э
+ —
(
\ n2
где
3>Ч ~ ду j т+1
(х,у)е £l1,t> 0,
,(xv)gQ,
Фг ~' /-\(л-. г). 11 \а.
. ф-т ^ й .
$7 “Г (х, v)e Г2,,/>0.
kfx(s)jfjj +kf2(s)/ju2, (x,v)e О, гг, (х, v)e ц \а Начальное условие
/?(x,j,0) = /?0, (х. г)С О,. Граничное условие
P(X,y,t) = Рс' (-Y. v)c dLi-.l >0.
(13)
(14)
(15)
(16)
(?)
Здесь, используя метод фиктивных областей, имеем следующее продолжение коэффициента
Уравнение (1) целесообразно с использованием функции Баклея-Леверетта и, пренебрегая малыми слагаемыми, привести к виду, также допускающему эффективную численную реализацию
m^=J&<!£+'*P'!£\ ф,.0 е Q.J \ От / > $, (17) Э/ \дх дх дуду)
где
<р--
л
Pi
(10)
, , р, =—, ^2 =
/1+//0/2 Pi
Dj = (х,у : 0 < х < /j, 0 < у < 12), е - достаточно малое действительное число, /Jx - коэффициент динамической вязкости воды, p(x,v,t) - давление, s(x,v,t) - водонасыщенность, с>(х-х,,у- V, ) -дельта-функция Дирака.
Аналогично записывается уравнение неразрывности для нефти
ds \-s др Э f ^ др)
- т(,----V тс,--------= — ф. •— +
0 Эt К, Эt Эх Г' Эх
1=^+1
В начальный момент времени распределение водонасыщенности задается в виде
л(х._гт0) = х(!. (18)
Вода в нагнетательные скважины закачивается с максимально допустимой водонасыщенностью 5 , т. е. граничные условия задаются в виде
я{х,у,() = !К. (X,зО€ О,./>0. ф,у,Г) = яо
(х. V)С: ЭЦ. ,/>(). (19)
Закон сохранения энергии с учетом подвижности флюидов записывается в виде уравнения теплопроводности с конвективными членами
cpK=±(jK)+jL
Г7Л'
(П)
(12)
Э t Эх у Эх) ду ^ ду
др дТ др дТ
+ Y-------+ —,
Эх Эх ду ду
(х,у)е 12-, / > 0.
+
(20)
Здесь
Ср = (1 - m)Csps + msCwpw + m{ 1 - s)C.opo Л = (1 - m)As + ms/I + m(l - s)A0,
соответственно, объемная теплоемкость и коэффициент неявную разностную схему на квадратной
теплопроводности системы пористая среда - вода - пространственной сетке, использующую поправочный
нефть, нижние индексы соответствуют коэффициент 1](х,у) в окрестности каждого
коллектору, воде и нефти; точечного источника (стока)
(Одвю/д +**%<»>//«,.<*.#• а Н-д-/ .мьг„„у-
Ш = { 0 у у I ,2
У [ Ог(х,_у)еЦШ, ^ Н
ч Р: - Р: 1 ; ,
| Я, (Х,у)бО Ф^112.}71(х^1/2^У]) ^0
_|- ф Т}\Х у )——
Распределение температуры в начальный момент уг^+1 12 2 5 ]+112 ^2
времени
Т(х,у,0) = Г0, (х,^)е Ц. (21) а г, ; *
Г раничное условие
Т(х,уА) = Тс,(х,у)еда1,1>0. (22) + ^чЛ^^ъу1~у1)
Для зависимости коэффициентов динамической
Г
N.
1=1
вязкости от температуры следует использовать / = 1,«1 — 1, у = 1,и, — 1. (23)
эмпирические формулы Сеточные аналоги граничных условий первого рода
1/(7 ) ~ 0.001 имеют вид
(0.11622 + 0.005(7’-273.15))- ри = р,й,7 = 0,, у = 0,и2 -1;/ = 0,
(24)
- формула В. Я. Булыгина [6] для нефти;
п2,1=\,п] -1;
//(У) =-----------------------^= %, / = 0. «1 - 7 :0./7л 1;0,
(1+0.033 X?1 - 273.15) + 0.00022 Ц7 - 273.15)) --------- (25)
%»г=1,»! -1;
- формула Пуазейля для воды, //п - коэффициент х = у / = 0 и / = 0 и - Г / = 0
^ 1 ’ 2 ’ (26)
динамической вязкости воды при температуре 273,15 К. Точное решение стационарной двумерной задачи о
%, / = 1Л} -1.
скважине имеет вид: Систему линейных алгебраических уравнений
Ч
относительно распределения давления (23) на каждом
1 = ре 4--------1п(г / Гс), временном слое целесообразно решать с помощью
итерационного метода сопряженных градиентов с где а - дебит скважины. использованием попеременно-треугольного метода [8].
Для повышения точности аппроксимации уравнения Распределение ^ во до насыщенности определяем с фильтрации можно пользоваться специальными помощью явной монотонной разностной схемы с разностными соотношениями, не требующими явного аппроксимацией конвективного члена против потока
задания контуров скважин. Пусть в области £21
5* —5*
введена квадратная сетка 0)ь с шагом /?. —-----— = —-
+ (7 Ф ~( *> . . а - а
'..У ‘,.1
„ , , t„ 2А 2 /?
Введем поправочные коэффициенты по аналогии с я
тем, что предложено А. Н. Чекалиным [7] Ь . — ю. \ т. —т . Ь — го. .1 т ,—ф
I I г_Л_ I ' 1' из~{ | г*3 | 1 ' *,./+1
/? 2 /?
71(х-хк,у-ук) =
0.5л:/\п(Ь/гс),(х-хк)2 +(у-ук)2 =1/4, У = 1,7?1-1, у = 1,7?2-1, (27)
1п(2.5)/1п(2), (х-хх.)2 +{у-ук)2 =5/4, где
1,(Х-%)2+0;-^)2 >=9/4Д=1,Ж . , ь = . Ри44 - А,;-1
"7 2/7 ’ ; 1 2/7
Численную реализацию построенной модели _____ ______
процесса разработки нефтяного месторождения I = 1,^ -1, ] = 1, »2 - 1.
осуществляем с учетом поправочных коэффициентов в Дискретные аналоги граничных и начальных условий
следующем порядке. Сначала уравнение (13) на каждом % _ у (х г )с £). (?8)
временном слое поставим в соответствие чисто
= 50, 7=1 ,п2,*и0 = хо,і = 0,пг,
*и = я0,і = 0,пг, ] = 0,п2,п = 0.
(29)
(30)
В начально-краевой задаче для уравнения теплопроводности (20)-(22) поставим в соответствие чисто неявную разностную схему с дискретизацией конвективных членов " против потока". Для численного решения полученной системы нелинейных алгебраических уравнений на каждом временном слое применяем метод простых итераций, каждый шаг которого осуществляется методом сопряженных градиентов.
Линеаризованную систему алгебраических уравнений представим в виде.
Лу = /.\щ
Первый шаг итерационного процесса осуществляется методом скорейшего спуска
вУк±1 ~ У к +Аук = уд = 0 Тк+1
Последующие шаги - методом сопряженных градиентов
ВУк+1 = ак+1 ~ Ткл\А)Ук + (1 _ ак+1Тк+\/’ к = .. ,
где
= }7к’Гк\’ к = Ш-, Гк = АУк Щ = Гк -
{Лм>к,м>к)
\ . ^к+х
1
4-1
Ри - Р: :
5 1-5
Чк7+1<7У
А / \ Р'+1-} Р-.І
- Фм1хМхмп>У /)------—- +
+ Фі-и2
Ри -Р.-и
, , . Рі.і
■ Фи+1 пМХі, / 2 )-----~2------ +
— Щ і У $ -У і)
1=1
где ? - легк0 обратимые операторы. Обращение
оператора В выполняется в два этапа. На первом этапе нужно решить уравнение
(Е + 0)11 )№* = /; .
На втором этапе - уравнение
Р + =%■.
Для первого этапа формулы для вычислений имеют вид
II,/) = СПУ( (/ -1,У) + 1%(/,7-1)+сгд.(/,./)
Где
а--
і = 1-я -1,./ = 1,«2 -!•
^ш?7(хм/2,.У,)
Г^-1 /2, А-112 ,Уі) + Фй і-ї / 2?7(*,, 7 /-1 / 2 ) +’
1г +60
+ т,
5 1-5
ЛІ12
уК, к2у
&:=-
Г Фи11% і■ / 2, V і ) + Фи 1-М 2П(Х:, V _!, 2 ) +'
+ О)
+ т,
• }' ' Г,,]
{ 5 1 — 5^1
Кг+ К2 )іп
с = -
{Фі■-*/ 2,/2,3;, ) + Ф, /-1 / Л/2) +
Таким образом, реализация метода сопряженных градиентов осуществляется в следующем порядке:
Ух = >7, -^1^0-А*-! = °%+1А + (1 — )Ук-1 — %+] М>к г
+ <»
+ /77,
V
к, + к
2 /
Для второго этапа:
Щ (/, /) = ста* (/ +1, ]) + Ът1 (/, /+1) + стк (/,}).
где
а
■ = 1, ;?1 1. / : : 1, ;?т - 1 .
^■+1/2,
Г ФыпЖхттУь) + & +
+ И-Ї,
' П,]-+1/
5 1-5Л
,К: + к2 у
6=-
фш / 2,у / 2, ;•/) + ФітїзЛ&і ’ у т 12) +
+ й>
+ /77,-
-2
7=1,77! -1, | = 1,% -1.
Для попеременно-треугольного метода оператор 5 представляется в виде произведения операторов В = (Е + (ОІІ )(Е + 0)Я2).
с =-
Й+1 / 2,/ /2 ^ ) + $,У+1 /2-^, ) +
+ Си
+ тг
I 5 1-5^
- + -
чК, К2У>Я
40000
30000
2» 20000 -
10000 =------------------------
о г I' I' I' ' I' 111 ' I 'I |1 I |1 |1 I I
О 10000 20000 30000 40000
X
Рис. 1. Графическая интерпретация конвейера
Для вычисления у?к (/, ]), i — \щ— 1, 7 = 1,п2 — 1 требуются его значения слева и ниже от нее. Если на внешнем цикле продвигаться вдоль полосы, а на внутреннем - поперек, то получится конвейер (рис. 1).
Вычисления ^(7,7) проводятся по аналогичному конвейеру справа налево и сверху вниз. На внутреннем цикле - поперек полосы, на внешнем - вдоль.
Результаты расчетов проведены при следующих значениях исходных данных:
р0 =20 МПа; /их = 10_3 Па*сек; //2 =2-10”' Па*сек; 5-0 =0.2; Т0 = 10°С ; т = 0.2; к = 1(Г12 мг; = 1.84 Дж/м3*К;
(\Р, =4.19 Дж/м3-К;
Соро = 0.9 Дж/м3,К; Л, =2.09 Вт/м«К; Ли. = 0.58 Вт/м*К; Ла = 0.5 Вт/м*К; в момент времени 10 лет, на пространственной сетке с количеством узлов 4097x4097,
Размеры области - 40 км. х: 40 км, Эксплуатационные и нагнетательные скважины расположены в шахматном порядке. Дебиты нагнетательных скважин одинаковы и равны 10"4 м/сек, дебиты эксплуатационных скважин равны 0.8"4 м3/сек. на 1 метр мощности коллектора.
На рис. 2-4 показаны распределения соответственно давления, температуры и водонасыщенности. В качестве языка программирования выбран язык С, а для реализации обмена между ядрами кластера использован интерфейс передачи данных МР1. На рис. 5 и 6 представлены результаты счета, полученные с помощью 32 ядер. Они показывают эффективность распараллеливания попеременно-треугольного метода.
Рис. 2. Распределение давления
Рис. 3. Распределение водонасыщенности
Рис. 4. Распределение температуры
I
600-
400 -
200-
°0 8 16 24 32 Р
Рис. 5. Зависимость времени расчета от количества ядер
и
24-
16-
8-
°0 "®Г 16 ' ' 24 32 Р
Рис. 6. Зависимость ускорения от количества ядер
Литература
1. Нигматулин Р. И. Динамика многофазных сред / Р. И. Нигматулин. - М. : Наука, 1987. - Ч. 2. - 360 с.
2. Васильев В. И. Специализированные симуляторы гидродинамического моделирования месторождений / В. И. Васильев, В. В. Попов, М. В. Васильева // Вестник ЦКР Роснедра. -2008. -№4. - С. 59-63.
3. Коновалов А. Н. Задача фильтрации многофазной несжимаемой жидкости / А. Н. Коновалов. - Новосибирск: Наука, 1988.-166 с.
4. Вабищевич П. Н. Метод фиктивных областей в задачах ма тематической физики / П, Н. Вабищевич. - М.: Изд. МГУ, 1991. -156 с.
5. Чарный И. А. Подземная гидрогазодинамика / И. А. Чарный. - М. : Гостоптехиздат, 1963. - 396 с.
6. Булыгин В. Я. Гидромеханика нефтяного пласта / В. Я. Булыгин. - М. : Недра, 1978. -232 с.
7. Чекалин А. Н. Численные решения задач фильтрации в водонефтяных пластах / А. Н. Чекалин. - Казань: Изд-во КГУ. 1982.-207 с.
8. Самарский А. А. Методы решения сеточных уравнений / А. А. Самарский, Е. С. Николаев. - М. : Наука, 1978. - 562 с.
УДК 53.043
П. В. Винокуров, Е. П. Неустроев, С. А. Смагулова, И. В. Антонова, М. С. Каган
ПЕРЕЗАРЯДКА ЛОКАЛИЗОВАННЫХ СОСТОЯНИЙ В СТРУКТУРАХ БДОЮЕ/БГ С КВАНТОВЫМИ ЯМАМИ В ЗАВИСИМОСТИ ОТ СОДЕРЖАНИЯ СхЕ
Методом зарядовой спектроскопии глубоких уровней исследованы напряженные структуры 81/810е/81 с квантовыми ямами. Определены энергии активации эмиссии носителей заряда с локализованных состояний в исследуемых структурах. Увеличение содержания германия в сплаве вЮе и, соответственно, рост механических напряжений в яме сопровождаются формированием дефектов, обладающих малыми временами перезарядки. В результате в структурах, в которых дислокации и другие крупные структурные нарушения еще не обнаруживаются, процессы обмена носителями между квантовой ямой и окружающей ее матрицей кардинально меняются. Уровни в квантовой яме проявляются только при низких температурах Т 120К'.
Ключевые слова: кремний, германий, гетероструктуры, квантовые ямы, метод зарядовой спектроскопии глубоких уровней, глубокие уровни, дефекты, напряженные структуры, локализованные состояния, перезарядка уровней.
P. V Vinokurov, E P. Neustroev, S. A. Smagulova, I. V. Antonova, M. S. Kagan
Recharging of localized states in Si/SiGe/Si structures with quantum wells in relation to Ge level
There were studied compressive structures of SiGe quantum wells by using charge deep level transient spectroscopy. Experimental data were used to determine activation energies of carrier’s emission from localized states. Increase of germanium level in SiGe alloy and growth of mechanical tensions in pit entailed formation of structural defects that had a short recharging time. Therefore, structures with low Ge content and without large-scale defects exhibited different behavior. Eevels in a quantum well were observed only at low temperatures T<120K.
Key words: Silicon, germanium, heterostructure, quantum wells, charge deep level transient spectroscopy, deep levels, defects, strength structures, confinement levels, recharging levels.
ВИНОКУРОВ Павел Васильевич - аспирант ФТИ СВФУ E-mail: yokkorashoMHnail :ru
НЕУСТРОЕВ Ефим Петрович - к.ф.-м.н., доцент ФТИ СВФУ
E-mail: nefittl(a>yandex.ru
СЛ-1АГУЛОВА Светлана Афанасьевна - к.ф.-м.н., доцент ФТИ СВФУ E-mail: smagulova@mail. ru
АНТОНОВА Ирина Вениалшновна - д.ф.-м.н.. Институт физики полупроводников им. А. В. Ржанова СО РАН E-mail: antonova@isp.nsc.ru
КАГАН Мирон Соломонович - д.ф.-м.н.. Институт радиотехники и электроники им. В. А. Котельникова РАН E-mail: kagan@,mail.cplire.ru
