CC BY
27
УДК 519.63
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ГЕОЭЛЕКТРИЧЕСКИХ ПОЛЕЙ В КУСОЧНО-АНИЗОТРОПНЫХ КВАЗИФРАКТАЛЬНЫХ СРЕДАХ ЖЮЛИА
P.P. Яматов, В.Н. Кризский
MATHEMATICAL MODELING OF GEOELECTRICAL FIELDS IN PIECEWISE ANISOTROPIC QUASIFRACTAL JULIA MEDIA
R.R. Yamatov, V.N. Krizsky
Рассмариваются алгоритмы компьютерного моделирования полей точечных источников постоянного электрического тока в кусочноанизотропных квазифрактальных средах Жюлиа, описывающих рудные и пористые нефтегазонасыщенные среды. На основе вариационных алгоритмов А.Н.Тихонова строятся процедуры решения обратных задач по определению параметров квазифрактальных сред.
Ключевые слова: квазифракталъные пористые кусочно-
анизотропные среды Жюлиа, поле постоянного тока, прямая и обратная задача.
The computer modeling algorithms for direct electric pointed source fields in piecewise quasifractal Julia media were build. Julia quasifractals are describe meaning and porous oil and gas media. On the base of A.N.Tikhonov variation type algorithms the solving procedures for inverse problems of determination of quasifractal media parameters are build.
Keywords: quasifractal porous piecewise anisotropic Julia media, direct
Введение
Характеристики нефтегазонасыщенных систем, представленных пористыми или трещиноватыми средами, в существенной мере определяются хаотическим распределением зерен породы, капилляров и трещин по форме и размерам. Как известно, пористые среды — среды фрактальной структуры. Важной эксплуатационной характеристикой подобных систем является коэффициент проницаемости, который можно рассчитать, если знать распределение поровых пустот (коэффициент пористости и связность поровых каналов), т.е. если известна структура среды, отвечающая стохастическому распределению поровых пустот, капилляров и трещинных каналов.
Большинство встречающихся в природе фрактальных структур являются квазифракталами, поскольку на некотором малом масштабе фрактальность исчезает. Квазифрактал отличается от идеальных абстрактных фракталов конечностью, неполнотой и неточностью повторений структуры. В силу этого, для расчета физических полей в квазифрактальных средах, представимых вкраплениями одной среды в другую (заполненных флюидом пор
в скелете, зерен одного вещества в другом и т.н.) могут быть использованы классические алгоритмы, учитывающие большое количество мелких анизотропных включений, генерируемых процедурами, реализующими при построении принцип фрактальности. В этой связи актуальной представляется проблема разработки алгоритмов решения прямых и обратных задач компьютерного моделирования геофизических (электрических, магнитных, тепловых, диффузионных и т.д.) полей в кусочно-анизотропных средах с квазифрактальными включениями.
Ранее в работах [1 - 3] на основе комбинации методов интегральных преобразований и интегральных уравнений были построены алгоритмы расчета стационарных и нестационарных физических полей в трехмерных кусочно-однородных средах с включениями. На их основе рассматривались обратные задачи определения границ и удельной проводимости включений в классе простых тел [4], тел вращения [5] и бесконечных цилиндрических тел [6, 7] с аппроксимированной сплайном границей. В [8, 9] результаты обобщены и распространены на кусочно-анизотропные среды.
В данной работе рассмотрены алгоритмы компьютерного моделирования стационарных электрических полей в квазифрактальных средах.
Для решения прямой задачи используются методы интегральных представлений с построением функции Грина вмещающего пространства и интегральных уравнений. На основе методов вариационного типа получены решения некорректной обратной задачи поиска параметров, положения и формы квазифрактального анизотропного тела.
Разработанные алгоритмы реализованы в комплексе программ и допускают распараллеливание при использовании суперкомпьютеров, многопроцессорных вычислительных комплексов или вычислительных кластеров.
1. Моделирование квазифрактальных сред Жюлиа
Рассмотрим множества Жюлиа (73^), которые в кубе стороны 4 с центром в начале координат, при делении его на N 3 элементарных объемов, образуют квазифрактальные включения в виде совокупности элементарных объектов. В качестве элементарного объекта среды удобно рассматривать «куб> со сглаженными углами и ребрами, с гранями, нарал-лельными координатным плоскостям и длинами ребер, много меньшими размеров всего фрактала, геоэлектрические свойства среды которого описываются симметричным постоянным тензором удельной электрической проводимости.
Для построения множества 73И рассмотрим последовательности гиперкомплексных чисел Ют(к), построенных для любого гиперкомплексного числа к по правилу:
Юо = к,'Шт+1 = Pnl('Wm)/Qn2('Wm) = 7 (ют )•
Здесь Рп 1 и <3п25 т>п2- гиперкомплексные полиномы степени п\ и П2 соответственно, а рациональное конформное отображение </(ги) - эндоморфизм.
Множеством Жюлиа 7(к) называется [10] замыкание множества точек притяжения бесконечности - таких гиперкомплексных чисел к, что последовательности уот{к) сходится к бесконечному пределу, при т стремящемся к бесконечности. Множество Жюлиа на гипер-комплексной плоскости х,у, г самоподобно. Множества Жюлиа делятся на два основных класса: связные и вполне несвязные. Во втором случае множество состоит из несчетного числа дискретных точек(«:ныль Фату> ).
Обозначим через <73-0 = С\7(к) - дополнение множества Жюлиа до гиперкомплексной плоскости. Множество 73И имеет конечное число компонент связности. Как множество 7(к), так и множество 73^, с точки зрения геофизики, интересны своими геологическими аналогами.
В частном наиболее исследованном случае, Pnl(u>)/Qn2(u>) = и>2 + 1, где 1 = а + гЬ + jc, а, Ь, с £ Я - произвольные гиперкомплексные числа. Множество Жюлиа 7(А’) и его дополнение — </31)(а, Ь, с) , при конечном т, будут квазифрактальными параметрическими множествами, зависящими от трех параметров а, Ь и с.
Квазифракталы, сгенерированные на основе 73^-алгоритма, обладают рядом свойств:
1. при а = Ъ = с = 0 множество <73-0(0, 0, 0) - единичный шар;
2. множество 13П(а, 0, 0) симметрично относительно осей у и г и при малых а — есть тело вращения относительно оси х;
3. если |(/| = л/а2 + Ь2 + с2 > 2 и х2 + у2 + г2 > д, то любая траектория отображения 7(х,у,г) стремится к бесконечности.
Перечисленные свойства позволяют оптимизировать алгоритм построения квазифракталь-ных включений в случае симметрии, ограничить вариации а, Ь и с, провести серию сравнительных вычислительных экспериментов при а = Ь = с = 0 в классе простых тел.
Рис. 1. Множества ,]ЗИ и их аксиальные сечения при некоторых значениях параметров а, Ь и с
На рис. 1 демонстрируются множества 73Р при различных значениях параметров а, Ь и с; Рп(ю) = Ю2 + 1.
Осуществляя перенос, растяжение/сжатие и поворот построенного множества, можно получать различные аналоги геологических объектов, подлежащих моделированию: «классические> включения (шар, линза, тело вращения); локальные включения, окруженные ореолом породы (зоной перемешивания сред); пористые среды.
При программной реализации алгоритм моделирования 73Р-квазифрактальных включений сталкивается с рядом проблем:
1. большое количество (103 — 104) элементарных включений, что ведет к росту размерности решаемых систем линейных алгебраических уравнений и необходимости учета мнимых границ при соприкосновении элементарных кубов при расчете физического поля;
2. большая глубина рекурсии используемых функций (большое число итераций т).
Все это влечет значительные затраты машинного времени на генерирование квазифракталь-ных включений и на расчеты физических полей в среде при их наличии. Для ускорения на этом этапе эффективно распараллеливание процедуры генерации.
2. Поле точечного источника постоянного тока
в кусочно-анизотропных квазифрактальных средах
Рассмотрим горизонтально-слоистую среду, разделенную гладкими параметрически заданными границами 71.0, 72.0>---> 7п.о на горизонтальные слои Г2о.о> ^1.о> ^2.о> • • •> ^п.о с удельными электрическими проводимостями, заданными симметричными тензорами (То.сь (Ті.о, 02.сь • • ч Оп.0) соответственно. Пусть каждый слой о содержит кі локальных квазифрактальных включений с границами 7^ и постоянными симметричными тензорами проводимости аj = 1 ,к{.
Пусть в точке А горизонтального слоя Г2т.о находится точечный источник постоянного электрического тока, с которого стекает ток силы I.
В произвольной точке Р ф А и{Р) - потенциал среды электрического поля, создаваемого источником постоянного тока, описывается следующей краевой задачей:
(1)
(Цу(а^Уи^(Р)) = -/¿.¿(Р),Р є = 0,n,j = 0,кі,
{о’і.а^и'і.а(Р), п) - ((п.^и^(Р),1у)\^_ = 0, £4о(Р) - ^(Р)!^., = 0,і = 0,п,і = 1 ,кі} (2)
4-і
((7г.0УС/г.0(Р),п) - ((Т*_1.0УС/*_1.0(Р), п) |^_ = 0, £Л.0(Р) - С/*-1.0(Р) п = 0,1 = 1, П, (3)
(4)
гДе /¿.і = { “umлl, 6(Р, А) - функция Дирака, п - вектор нормали, (2), (3)
17г.о " ’ ' > 17г.0
иг.о(Р) ^ 0,Р(х,у, г) ^ ж,г = 0,п,
_ ЩР,А),Р€Пт. о,
^ \0,Р£Пт.0,
- условия непрерывности потенциала и плотности тока на границах контактов сред, (4) -условие регулярности решения на бесконечности.
Для решения задачи воспользуемся методом интегральных представлений, поэтапно понижая геометрическую сложность задачи [3, 8, 9].
Рассмотрим вспомогательную задачу для функции Грина вмещающего пространства (без квазифрактальных включений в последнем слое) с точечным источником единичной интенсивности, находящемся в произвольной точке Q(хq,уд, ):
(Іт^і^С^ЛР,(5)) = -¿(Р, <5),Р Є = 0,n,j = 0,кі
(5)
КоУС^0(р,д),п) - (аг.^сыр,(Э), п)1 ., = о, сир, д) - с?.(р,д)I , = о, (6)
'7г.;/
'7г.;/
г = 0,п - 1,^ = 1, кг,
{<Ti.oVG2.oiP,Я),п) - (^_1.оУС?_1.0(р,д),п)|74о = о, С£0(р,д) - с?_1.о(р,д)|740 = о, (?)
г = 1, п,
Щ.о(Р) ^ 0,Р(х,у, г) ^ ж,г = 0,п.
Интегральное представление решения задачи (1) - (4) имеет вид: р
и(р) = Е и((^-0 - (Тп.з)УС^(Р, д), пд) с11п.]д + Сп(р, А),
«/'Уг, *
(8)
(9)
^’=1 7>г.]
из которого следует, что решение исходной задачи (1) - (4) может быть получено в любой точке Р кусочно-анизотронной области, если будет определено решение задачи (5) - (8) -функция Грина Сп (Р^), и будут известны граничные значения потенциала на границах подобластей, не вошедших в задачу для функции Грина. Здесь п вектор внешней нормали в точке Q.
Опуская в (9) точку Р на каждую из границ 7пл, I = 1, кп и учитывая условия (2), получим систему интегральных уравнений Фредгольма второго рода относительно неизвестных граничных значений потенциала вида: р
ипЛ{Р)-Е ип.,{я) (к. о - <7П.д-)УС^(рд), пд) = СП(Р, А), Ре ъ,1- (10)
j=l ^7п.3'
Для решения задачи (5) - (8) снова применим метод интегральных представлений и интегральных уравнений. Рассмотрим вспомогательную задачу для функции Грина вмещающего пространства, считая включением последний слой:
¿^{а^УС^ЛР, О)) = -¿(Р, д),Р е Г^,г = 0,п- 1,^' = 0, кг
{аг.оЧ&1о{р, д), п) - {Ог.^сир, д), п) = о, С^0(Р, д) - СЫР, д) = о, (12)
7г.;
7г.;
г = 0,П - 1,^ = 1, кг,
К0ус£0(р, д), п) - (^_1.0ус?_1.о(р, д), п)
= 0,
7г. 0
^0(р,д)-с?_1.0(р,д)
= 0, г = 1, п — 1,
7г. 0
(11)
(13)
<Зг.о(Р, д) -> о, Р{х, у, г) ->• 00, г = 0, п - 1. Интегральное представление решения задачи (5) - (8) будет иметь вид:
(14)
сп{р,я) = сп{я 1,д) (к_1.о -<7П.о)УСп(р,д1),пд1] +сп{р,я), (15)
*1п.О
а значения Сп(Р, (5) на границе 7п.о ~ определяются из интегрального уравнения
' "Уп.О
СП(Р,Я) - сп(дьд) (<7„_1.0 - <7„.0)усп(р, до, пд1 £г7п.0д1 = сп(р,д), (16)
где Р е 7п.0.
Таким образом, переход от исходной задачи (1) - (4) к аналогичной задаче (11) - (14) упрощает геометрию области решения - в среде стало на один слой с включениями меньше.
К задаче (11) - (14), следовательно, снова можно применить вышеописанный способ, вплоть до среды, состоящей из одного анизотропного СЛОЯ По.О- Если при этом тензор его удельной электрической проводимости задан диагональной матрицей <то.о = 011 0 0 \
О (Т22 0 , то решение последней задачи можно построить аналитически:
О 0 озз /
-1
с0(р, Я) =
і
47Г^(Тіі(Т22<7зз
Схр ~ хч)2 (ур - учУ
(У И
<722
+
(%> + ^)2
<^33
где Р(хр,ур, гр), Я(хд,уд, гд) - точки, содержащие приемник и источник тока.
В табл. 1 приведены результаты сопоставления значений аномального потенциала постоянного электрического тока системы двух точечных источников силы /х = +1 А ,/2 = —1 А с координатами Ах(+40, +40, 0) м, Аг(—40, —40, 0) м, соответственно, вычисленного но алгоритму, изложенному выше, при различном порядке N детализации квазифрактала 7ЗИ(0, 0, 0) и шара с центром в 0(0, 0, 6) м, радиуса И = 2 м и удельной проводимостью <70.1 = 0,02 См/м. Система приемников расположена на площадке Д = {1 £ [—20, 20], у £ [—20, 20], г = 0)} «дневной» поверхности однородного полупространства с удельной электрической проводимостью (То.О = 0, 01 См/м.
Таблица 1
Сопоставление значений аномального потенциала в однородном полупространстве и]Ы в присутствии квазифрактала /3-0(0, 0, 0) и 11зрг - шара, - количество экспериментальных данных
N 1 \Uspr-Uj3d\i (У/_ N ' ' Л 5 ^ і*гг . м г=1 1 Ыи 0Е(^рт.-ад2, в тах \Uapr - и^і, В гєі,^ц
20 5,13 0,00017722 3,18е-06
30 1,29 0,00006944 1,03е-06
40 1,00 0,00006117 8,89е-07
50 0,69 0,00006047 5,61е-07
Относительная погрешность вычислений составила менее 6 %, что свидетельствует о достаточной хорошей точности расчетов. Увеличение N количества элементарных «кубов», входящих в квазифрактальное множество, влечет повышение точности расчетов.
На рис. 2а приведены изолинии аномального потенциала постоянного электрического тока системы источников и приемников рассмотренной для табл. 1, в случае однородного полупространства с удельной электрической проводимостью <то.о = 0, 01 См/м в присутствии квазифрактальных однородных включений 7ЗИ(а,Ь, с), построенных в шаре радиуса 2 м с центром О = (0, 0,—4) м и удельной электрической проводимостью <то.1 = 0,1 См/м, для значений параметров а, Ь и с, приведенных на рис. 1. На рис. 2б приведены аномальные поля, когда те же квазифрактальные включения 7ЗИ(а, Ь, с) в однородном полупространстве анизотропны, с ненулевыми коэффициентами удельной электрической проводимости вдоль осей системы координат 0, 2, 0, 05 и 0, і См/м соответственно.
а=-0.7120, Ь=-0.1300, с=0.1280 а=0.3660, Ь=0.0860, с=0.3000 а=-0.3820, Ь=0.5960, с=-0.1120
а - квазифрактальные включения - однородны
а=-0.7120, Ь=-0.1300, с=0.1280 а=0.3660, Ь=0.0860, с=0.3000 а=-0.3820, Ь=0.5960, с=-0.1120
б - квазифрактальные включения - анизотропны
Рис. 2. Изолинии аномального потенциала постоянного электрического тока, В
3. Восстановление параметров квазифрактальных включений
Задача определения границ и параметров квазифрактальных включений по наблюдаемым значениям потенциала электрического тока относится к классу обратных задач электроразведки. Вследствие неединственности и неустойчивости решения задача некорректна.
Будем искать параметры квазифрактальных включений, как экстремальные решения регуляризирующего функционала А.Н. Тихонова.
Зададим форму і.]-го квазифрактального включения
а
•г-0
векторами
0,.;{ Н,.г(),.г а,.г Ь,.г с;.г Д',.;). где = { іі11. И!- г ^) - вектор масштабирования с
коэффициентами сжатия/растяжения куба стороны 4 м по осям х, у и г соответственно. Оі.і = I ()!' і. (■.()] і і - смещенная от точки (0,0,0) координата центра, г/,.;. с,.; - пара-
метры JЗD функции, N-1^ - параметр детализации фрактала, определяющий количество элементарных «кубов».
Введем В рассмотрение конечномерный вектор 5 = (фол, Фо.2, • • • ) фі,і, • • • ) Фп.кп ) і определяющий вид всех квазифрактальных включений среды. Будем искать его как решение, минимизирующее регуляризирующий функционал А.Н. Тихонова вида:
Р(з) = \\ие-ит(з)\\ЫВхГ)) +
й1»*
(17)
где 1/е, ит - соответственно экспериментальные и модельные (как решение задачи (1) - (4)) значения потенциалов между датчиками, расположенными в узлах сеточного множества Т) х И приемников/источников ноля, а — параметр регуляризации, 5° — А’-компонентный вектор (описывающий к включений) опорной модели, который строится с учетом всей априорно известной информации о структуре исследуемой среды (например, по данным сейсмо- и/или гравиразведки).
Ограничивая вариацию компонент конечномерного вектора в, получим компактное множество корректности А.Н. Тихонова, на котором существует единственное квазирешение, определяющее форму квазифрактальных включений.
Поиск экстремали регуляризирующего функционала (искомые параметры включений) проводится методом Хука - Дживса поиска минимума сильно-овражной функции.
В табл. 2 приведены результаты восстановления параметров квазифрак-тального включения </3-0( —0,3820; 0,5960; —0,1120), смещенного в точку Оод =
(1,1230;—1, 3453;—6, 2343) м относительно начала координат, с вектором масштабирования До.1 = (2, 2340; 2, 2340; 2, 2340), по значению аномального потенциала постоянного электрического тока от системы источников силы /1 = /3 = +1А, /2 = /4 = —1А с координатами А\ = (40,40,0) м, А2 = (-40,-40,0) м, А% = (40,0,0) м, А4 = (—30,0,0) м соответственно, на площадке В = {х £ [—20, 20], у £ [—20, 20], г = 0)} «дневной» поверхности однородного полупространства с удельной электрической проводимостью <то.о = 0, 01 См/м. Удельная электрическая проводимость включения 00.1 = 0,1 См/м в изотропном и с ненулевыми коэффициентами проводимости вдоль осей системы координат 0, 2 См/м, 0, 05 См/м и 0,1 См/м, соответственно, анизотропном случае. Параметр детализации N0.1 = 40. В качестве априорно известной информации рассматривались значения: Оо.г = (0; 0; —10), До.1 = (3, 3, 3), а = 0, Ь = 0, с = 0.
Таблица 2
Восстановление параметров квазифрактального включения при различном количестве источников
Система источников о ох м о о М о0.і, м Ко. і &о.і Ьо.і Со.і Улзв,%
1,123 -1,345 -6,234 2,234 -0,382 0,596 -0,112 4,352
включение однородное
1,225 -1,468 -6,550 1,975 -0,376 0,590 -0,120 8,458
5, % 9,08 9,12 5,06 11,59 1,57 1,02 7,14 94,36
¿1, ¿2 1,075 -1,413 -6,450 2,375 -0,420 0,575 -0,124 6,965
5, % 4,27 4,99 3,46 6,31 9,95 3,52 10,71 60,05
¿Ъ ¿2, ¿3, ¿4 1,144 -1,325 -6,325 2,275 -0,420 0,564 -0,120 4,470
5, % 1,87 1,51 1,45 1,84 4,71 1,85 7,14 2,80
включение анизотропное
¿1 1,288 -1,125 -6,650 2,025 -0,420 0,564 -0,140 7,640
5, % 14,65 16,38 6,67 9,36 9,95 5,37 25,00 75,55
¿Ъ ¿2 1,250 -1,475 -6,525 2,425 -0,378 0,600 -0,120 7,594
5, % 11,31 9,64 4,66 8,55 1,05 0,67 7,14 74,50
¿Ъ ¿2, ¿3, ¿4 1,155 -1,225 -6,375 2,344 -0,390 0,592 -0,116 4,986
5, % 2,85 8,94 2,26 4,91 2,09 0,67 3,57 14,59
В табл. 2 е - относительная погрешность, значение = -^шэ/А^о 1 х 100% - аналог
коэффициента пористости, где Л^/зд - число элементарных объектов, образующих включение. Для шара V]зд = 53, 846, при параметре детализации N = 40.
Вычислительные эксперименты показали, что, наиболее эффективной стратегией поиска решения задачи определения параметров квазифрактального включения является следующая: на нервом шаге определяется местоположение (координаты центра), на втором —
вектор масштабирования, на третьем - выполняется уточнение коэффициентов генератора фрактала. Данной стратегией определялась последовательность варьирования параметров в методе Хука - Дживса. Найденные параметры квазифрактала позволяют строить оценки коэффициента пористости среды (см. последний столбец табл. 2), а последующие исследования квазифрактала на связность (наличие капиллярных каналов) и флюидодинамику позволят оценить коэффициент проницаемости нефтегазонасыщенных систем.
Литература
1. Иванов, В.Т. Методы решения прямых и обратных задач электрокартажа / В.Т. Иванов,
2. Кризский, В.Н. Математическое моделирование потенциальных геоэлектрических по-
s. Кризский, В.Н. О способе вычисления геофизических полей в кусочно-однородных сре-
4. Герасимов, И.А. Математическое моделирование геоэлектрических полей в осесимметричных кусочно-однородных средах: дис. ... канд. физ.-мат. наук / И.А. Герасимов. — Стерлитамак, 2004.
5. Викторов, С.В. Математическое моделирование геоэлектрических полей в осесимметричных средах со сплайн-аппроксимацией границ: дис. ... канд. физ.-мат. наук /
6. Беляева, М.Б. Математическое моделирование электрических полей в цилиндрических кусочно-однородных средах со сплайн-аппроксимацией границ: дис. ... канд. физ.-мат.
7. Трегубов, Н.В. Программно-алгоритмическое обеспечение навигации бурения горизонтальных скважин / Н.В. Трегубов, В.Н. Кризский // Системы управления и информа-
8. Кризский, В.Н. О способе вычисления физических полей в кусочно-анизотропных средах. Ч. I. Стационарные ноля / В.Н. Кризский // Вестн. Башк. гос. ун-та. - 2009. -
9. Кризский, В.Н. О способе вычисления физических полей в кусочно-анизотропных средах. Ч. II. Нестационарные ноля / В.Н. Кризский // Вестн. Башк. гос. ун-та. - 2009, -
10. Морозов, АД. Введение в теорию фракталов / АД. Морозов. - М.; Ижевск: Ин-т компьютер, исследований, 2002. - 160 с.
Кризский Владимир Николаевич, доктор физико-математических наук, профессор, кафедра математического моделирования, Стерлитамакская государственная педагогическая академия, Krizsky@ rambler .ru.
Яматов Рим Рафикович, кафедра математического моделирования, Стерлитамакская государственная педагогическая академия, [email protected].
Поступила в редакцию 15 апреля 2011 г.
