Математическое моделирование нестационарного электромагнитного поля дефектоскопа обсадных колонн Текст научной статьи по специальности «Физика»

Вычислительные технологии

Том 7, № 6, 2002

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ НЕСТАЦИОНАРНОГО ЭЛЕКТРОМАГНИТНОГО ПОЛЯ ДЕФЕКТОСКОПА ОБСАДНЫХ КОЛОНН *

Э. П. ШУРИНА, А. В. ГЕЛЬБЕР, М. А. ГЕЛЬБЕР Новосибирский государственный технический университет, Россия М. И. Эпов Институт геофизики СО РАН, Новосибирск, Россия e-mail: shurina@fpm.nstu.nsk.su

This work is devoted to modelling of non-stationary electromagnetic fields which arising from defectoscopy of oil-and-gas casings. Modelling is carried out in the terms of vector potential in axial-symmetric case. The method of account of electromagnetic field’s discontinuous on the boundary of different materials is proposed.

Введение

Электромагнитные методы исследования получили широкое распространение в геофизике [1, 2]. Одна из актуальных задач инженерной промысловой геофизики состоит в эффективном контроле за техническим состоянием и толщиной обсадных колонн скважин, подверженных техногенным и геологическим воздействиям. Для ее решения созданы электромагнитные дефектоскопы, практическое применение которых ориентировано на выявление дефектов, обусловленных механическими повреждениями и вызывающих изменение физических свойств материала колонны.

Целью работы является моделирование нестационарных электромагнитных полей, возникающих при работе дефектоскопа нефтегазовых труб, разработанного в Институте геофизики СО РАН и НПП ГА “Луч” [2].

Математической моделью, описывающей поведение электромагнитного поля, является система уравнений Максвелла. В зависимости от типа получаемых дифференциальных уравнений, методы ее решения могут быть условно разделены на три группы [3-6]:

— методы, основанные на решении уравнений первого порядка;

— методы, основанные на введении скалярного и векторного потенциалов;

— методы, основанные на решении уравнений второго порядка (в естественных переменных) .

В первом случае переходят к представлению уравнений в виде закона сохранения. Для их решения применяется метод конечных объемов [7, 8] либо формулируется специальная вариационная постановка [9].

* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант №00-01-00-899, и интеграционной программы СО РАН, грант №2000-1.

© Э.П. Шурина, А. В. Гельбер, М. А. Гельбер, М. И. Эпов, 2002.

Классическим способом решения задач электромагнетизма является использование скалярного или векторного потенциалов [9-11]. Они широко используются при вычислении стационарных и нестационарных электромагнитных полей, а также в задачах с вихревыми токами. Кроме того, потенциалы широко применяются и для решения задач моделирования гармонических полей (системы уравнений Гельмгольца). Аппроксимации с применением потенциалов позволяют избегать ложных или так называемых паразитных решений [9]. Это связано с тем, что ограничения на дивергенцию поля обеспечиваются введением специальных условий калибровки [10, 11]. Кроме того, применение аппарата векторного потенциала позволяет эффективно решать задачи в областях со скачкообразным изменением свойств материалов. Среди преимуществ соответствующих численных процедур следует отметить также относительно низкие вычислительные затраты [12].

Так, одним из наиболее эффективных методов решения задач магнитостатики является использование двух потенциалов: полного - в областях с ферромагнитными материалами, и неполного — в областях, содержащих токовые обмотки [10]. Примером эффективного использования векторного и скалярного потенциалов может также служить работа [11]. Авторами используется потенциально-токовая формулировка с калибровкой Кулона. Для решения полученных уравнений предлагается шахматная (staggered) схема, которая естественным образом учитывает осреднение проводимости внутри ячеек сетки. Также в

[11] исследуются способы решения получаемых систем линейных алгебраических уравнений (СЛАУ).

Альтернативой использованию потенциалов является моделирование задач электромагнетизма в естественных переменных [12-14], т. е. формулировка уравнений относительно напряженностей электрического и магнитного полей. Уравнения Максвелла в этом случае записываются в форме дифференциальных уравнений в частных производных второго порядка. Они могут быть гиперболического или параболического типов в зависимости от физических свойств среды. Несмотря на ряд преимуществ, такие постановки приводят к решению задачи о седловой точке, т. е. к необходимости реализации алгоритмов Удзавы или penalty-метода [15]. Кроме того, возникают некоторые сложности при учете условий непрерывности решения на границах материалов и проводящих границах [9].

Критика постановок в терминах векторных и скалярных потенциалов обычно связана с потерей точности и непрерывности при численном дифференцировании для нахождения векторов напряженности электрического поля и магнитной индукции или поля. Целью моделирования является построение переходной кривой ЭДС E в измерительном контуре, поэтому в численном дифференцировании потенциалов по пространственным переменным нет необходимости.

Для дискретизации по времени мы использовали трехслойную полностью неявную схему. Преимуществом явной схемы [16] является возможность избежать решения СЛАУ на каждом временном слое при специальным образом построенной аппроксимации. Однако требования, предъявляемые к величине шагов дискретизации по времени и пространству, становятся очень жесткими [16, 17]. Для пространственной дискретизации был использован метод конечных элементов (МКЭ) с применением билинейных базисных функций [17].

1. Постановка задачи

1.1. Модель среды

Опишем модель среды, соответствующую реальной ситуации при дефектоскопии обсадных колонн. Обсадная колонна представляет собой полый проводящий цилиндр, его удельная электрическая проводимость составляет примерно (0.77... 5.71) • 107 См/м, а относительная магнитная проницаемость может изменяться в пределах 1... 100 ^0. Толщина стенки колонн редко превосходит 0.2 м при внутреннем ее радиусе 0.05... 0.3 м.

В процессе эксплуатации скважин происходят коррозия обсадной колонны и уменьшение ее толщины. Под влиянием механических напряжений в обсадной колонне могут возникать разнонаправленные “слепые” и проникающие трещины. При соединении отдельных секций обсадной трубы могут оставаться тонкие кольцеобразные области, связанные с их неполным резьбовым соединением. При увеличении механических напряжений вокруг скважины упругие деформации колонн могут переходить в пластические с потерей магнитных свойств. Выявление этих изменений, их идентификация и определение размеров и являются одной из главных целей дефектоскопии.

Область внутри обсадной колонны заполнена либо водонефтяной смесью, либо водой. При лабораторном моделировании это может быть воздух. Эта область не содержит ферромагнитных материалов, и ее удельная электрическая проводимость не может превышать единиц сименс.

Внешняя область представлена прилегающим к колонне квазицилиндрическим слоем цементного камня и горными породами. Эта область заполнена немагнитным материалом с удельной электропроводностью, редко превосходящей единицы сименс.

1.2. Модель зондирующей установки

Рассматриваемый дефектоскоп — это зонд, движущийся по исследуемой скважине. Зонд представляет собой непроводящий цилиндр. На нем соосно расположены две генераторных и одна приемная катушки одинакового радиуса (0.024 м). Последняя находится на равном расстоянии между генераторными катушками. Каждая из катушек представляет собой медный цилиндр (а = 5.71 • 107 См/м, ^ = ^0) с поперечным сечением 6 • 10-4 м и высотами 10-2 м (приемная катушка) и 4• 10-2 м (генераторные катушки). Величина питающего тока в дальнейшем принята равной 1 А. Расстояние между катушками меньше радиуса трубы.

Приемный и генераторные контуры — это реальные катушки с определенными геометрическими размерами. Условной замены контуров магнитными диполями нет. Прибор относится к классу электромагнитных систем, работающих во временной области (ВЭМ-систем) [3-6]. В системах такого типа электромагнитный импульс формируется посредством резкого изменения значения постоянного тока в генераторном контуре. Электромагнитным откликом является ЭДС наведенного тока в измерительном контуре. Одним из достоинств ВЭМ-систем является отсутствие первичного поля при измерениях, что увеличивает их разрешающую способность [3 - 6].

Работу дефектоскопа можно условно разделить на два этапа. На первом этапе по генераторным катушкам протекает сторонний ток ^ = 0. Затем в момент времени £0 этот ток отключается. На втором этапе работы происходит регистрация показаний ЭДС в измерительном контуре в зависимости от времени после выключения тока. При этом

(0, ^0^) 0)Т) £ £ [Мо], 0,£ > £0,

где J0 у = const. Отметим, что в данном устройстве могут быть использованы полный и дифференциальный сигналы, т. е. токи в генераторных контурах могут быть как одинаково, так и противоположно направленными.

1.3. Формулировка уравнений Максвелла в терминах вектор-потенциала

В векторно-дифференциальной форме система уравнений Максвелла имеет вид

dB

— + Vx E = 0, (1)

д D

Vx H = J + —, (2)

V ■ D = Р, (3)

V ■ B = 0, (4)

D = е E, B = ^ H, (5)

где E — напряженность электрического поля; B — магнитная индукция; D — электрическое смещение; H — напряженность магнитного поля; р — плотность зарядов; ^ = ^0^*, ^0 = 4п ■ 10-7 Гн/м, — относительная магнитная проницаемость; е = е0е*, е0 = 8.85 ■ 10-12 Ф/м, е* — относительная диэлектрическая проницаемость.

Плотность тока определяется соотношением

J = J0 + стЕ, (6)

где J0 — плотность стороннего тока; ст — удельная электропроводность.

На границе Г материалов с различными свойствами выполняются следующие условия непрерывности:

Е x n |г+ = Е x n |г- , B ■ n |г+ = B ■ n |г- . (7)

Если один из материалов является проводником, то на границах “проводник — диэлек-

трик” имеют место условия

H x n |г+ — H x n |г- = Jr, D ■ n |r+ — D ■ n |г- = рг . (8)

Здесь n — вектор внешней нормали к границе Г; Jr — плотность поверхностного тока; рг

— плотность поверхностного заряда.

Плотность тока Jr в случае отсутствия сторонних поверхностных токов равна 0. Однако при решении нестационарной задачи в окрестностях границы раздела сред с разными удельными электропроводностями в высокопроводящей области возникает приграничный ток. Толщина слоя, по которому протекает этот ток, мала, и его расчет вызывает большие вычислительные сложности. Поэтому предлагается следующая схема его учета в конечноэлементной аппроксимации — вводится плотность поверхностного тока Jr = ст2 Е h/2, где ст2 — удельная электропроводность стали, h — размер (по г) ячейки конечно-элементной дискретизации на границе Г. Коэффициент 1/2 связан с типом используемых конечных элементов.

Для задачи (1) - (6) граничные условия задаются следующим образом:

E |дп — 0, H |дп — 0, B |дп — 0, D |дп — 0,

где дП — Липшиц-непрерывная граница области П. Размеры области определяется условиями “большого бака”.

Введем векторный потенциал А:

В = Ух А.

(9)

Подставляя (9) в (1), (2) и учитывая (6), получаем в П:

-2А -А ^ т

£ -^2 + ^~дї + У Х У Х =

(10)

В осесимметричном случае, не нарушая общности, можно переформулировать задачу в цилиндрических координатах [14]. Тогда

А = (0,Ау, 0)т, Ау = Ау(г,г),

и уравнение (10) принимает вид

д2А р----

д£2 ' ~ д£

В осесимметричном случае для гг-координат граничные условия имеют вид

|дП = 0.

-А 1

у + - У(^-1УАу) + —Ау = ^.

г2^

(11)

(12)

Отметим, что для осесимметричных задач соотношение (4) выполняется автоматически. Расчетная область П в гг-координатах схематически показана на рис. 1.

2. Этапы моделирования

Первый этап работы дефектоскопа (£ < £0) можно описать стационарными уравнениями Максвелла ( ^ А? = 0, ^ у =0 |:

-І2

-І

-У^" :УАу) + Ау = Зоу

-1

(13)

с граничными условиями (12). Особенностью является учет условий (7), (8) на границах Г12 или Г23 для сред с разными магнитными свойствами

Н х п

+

- Н х п

Г12

Г12

= 12

В ■ п

+

Г12

В ■ п

Г12

Здесь J12 — поверхностный ток. В стационарном режиме ^2 = 0, следовательно,

+ +

Нх п =Нхп 2 , В ■ п = В ■ п 2

Г12 Г12 Г12 Г12

(14)

(15)

в терминах вектор-потенциала

+

^і (Ух А) х п = (Ух А) х п

(Ух А) ■ п

Г12

+

Г12

Г12

(Ух А) ■ п

Г12

Границы Г12 и Г23 — коаксиально-цилиндрические (г = г12, г = г23, г € [а, Ь]), и в связи с этим в гг-координатах (15) преобразуется к виду

+

і1

^1

іА + -Ау

Ау + ГУ

г -г

-А,

-г

Г12

+

Г12

і1

^2

-А

-г

1 а + -Ау

А у + Г) г -г

Г12

(16)

(17)

Второй этап работы прибора (£ > £0) описывается уравнением (11) с Т0- = 0 и граничными (12) и начальными условиями

|+ * = -

у I ¿=¿0

-А у '

ір)

0.

¿=¿0

На границах Г 2 и Г23 для нестационарной задачи условия сопряжения имеют вид

+

і А + -Ау

г у+ -г

Г12

-А,,

-г

+

Г12

г

-А,

1А + -Ау

у+ -г

Г12

-г

(18)

(19)

(20) (21)

Г12

Для определения Т1 2 в (20) необходимо вычислить объемную плотность наведенного тока 1 в стальной трубе

Т = Е = дА-

1 = а2Е - = - °2—о^,

по которой в дальнейшем вычисляется плотность поверхностного тока Т12 и Т23. В измерительном контуре рассчитываем плотность наведенного тока

(22)

-1 д£

Для вычисления наведенной ЭДС в измерительной катушке использовался закон Ома.

Г

2.1. Дискретизация по времени нестационарной задачи

Дискретизация по времени построена для трехслойной полностью неявной схемы [9, 10]. Для (11) на п-м слое (* = *п) получаем

А- - 2А”;-1 + А^-2 А- - А^-1 1

р---------Д2---— + - - д* - - V(^-1VA-) + —А- = 0) (23)

где Д* — шаг по времени. Аналоги начальных условий (18), (19):

Ап-1 = Ап-2 = А А- А- А--Для каждого временного слоя получаем уравнение

11 р д*2 + —д

А” - V(м-1VA”) + А” = Р (А”"1, А“"2), (24)

где

2 Ап"1 Ап"2 Ап"1

Р(А“-1, А“-2) = р “ д *2 “ + —-Ду. (25)

3. Вариационные формулировки

Введем гильбертовы пространства Ь2(П), Н 1(П), Нд(П), где

Н 1(П) = {V € Ь2(П) | ^ € Ь2(П)3},

Н1(П) = {V € Н 1(П) | V = 0 дП}.

Сформулируем вариационную постановку для стационарной задачи.

Для заданной плотности тока Т0 “ найти А” € Щ, такое, что для любого V € Н1 выполняется:

У ^^А^^П + У ^ А-^П + У А-^5 -у ^ А-^5 +

П П Г12 Г12

+ У А-^5 - I А-^5 =у То-V ^П. (26)

Г23 Г23 п

Для нестационарной задачи вариационная постановка на каждом временном слое п имеет вид:

Для заданной функции Р(А--1, А--1) найти А-, такое, что для любого V € Н выполняется:

11 £——г + — ——

д*2 д*

+ У А-^5 - I А-^5- I Тг12 V ^5 - I Тг23 ^5 +

Г12 Г12 Г12 Г23

+ У А-^5 - I А-^5 = I F(A“-1)A“-2)vdП) (27)

Г23 Г23 п

где Р(А- 1, А- 1) определяется (25).

П

3.1. Построение дискретных аналогов вариационных формулировок

Введем подпространство Vн С Н0 билинейных базисных функций, определенных на прямоугольной сетке П = и Пе. Представим приближенное решение А- задач (26), (27) в

е -

виде линейной комбинации базисных функций:

N

АУ = £ Ау,^,,

,= 1

где N — число узлов дискретизации конечно-элементной сетки. Тогда дискретный аналог формулировки (26):

Найти А- € Vй, такое, что для любого ^ € Vh выполняется

I ^-1УАу -У^к dП + J ^Ау^к dП + J ^Ау^к -^ ^Ау^ +

П П Г12 Г12

+ У ^Ау^к dS - I ^Ау^к dS = I dП, (28)

Г23 Г 23 П

а дискретный аналог формулировки (27) —

Найти А-™, такое, что для любого ^ ^ выполняется

11 є——- + о ——

Аі2 Аі

Ауп^к dS + / ^-1УАупУ^к dП+ / ^Ауп^к dП +

П

П

+ [ ^ Ауп^к dS - / ^ Ауп^к dS- / ІГ12 ^к dS - / ^к dS +

г

г

Г12 +

Г12

Г12

^2-Ауп^ dS - / ^Ауп^к dS = / Р(А^\А^2)^к dП.

Г23

\}т— 1 А]т—2\

г

г

Г23

Г23

П

(29)

П

Дискретные аналоги (28), (29) приводят к системам линейных алгебраических уравнений. Матрицы СЛАУ для задач (28) и (29) имеют разреженную структуру, их портреты идентичны, если расчеты выполняются на одной и той же сетке.

3.2. Вычислительная схема

1. Решение стационарной задачи (вариационная формулировка (26)).

2. На каждом шаге по времени:

а) расчет наведенных токов Т12 и Т23 в соответствии с (22);

б) решение нестационарной задачи (вариационная формулировка (27));

в) расчет наведенного в измерительном контуре тока;

г) расчет наведенной в измерительном контуре ЭДС.

4. Результаты вычислительного эксперимента

Дефектоскоп представляет собой зонд, движущийся вдоль трубы. Определим понятие “центр зонда”. В рассматриваемом случае центр зонда находится на оси симметрии трубы (г = 0), поэтому в дальнейшем мы будем указывать только его z-координату. Центр зонда по оси z — положение центра измерительной катушки.

Тестирование. Для исследования точности предложенной численной схемы проведено сравнение расчетов, полученных на двух вложенных сетках (7560 и 29 887 узлов). Параметры среды следующие: внутренний радиус трубы без дефектов 0.073 м, внешний

Компонента Нг на двух сетках. Приведено сечение г = 0 вблизи трубы

г, м Hr, А/м (7560) Hr, А/м (29 887) Абсолютная разница, А/м Относительная разность

6.10e - 02 1.7486e - 02 1.7496e - 02 -1.0552e - 05 -6.0349e - 04

6.50e - 02 1.4943e - 02 1.4953e - 02 -1.0185e - 05 -6.8161e - 04

6.99e - 02 1.2922e - 02 1.2931e - 02 -9.3892e - 06 -7.2663e - 04

7.25e - 02 1.2259e - 02 1.2268e - 02 -8.9952e - 06 -7.3376e - 04

7.30e - 02 1.2801e - 04 1.2810e - 04 -9.3954e - 08 -7.3395e - 04

7.35e - 02 1.1696e - 04 1.1704e - 04 -8.5873e - 08 -7.3423e - 04

7.40e - 02 1.0606e - 04 1.0614e - 04 -7.7902e - 08 -7.3452e - 04

7.50e - 02 8.4692e - 05 8.4754e - 05 -6.2277e - 08 -7.3534e - 04

7.60e - 02 6.3840e - 05 6.3887e - 05 -4.6976e - 08 -7.3584e - 04

7.70e - 02 4.3434e - 05 4.3466e - 05 -3.1932e - 08 -7.3519e - 04

7.80e - 02 2.3409e - 05 2.3426e - 05 -1.7127e - 08 -7.3165e - 04

7.90e - 02 3.5166e - 04 3.5190e - 04 -2.4111e - 07 -6.8562e - 04

8.20e - 02 2.8441e - 04 2.8460e - 04 -1.9220e - 07 -6.7579e - 04

8.60e - 02 2.0694e - 04 2.0708e - 04 -1.3710e - 07 -6.6253e - 04

8.84e - 02 1.6570e - 04 1.6581e - 04 -1.0860e - 07 -6.5539e - 04

z.

0.004

0.002

0

-0.002

-0.004

0.02 0.025 0.03 Г’М

Рис. 2. Фрагмент конечно-элементной сетки, содержащей 29 887 узлов вблизи зонда.

радиус трубы без дефектов 0.079 м, = 95, а = 0.77 См/м. Вычислялся полный сигнал от двух генераторных контуров, центр зонда z = 0. Результаты расчета горизонтальной компоненты стационарного магнитного поля приведены в таблице. Фрагмент сетки, содержащей 29 887 узлов, представлен на рис. 2. Пространственное распределение компонент нестационарного магнитного поля приведено на рис. 3.

Исследование решения в зависимости от наличия дефектов в трубе. Абсолютная разность ЭДС, наведенной в измерительном контуре в различные моменты времени для труб с дефектами и бездефектной трубы, приведено на рис. 4. Использовался полный сигнал генераторных контуров. В данной работе в качестве дефектов рассматриваются кольцевые углубления на внутренней поверхности трубы. В трубе имеют место различные

Я,-ЮЛ A/M Д,-ЮЛ A/M

а б

Рис. 3. Компоненты магнитного поля Нг (а) и Нг (б) на вложенных сетках (7560 узлов — сплошная линия и 29 887 узлов — линия с маркерами). Сечение г = 0, Ь = 0.51 мс.

ЕЛОг6, В

Рис. 4. Абсолютная разность полных ЭДС для дефектной и бездефектной труб.

комбинации следующих геометрических дефектов (тип дефекта указан на рисунке):

— дефект 1 0.073 м < г < 0.0735 м -0.001 м < г < 0.001 м;

— дефект 2 0.073 м < г < 0.0735 м -0.015 м < г < -0.013 м;

— дефект 3 0.073 м < г < 0.0735 м -0.018 м < г < -0.016 м;

— дефект 4 0.073 м < г < 0.0735 м 0.014 м < г < 0.016 м;

— дефект 5 0.073 м < г < 0.0735 м 0.0165 м < г < 0.0185 м;

— дефект 6 0.073 м < г < 0.0735 м 0.019 м < г < 0.021 м.

Анализ влияния соединительной муфты. Проведено исследование кривой ЭДС, наведенной в измерительном контуре, получено распределение плотности вихревого тока в трубе в различные моменты времени при наличии соединительной муфты. Соединительная муфта представляет собой полый цилиндр, внутренний радиус которого равен внешнему радиусу соединяемых труб (внутренний радиус трубы 0.065 м, внешний радиус трубы 0.073 м). Параметры муфты: 0.073 < г < 0.083 м, -0.09 < г < 0.09 м, р,м = 95^0, а = 0.77 ■ 107 См/ м. В данном случае предполагается, что между соединяемыми трубами есть зазор. Параметры зазора: 0.064 < г < 0.073 м, -0.001 < г < 0.001 м. Режим функционирования установки — измерение полного и дифференциального сигналов. Рассматривались несколько положений зонда. Результаты приведены на рис. 5-10.

Отметим, что основными источниками измеряемого сигнала являются вихревые токи, текущие в трубе. Поэтому характеристики пространственного разрешения зондирующей системы непосредственно связаны с временной эволюцией этих токов в проводящем материале. Пространственно-временное распределение вихревых токов показано на рис. 8, 9.

ЕЛО-*, В

ЕЛЪ*, В о

1\

-2 -

-4 -

-6 -

К'

/ /

У /

/у /' /!

/ ' //

-• / //

\ ' //

\\ / // 2 ~ ®

•Л / / і ------------------------------- 0.01 м

' \ / / ,

\ V - ' / / ---------- 0.02 м

\

\\ \

'л\ \

Л \

0.04 м 0.05 м

Т—I—I—I—|—I—I—I—I—|—I—I—I—I—|—I—I—I—I—|—I—I—I—I

0.2 0.4 0.6 0.8 МО“3, с

Рис. 5. Абсолютная разность полных (а) и дифференциальных (б) ЭДС для труб с муфтой при

различных положениях зонда.

а

Сравнение с результатами физического эксперимента. Проведено сравнение вычисленных и физически измеренных значений ЭДС. Параметры установки следующие: внутренний радиус трубы 0.073 м, внешний радиус трубы 0.079 м, = 95^0, труба без дефектов, полный сигнал генераторных контуров, центр зонда г = 0. Результаты приведены на рис. 10.

ЛЮ2, А/м2

Л-Ю2, А/м2

Рис. 6. Плотность полного наведенного тока (сечение г = 0.065 м) в разные моменты времени: а

— труба с муфтой, б — труба без муфты; центр зонда г = 0.02 м.

а

Л-Ю2, А/м2

А

0.21 мс \\ \\

0.51 мс \\

0.99 мс V ' \ V

/ \\\

\\\

¡1 /' -ч\\

А4 '

\\\ /7

\\ Ч

- \\ \у \\ 1

\\ 1 \

_ V/

-0.1 0 0.1 г,

а

Л-102, А/м2

Рис. 7. Плотность дифференциального наведенного тока Л (сечение г = 0.065 м) в разные моменты времени: а — труба без муфты, б — труба с муфтой; центр зонда г = 0.02 м.

г-10-2, м

0.5

г-10"2, м

Рис. 8. Распределение плотности полного тока в трубе с муфтой: центр зонда 0.02 м;

Ь = 0.21 мс (а) и Ь = 0.51 мс (б).

а

аб Рис. 9. Распределение плотности дифференциального тока в трубе с муфтой: центр зонда 0.02 м; Ь = 0.21 мс (а) и Ь = 0.51 мс (б).

Выводы по результатам вычислительного эксперимента. Сравнение решений, полученных на вложенных сетках, показало приемлемую точность предложенных вычислительных схем. Это говорит об их применимости для моделирования электромагнитных полей при дефектоскопии нефтегазовых труб.

Сравнение ЭДС, рассчитанных для дефектной и бездефектной труб (см. рис. 4), позволяет сделать вывод о том, что наличие большего числа дефектов увеличивает абсолютную разность ЭДС, но сохраняет характер зависимости, увеличивая лишь ее максимальные и минимальные значения.

£10ЛВ

Рис. 10. Расчетная ЭДС с учетом реальных размеров катушек (сплошная линия), экспериментальная ЭДС (пунктирная линия) и расчетная ЭДС без учета реальных размеров катушек (штрихпунктирная линия) для трубы без дефектов.

Исследование влияния муфты (см. рис. 5) на поведение ЭДС показали, что при использовании полного сигнала генераторных контуров максимальное влияние муфты имеет место в положении центра зонда г = 0, а затем, по мере удаления центра зонда от зазора между соединяемыми трубами, убывает. При дифференциальном сигнале генераторных контуров при положении центра зонда г = 0 в силу симметричности конструкции влияние муфты минимально. Максимальное влияние муфты отмечается при некотором небольшом смещении центра зонда от зазора, и это влияние убывает по мере дальнейшего смещения центра зонда. Совместное использование полного и дифференциального сигналов генераторных контуров дает возможность исключить влияние крепежных зазоров муфты на результаты дефектоскопии скважин.

Заключение

В работе предложена вычислительная схема для моделирования электромагнитных полей при дефектоскопии нефтегазовых скважин.

Моделирование проводилось в осесимметричном случае с использованием векторного потенциала. Применение векторного потенциала в данном случае оправдано тем, что основной задачей моделирования было определение значений ЭДС, которое не требует операции пространственного дифференцирования, и использование векторного потенциала оказывается наиболее экономичным.

Первый этап моделирования состоит в решении стационарной задачи при наличии постоянного тока в генераторных контурах, второй — в решении нестационарной задачи после отключения тока в контурах. Для дискретизации по пространству использован метод конечных элементов на билинейных базисных функциях. Предложен способ учета скачков магнитного поля на границах материалов с различными физическими свойствами. Трехслойная полностью неявная схема применена для дискретизации по времени.

Проведен ряд вычислительных экспериментов. Показана хорошая точность предложенной схемы. Сравнение с физическим экспериментом показало приемлемое моделирование поведения кривой ЭДС для исследования влияния особенностей трубы. Проведен анализ влияния групп дефектов и соединительной муфты на поведение кривой ЭДС. Анализ проводился как для полного, так и для дифференциального сигналов генераторных контуров.

В дальнейшем планируется построение аналогичных вычислительных схем для моделирования электромагнитных полей в трехмерном случае.

Список литературы

[1] ДВОРЕЦКИЙ П. И., ЯРМАХОВ И. Г. Электромагнитные и гидродинамические методы при освоении нефтегазовых месторождений. М.: Недра, 1998. 318 с.

[2] ТЕХНОЛОГИЯ исследования нефтегазовых скважин на основе ВИКИЗ: Метод. руководство / Под ред. М. И. Эпова, Ю. Н. Антонова. Новосибирск: НИЦ ОИГГМ СО РАН, Изд-во СО РАН, 2000. 121 с.

[3] NEWMAN G. A., Hoversten G.M. Solution strategies for two- and three-dimensional electromagnetic inverse problems // Inverse Probl. 2000 №27. P. 1357-1381.

[4] BOURGEOUS B., SuiGNARD K., PERRUSSON G. Electric and magnetic dipoles for geometric interpretation of three-component electromagnetic data in geophysics // Inverse Probl. 2000. №27. P. 1209-1225.

[5] Dorn O., Bertete-Aguirre H., Berryman J. G., Papanicolaou G. C. A nonlinear inversion method for 3D electromagnetic imaging using adjoint fields // Inverse Probl. 1999. №15. P. 1523-1558.

[6] Hjelt S.-E., PiRTTIJARVI M. Some characteristics of the conducting plate model in the inversion of geophysical electromagnetic data // Inverse Probl. 2000. №27. P. 1209-1225.

[7] ClONI J.-P., FEZOUI L., ISSAUTIER D. High order upwind schemes for solving timedomain Maxwell equations //La Recherche Aeropatiale. 1994. №5. P. 319-328.

[8] Munz с.-D., Schneider R., Voss U. A Finite-volume Method for the Maxwell Equations in the Time Domain. Forschungszentrum Karlsruhe Technik und Umwelt Interner Bericht, 51.03.01, 1996.

[9] JiANG B., Wu J., POVINELLI L. A. The origin of spurious solutions in computational electromagnetics // J. Comput. Phys. 1996. Vol. 125. P. 104-123.

[10] ШУРИНА Э. П., СОЛОВЕЙЧИК Ю. Г., РОЯК М. Э. Решение трехмерных нелинейных магнитостатических задач с использованием двух потенциалов. Новосибирск, 1996 (Препр. СО РАН; ВЦ. №1070).

[11] HABER E., Ascher U.M. ET AL., Fast Simulation of 3D Electromagnetic Problems Using Potentials // J. Comput. Phys. 2000. Vol. №163. P. 150-171.

[12] ASSOUS F., DEGOND P., HEINTZE E., ET AL. On a finite-element method for solving the three-dimensional Maxwell equations // J. Comput. Phys. 1993. Vol. 109. P. 222-237.

[13] ASSOUS F., DEGOND P., SEGRE J. Numerical approximation of the Maxwell equations in inhomogeneous media by a conforming finite element method // J. Comput. Phys. 1996. Vol. 128. P. 363-380.

[14] LACOSTE P. Solution of Maxwell equation in axisymmetric geometry by Fourier series decomposition and by use of conforming finite element // Numer. Math. 2000. Vol. 84. P. 577-609.

[15] COULLIETTE D.L., Koch M. On the difficulties and remedies in enforcing the div=0 condition in the finite element analisys of thermal plumes with strongly temperature-dependent Viscosity // Intern. J. Num. Methods in Fluids. 1994. Vol. 18. P. 189-214.

[16] САМАРСКИЙ А. А., ГУЛИН А. В. Численные методы. М.: Наука, 1989. 432 с.

[17] РОЯК М. Э., СОЛОВЕЙЧИК Ю. Г., ШУРИНА Э. П. Сеточные методы решения краевых задач математической физики. Новосибирск: Изд-во НГТУ, 1998. 120 с.

Поступила в редакцию 17 июня 2002 г., в переработанном виде — 23 июля 2002 г.