УДК 539.3, 550.3, 517.968.28, 517.956.224
ИСПОЛЬЗОВАНИЕ МЕТОДА ГРАНИЧНЫХ ЭЛЕМЕНТОВ
ДЛЯ МОДЕЛИРОВАНИЯ ОТРАЖЕНИЯ ОТ ШЕРОХОВАТЫХ ГРАНИЦ
Егор Борисович Сибиряков
Институт нефтегазовой геологии и геофизики им. А. А. Трофимука СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Коптюга, 3, кандидат физико-математических наук, старший научный сотрудник, тел. (383)330-90-02, e-mail: [email protected]
В данной работе использовался модифицированный метод граничных элементов для решения краевых задач третьего рода в случае упругих стационарных колебаний. В широком диапазоне частот вычислялись перемещения двухслойной среды под действием сосредоточенного поверхностного источника в случае гладкой, а также малоамплитудной шероховатой границы. Показано, что сигнал, порожденный осесимметричным источником, после отражения от шероховатой границы может получить угловые компоненты, которые сравнительно медленно затухают при удалении от источника.
Ключевые слова: метод граничных элементов, упругие стационарные колебания, шероховатые границы.
METHOD OF BOUNDARY ELEMENTS FOR MODELING OF THE ROUGH INTERFACE REFLECTION
Egor B. Sibiriakov
Trofimuk Institute of Petroleum Geology and Geophysics SB RAS, 630090, Russia, Novosibirsk, Koptyug Prospect 3, Ph. D., Senior Researcher, tel. (383)330-90-02, e-mail: [email protected]
A new method of boundary value problem solution for the elastic stationary oscillation was used for calculation of the full displacement field caused by two-layered medium with rough and smooth interfaces in wide enough frequency domain. It was established the signal caused by the axially symmetrical source could contain transvers components after the reflection from the rough interface. Mentioned component caused by diffraction from rough surface only did not attenuate in some frequency domain and can be detected far from the source.
Key words: method of boundary elements, rough surface, elastic stationary oscillation.
ВВЕДЕНИЕ
Главная задача этой работы - найти не только количественные, но и качественные различия полного поля перемещений при отражении от гладкой и шероховатой поверхностей. Проблема отражения от подобных поверхностей интересна для сейсмологии (для понимания природы глубоко залегающих границ), сейсморазведки (выявления эродированных границ), а также для акустики (NDE, выявление дефектов неразрушающими методами).
Лучевое разложение в случае малоамплитудной шероховатости дает очень большие погрешности. Аналитические методы [1], использующие метод малых возмущений (как правило, метод Рэлея, Кирхгофа или их модификации), предполагают предельно простую геометрию поверхности раздела. В случае, если
вектор нормали к поверхности изменяется достаточно быстро, они дают большие ошибки. Использование метода конечных элементов для подобных границ раздела крайне затруднительно, поскольку триангуляция подобных поверхностей становится отдельной проблемой [2]. Кроме того, число обусловленности после проведения триангуляции становится настолько большим, что численное решение нельзя считать достоверным даже после проведения регуляризацион-ных процедур.
В задачах, где имеется как сосредоточенный поверхностный источник (вектор нагрузок изменяется очень быстро), так и малоамплитудная шероховатая граница раздела (быстро изменяется вектор нормали), лучше всего использовать модифицированный метод граничных элементов [3].
ПОСТРОЕНИЕ НОВЫХ ЯДЕР
Суть метода граничных интегральных уравнений в том, что решение в любой точке объема ищется в виде интеграла по поверхности:
При этом точка х0 фиксирована, точка х бежит по поверхности и по ней проводится интегрирование. Тензор М1к(х0,х) удовлетворяет уравнению упругих стационарных колебаний автоматически. Вектор Р^х) вычисляется (подбирается) таким образом, чтобы удовлетворить граничным условиям. Обозначения. Есть фиксированная и бегущая точка поверхности, а также вектор, проведенный из первой во вторую. На любой поверхности есть одно выделенное направление - вектор внешней нормали. По этой причине в качестве координаты хг = (г,п) выбирается скалярное произведение этого вектора на нормаль к поверхности в бегущей точке. Также на поверхности есть два взаимно ортогональных касательных направления: х1,х2. Обозначим х2 = (г,^), х3 =
= (г, т2), гт = у/х% + х\, х2 = гтС05ф, х3 = г^соБф. Лучшие с точки зрения обусловленности ядра для статических задач дают тензор Грина для полупространства [3]. На первый взгляд, в случае стационарных колебаний также нужно использовать тензор Грина для полупространства с использованием цилиндрических координат. Однако решение достаточно простой задачи Лэмба о нормальном ударе по полупространству в случае стационарных колебаний имеет серьезные недостатки. Упомянутое решение имеет вид [4]:
(1)
00
ип — J /оС^г^т)
(2кг2 ехр(—|л:1|у5) — (2к? — к2) ехр(—|*1|ур)) с1к.
г
О
00
и г = 11'о(кггт)
2п\1К(кг, к)
(—2у5ур ехр(—1*1^5) 4- (2к% ~ ехр(— (1к
г
О
2 _
к2 = —, = — к2, у.р = д/^г _ У2к2, кг есть параметр интегрирования,
Я = 4к2^к$ — к2-^^ — у2к2 — (2к2 — к2)2 - знаменатель Рэлея, у2 =
/о(^ггт) - функция Бесселя первого рода. Недостатки. Первый - интервал интегрирования бесконечен (потеря точности). Второй - в некоторой точке знаменатель Рэлея обращается в ноль (еще большая потеря точности). Первый этап модернизации ядер - конечный предел интегрирования. В итоге это даст феноменальное повышение точности. Второй этап заключается в изменении представления дельта-функции, т. е. (2) есть отклик по перемещениям на нормальную нагрузку:
8(5) = -¿.С К! о (кггх)(1кг (3)
Вместо (3) лучше использовать смещенный конечный аналог:
8^5)
N
= ~2кг1° + 2к2) &кг (4)
о
Это дает возможность изменить (2) и получить два из пяти членов матрицы ядер.
N
кгу2
и, ' * 4 г 2
П = I/оОз^т) ехРН*1К) ~ (2кг + зк2) ехр(-|*11^2))
о
N
Г к V
иг = ] /оОз^т) 2^1(1 к) (~2УгУ2 ехРС—+ (2к2 + Зк2) ехр(-|*1|у2)) &к,
о ' ___
V! = у/к2 + к2, у2=у/к2 + к2(2-у21 у3 = 4к2 + 2к2,
_ _
= 4(к2 + 2к2)у/к2 + к^к2 +к2(2- у2) - (2к2 + Зк2)2 (6)
Тензор М1к(х0,х) имеет только пять независимых ненулевых компонент в координатах, связанных с направлениями п,гх и ср. Два из них (Мппи Мгп) есть (5). Чтобы получить остальные три, необходимо решить задачу о касательном ударе по полупространству. Для этого нужно повторить изящные вычисления [5], принимая во внимание различие между (4) и (3). Последнее усовершенствование заключается в замене интеграла (1) конечной суммой. Таким образом, решение ищется в виде и^х0) = М1к(х0, х)¥к(х). Вектор нагрузок можно вычислить аналитически с помощью закон Гука (Р[(х0) = щкпок)\ р^х0) = = —Т,Рцс(х0,х)¥к(х). Все необходимые проекции М^(л:0,л:) могут быть вычислены как компоненты вектора после вращения. Например, ип0 = = ип(п0,п) + иГг(п0,еГт) + и^щ, еф).
Е Мп0п (х0, х) Рп (х) + Мп0т1 (х0, х) Е^ (х) + Мп0т2 (х0, х) Рт2 (х) (7)
ФОРМУЛИРОВКА ЗАДАЧ
Верхняя свободная поверхность - плоскость (7=0). Граница раздела - два варианта: плоскость (7= - 0.5) либо шероховатая поверхность. Ее уравнение есть z = —0.5 + 0.01 * (со5(10лх) + со5(107гу)) (рис.1). Параметр к изменяется от 0.1 до 20, что соответствует эффективной длине Б-волны отот 20п до Упомянутая частотная область покрывает диапазон от квазистатики до лучевой динамики. Нагрузка на верхней поверхности осесимметрична и имеет только нормальную компоненту, равную рп = ехр (—ЮООхг2). гЕ [0,1]. Упругие константы верхнего слоя: X = ц = д = 1. Нижний слой имеет такие же параметры, за исключением \х~ = 1.2. Условия на границе раздела - непрерывность вектора перемещений и нулевая векторная сумма вектора нагрузок. Использовались цилиндрические координаты с равномерной сеткой 25*25 на каждой из поверхностей. На верхней поверхности вычислялся вектор перемещений, порожденный нагрузкой и двухслойной средой при различных частотах к.
Рис. 1. Шероховатая поверхность раздела. г = —0.5 + 0.01 * (собСЮях) + +со5(1071у)), х и у изменяются от 0 до 1
РЕЗУЛЬТАТЫ
Отражение осесимметричного сигнала от шероховатой поверхности дает угловую компоненту на всех частотах. Однако на достаточно коротких волнах амплитуда этой компоненты слишком мала по сравнению с поверхностной волной. Угловая компонента представлена на рис. 2 при к=10. К сожалению, амплитуда этой компоненты уменьшается при удалении от источника. Разность радиальных компонент, порожденных глажкой и шероховатой границами раздела (к=10), представлена на рис. 3. На высоких частотах шероховатость может быть выделена вычитанием радиальных компонент, полученных при отражениях от гладкой и шероховатой. Эта разность есть квазипериодический сигнал с незначительным амплитудным затуханием.
0 0
Рис. 2. Угловая компонента (еф) вектора перемещений на верхней границе, порожденная шероховатой границей раздела (к=10)
-9
Рис. 3. Зависимость разности между радиальными компонентами, порожденными отражениями от плоской и шероховатой границ (к=10)
Угловая компонента (е^) вектора перемещений на свободной поверхности, порожденного отражением от шероховатой поверхности, представлена на рис. 4 при к=0.1. Это означает, что, если волна достаточно длинная (квазистатика), затухание угловой компоненты становится незначительным и ее амплитуда становится сравнимой с амплитудой поверхностных волн.
О о
Рис. 4. Угловая компонента (еф) вектора перемещений на верхней поверхности, порожденная отражением от шероховатой поверхности (к=0.1)
ВЫВОДЫ
1. Шероховатая граница раздела порождает угловую компоненту вектора перемещений на свободной границе при осесимметричном источнике. В случае, если граница раздела плоская, эта компонента отсутствует.
2. Амплитуда этой компоненты растет с ростом эффективной длины волны. Эта компонента затухает существенно слабее, чем остальные компоненты. Это дает возможность выделить ее на большом расстоянии от источника.
3. Если длина волны достаточно мала (по сравнению с расстоянием между слоями), то тогда обнаружить шероховатость можно с помощью вычитания радиальных компонент, порожденных плоской и шероховатой границами раздела.
4. Если длина волны достаточно велика (квазистатика), информация о шероховатых свойствах границы может быть получена с помощью угловой компоненты вектора перемещений на свободной поверхности.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. B.L.N. Kennett // Geophys. J R. astr. Soc. - 1972. - Vol. 28. - P. 249-266.
2. James R. Pettit, Anthony E. Walker, Michael J.S. Lowe. Improved Detection of Rough Defects for Ultrasonic Nondestructive Evaluation Inspections Based on Finite Element Modelling of Elastic Wave Scattering. IEEE Transactions on ultrasonics, ferroelectrics, and frequency control. -2015. - Vol. 62. - N 10. - P. 1797-1808.
3. Сибиряков Б.П., Сибиряков Е.Б. Области локального понижения давлений как вероятные аккумуляторы флюидов в геологических структурах // Геология и геофизика. - 2015. -Т. 56. - № 7. - С. 1391-1397.
4. Новацкий В. Теория упругости. - М.: Издательство «Мир», 1975. - 872 с.
5. Ziatdinov S.R.; Chestnut B.M. Extrinsic components of the Rayleigh wave // Geophysical Questions. - 2005. - Vol. 38. - P. 46-55.
© Е. Б. Сибиряков, 2016