Научная статья на тему 'Моделирование динамики частично насыщенных пороупругих тел на основе метода гранично-временных элементов'

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

CC BY
231
53
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
ПОРИСТАЯ СРЕДА / МОДЕЛЬ БИО / ПРЕОБРАЗОВАНИЕ ЛАПЛАСА / ГРАНИЧНЫЕ ИНТЕГРАЛЬНЫЕ УРАВНЕНИЯ / МЕТОД ГРАНИЧНЫХ ЭЛЕМЕНТОВ / ШАГОВАЯ СХЕМА / ПРИЗМАТИЧЕСКОЕ ТЕЛО / АНАЛИТИЧЕСКОЕ РЕШЕНИЕ / OROUS MEDIUM / MODEL BIOT'S / LAPLACE TRANSFORM / BOUNDARY INTEGRAL EQUATION / BOUNDARY ELEMENT METHOD / STEP SCHEME / PRISMATIC BODY / ANALYTIC SOLUTION

Аннотация научной статьи по физике, автор научной работы — Игумнов Л.А., Петров А.Н.

Дано описание математической модели частично насыщенной пористой среды, основанной на модели Био, с пятью базовыми функциями описания волнового процесса. В изображениях по Лапласу приведена математическая модель краевой задачи трехмерной динамической теории пороупругости. На основе теоремы операционного исчисления об интегрировании оригинала представлен шаговый метод численного обращения преобразования Лапласа. В качестве метода решения краевых задач трехмерной динамической теории пороупругости выбран прямой вариант метода граничных интегральных уравнений и приведено соответствующее граничное интегральное уравнение. Выписаны соответствующие матрицы фундаментальных и сингулярных решений трехмерной динамической теории пороупругости. Представлено краткое описание гранично-элементной дискретизации. Методическое обеспечение опирается на использование регуляризованного граничного интегрального уравнения, записанного с учетом преобразований симметрии задачи. Граничная поверхность исследуемого тела разбивается обобщенными восьмиузловыми четырехугольными элементами. Применяется согласованная поэлементная аппроксимация. Коллокационные точки решения граничного интегрального уравнения совпадают с узлами интерполяции неизвестных граничных функций. Для повышения точности интегрирования по элементу, не содержащему коллокационную точку, кроме формул интегрирования Гаусса, применяется также иерархический алгоритм интегрирования. Возникающие дискретные аналоги решаются методом Гаусса на основе шагового процесса получения значений граничных функций. Шаговый процесс определяется шаговым алгоритмом численного обращения преобразования Лапласа. Рассмотрена задача о скачке единичной поверхностной силы на торце призматического частично насыщенного пороупругого тела. В качестве пористого материала выбран песчаник. Для верификации гранично-элементной модели используется аналитическое решение соответствующей одномерной задачи. Проведено исследование решения задачи на сходимость по расчетной сетке, а также рассмотрено влияние параметра шаговой схемы на решение.

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Dynamics of partially saturated poroelastic solids by boundary-element method

The paper describes the mathematical model for partially saturated porous media based on the Biot’s model with five basic functions to characterize wave process. The mathematical model of the boundary value problem for the three-dimension dynamic theory of poroelasticity is given in Laplace transform. On the basis of the operational calculus theorem about the original integration, the step method of numerical inversion of the Laplace transform is presented. The direct method of boundary integral equations is selected to solve value problems of the three-dimensional dynamic poroelasticity theory, the corresponding boundary integral equation is given. The corresponding matrices of fundamental and singular solutions of the three-dimensional dynamic poroelastic theory are given. A brief description of the boundary-element discretization is presented. Methodological assurance is based on the regularized boundary integral equation usage. The regularized boundary integral equations are written considering the problem of symmetry. The boundary surface of the investigated solid is divided by generalized eight-node quadrangular elements. The consistent elementwise approximation is used. Collocation solution points of the boundary integral equation coincide with the interpolation nodes of unknown boundary functions. To increase the integration accuracy for an element not containing the collocation point, the hierarchical integration algorithm and Gauss integration formulas are applied. Arising discrete analogs are solved by the Gauss-based step process to obtain the values of the boundary functions. Step process is determined by a step algorithm of the numerical Laplace transform inversion. The problem of a unit surface force shock on the free end of the prismatic partially saturated poroelastic solid is considered. Sandstone is selected as the porous material. An analytical solution of the corresponding one-dimensional problem is used for the boundary element model verification. The dependence of computational grid on convergence of the problem solution is investigated; and the impact of step scheme parameter on the solution is examined.

Текст научной работы на тему «Моделирование динамики частично насыщенных пороупругих тел на основе метода гранично-временных элементов»

Игумнов Л.А., Петров А.Н. Моделирование динамики частично насыщенных пороупругих тел на основе метода гранично-временных элементов // Вестник Пермского национального исследовательского политехнического университета. Механика.- 2016. - № 3. - С. 47-61. DOI: 10.15593/perm.mech/2016.3.03

Igumnov L.A., Petrov A.N. Dynamics of partially saturated poroelastic solids by boundary-element method. PNRPU Mechanics Bulletin. 2016. No. 3. Рр. 47-61. DOI: 10.15593/perm.mech/2016.3.03

ВЕСТНИК ПНИПУ. МЕХАНИКА

№ 3,2016 PNRPU MECHANICS BULLETIN

http ://vestnik.pstu. ru/mechanics/about/inf/

001 10.15593/регш.шесЬ/2016.3.03 УДК 539.3

МОДЕЛИРОВАНИЕ ДИНАМИКИ ЧАСТИЧНО НАСЫЩЕННЫХ ПОРОУПРУГИХ ТЕЛ НА ОСНОВЕ МЕТОДА ГРАНИЧНО-ВРЕМЕННЫХ ЭЛЕМЕНТОВ

1 2 Л.А. Игумнов', А.Н. Петров2

Научно-исследовательским институт механики Национального исследовательского Нижегородского государственного университета им. Н.И. Лобачевского, Нижний Новгород, Россия 2Донской государственный технический университет, Ростов-на-Дону, Россия

О СТАТЬЕ

АННОТАЦИЯ

Получена: 3 июня 2016 г. Принята: 10 сентября 2016 г. Опубликована: 30 сентября 2016 г.

Ключевые слова:

пористая среда, модель Био, преобразование Лапласа, граничные интегральные уравнения, метод граничных элементов, шаговая схема, призматическое тело, аналитическое решение

Дано описание математической модели частично насыщенной пористой среды, основанной на модели Био, с пятью базовыми функциями описания волнового процесса. В изображениях по Лапласу приведена математическая модель краевой задачи трехмерной динамической теории пороупругости. На основе теоремы операционного исчисления об интегрировании оригинала представлен шаговый метод численного обращения преобразования Лапласа. В качестве метода решения краевых задач трехмерной динамической теории пороупругости выбран прямой вариант метода граничных интегральных уравнений и приведено соответствующее граничное интегральное уравнение. Выписаны соответствующие матрицы фундаментальных и сингулярных решений трехмерной динамической теории пороупру-гости. Представлено краткое описание гранично-элементной дискретизации. Методическое обеспечение опирается на использование регуляризованного граничного интегрального уравнения, записанного с учетом преобразований симметрии задачи. Граничная поверхность исследуемого тела разбивается обобщенными восьми-узловыми четырехугольными элементами. Применяется согласованная поэлементная аппроксимация. Коллокационные точки решения граничного интегрального уравнения совпадают с узлами интерполяции неизвестных граничных функций. Для повышения точности интегрирования по элементу, не содержащему коллока-ционную точку, кроме формул интегрирования Гаусса, применяется также иерархический алгоритм интегрирования. Возникающие дискретные аналоги решаются методом Гаусса на основе шагового процесса получения значений граничных функций. Шаговый процесс определяется шаговым алгоритмом численного обращения преобразования Лапласа. Рассмотрена задача о скачке единичной поверхностной силы на торце призматического частично насыщенного пороупругого тела. В качестве пористого материала выбран песчаник. Для верификации гранично-элементной модели используется аналитическое решение соответствующей одномерной задачи. Проведено исследование решения задачи на сходимость по расчетной сетке, а также рассмотрено влияние параметра шаговой схемы на решение.

© ПНИПУ

© Игумнов Леонид Александрович - доктор физико-математических наук, профессор, e-mail: igumnov@mech.unn.ru Петров Андрей Николаевич - кандидат физико-математических наук, старший научный сотрудник, e-mail: andrey.petrov@mech.unn.ru

Leonid A. Igumnov - Doctor of Physical and Mathematical Sciences, Professor, e-mail: igumnov@mech.unn.ru Andrey N. Petrov - PhD in Physical and Mathematical Sciences, Senior Researcher, e-mail: andrey.petrov@mech.unn.ru

DYNAMICS OF PARTIALLY SATURATED POROELASTIC SOLIDS BY BOUNDARY-ELEMENT METHOD

L.A. Igumnov1, A.N. Petrov2

1

Research Institute of Mechanics, National Research Lobachevsky State University of Nizhny Novgorod, Nizhny Novgorod, Russian Federation

Don State Technical University, Rostov-on-Don, Russian Federation

ABSTRACT

The paper describes the mathematical model for partially saturated porous media based on the Biot's model with five basic functions to characterize wave process. The mathematical model of the boundary value problem for the three-dimension dynamic theory of poroelasticity is given in Laplace transform. On the basis of the operational calculus theorem about the original integration, the step method of numerical inversion of the Laplace transform is presented. The direct method of boundary integral equations is selected to solve value problems of the three-dimensional dynamic poroelasticity theory, the corresponding boundary integral equation is given. The corresponding matrices of fundamental and singular solutions of the three-dimensional dynamic poroelastic theory are given. A brief description of the boundary-element discretization is presented. Methodological assurance is based on the regularized boundary integral equation usage. The regularized boundary integral equations are written considering the problem of symmetry. The boundary surface of the investigated solid is divided by generalized eight-node quadrangular elements. The consistent elementwise approximation is used. Collocation solution points of the boundary integral equation coincide with the interpolation nodes of unknown boundary functions. To increase the integration accuracy for an element not containing the collocation point, the hierarchical integration algorithm and Gauss integration formulas are applied. Arising discrete analogs are solved by the Gauss-based step process to obtain the values of the boundary functions. Step process is determined by a step algorithm of the numerical Laplace transform inversion. The problem of a unit surface force shock on the free end of the prismatic partially saturated poroelastic solid is considered. Sandstone is selected as the porous material. An analytical solution of the corresponding one-dimensional problem is used for the boundary element model verification. The dependence of computational grid on convergence of the problem solution is investigated; and the impact of step scheme parameter on the solution is examined.

© PNRPU

Введение

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

Началом исследований волновых процессов в насыщенных пористых средах послужила работа Я.И. Френкеля [1]. Л.Я. Косачевский [2] показал, что теория М. Био [3,4], как и подход Я.И. Френкеля, опирается на те же соотношения между напряжениями и деформациями, но отличается большей общностью. Естественным продолжением этих работ стали попытки различных авторов обобщить теорию на случай пористых сред, насыщенных двумя несме-шивающимися текучими наполнителями. Принципиальным отличием в данном случае является необходимость учитывать поверхностные натяжения и капиллярные эффекты, возникающие в результате взаимодействия наполнителей. ВгШБаег! [5] первым расширил теорию Био на случай трехфазной среды, однако исключив из рассмотрения капиллярные эффекты.

ARTICLE INFO

Received: 3 June 2016 Accepted: 10 September 2016 Published: 30 September 2016

Keywords:

porous medium, model Biot's, Laplace transform, boundary integral equation, boundary element method, step scheme, prismatic body, analytic solution

Общая система уравнений многофазной пористой среды была предложена Николаевским [б]. Работы Garg and Nayfeh [7], Santos и др. [8, 9], Tuncay и Corapcioglu [lO, ll], Wei и Muraleetharan [l2], Lo и др. [l3], Lu и др. [l4] по теории пористых сред с двумя жидкостями принадлежат концептуально к общему подходу, отличаясь в деталях реализации. Метод в целом состоял в разработке соотношений между напряжениями и деформациями для составляющих фаз с последующим использованием их в уравнении баланса в рамках стандартной механики сплошных сред. Berryman [l5] вывел уравнения для частично насыщенных сред с использованием предположения о незначительном изменении капиллярного давления при прохождении через них акустического сигнала. Smeulders и др. [l6] изучали распространение волн сжатия в пористой среде, заполненной жидкостью с небольшой объемной долей газа. В работах Thomas [l7, l8] предложена конечно-элементная модель для анализа волновых полей напряжений, перемещений и порового давления в ненасыщенном морском дне. Мейснер и Беккер [l9] установили влияние капиллярных эффектов на динамическое поведение частично насыщенного песка. Обзор доступных аналитических решений задач о распространении волн в насыщенных пористых средах можно найти в статье Schanz [20]. Однако для ненасыщенных моделей известно немного существующих аналитических решений. Можно отметить решение, полученное в работе Li и Schanz [2l], для одномерного частично насыщенного пороупругого стержня. Поскольку получение аналитического решения часто возможно лишь в некоторых специальных случаях и при определенных видах граничных условий, решение задач динамики частично насыщенных пористых сред требует привлечения развитых численных методов, таких как МКЭ и МГЭ. В случаях бесконечных и полубесконечных тел и сред метод МГЭ является наиболее предпочтительным универсальным подходом к исследованию волн. Для рассмотрения задачи о распространении волны необходимы динамические формулировки исходной системы дифференциальных уравнений и соответствующие ей ГИУ, фундаментальные и сингулярные решения. В работах Gatmiri и Jabbari [22, 23] получены фундаментальные решения для статических и квазистатических задач теории ненасыщенных грунтов в двух- и трехмерной постановке. Соответствующие термо-упругие фундаментальные решения двух и трехмерных статических и квазистатических задач для частично насыщенных сред представлены в работах Jabbari и Gatmiri [24], Gatmiri и др. [25], Maghoul и др. [26]. Динамические формулировки ГИУ и фундаментальные решения для двухмерных задач представлены Maghoul и др. в [27], для трехмерных Ashayeri и др. в [28], Li в [29]. Авторами рассматривается трехкомпонентная пороупругая среда Био с пятью базовыми функциями описания волнового процесса - перемещениями скелета и поровым давлением жидкости и газа. ГИУ получены на основе метода взвешенных невязок, фундаментальные решения представлены в терминах преобразования Лапласа по переменной времени, при построении шаговой гранично-элементной схемы используется метод квадратур свёрток.

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

Математическая модель

Рассматривается частично насыщенный пористый материал, который представляет собой упругий скелет с порами, заполненными жидкостью и газом, при этом исключается испарение жидкости и растворение газа [29]. Кроме того, предполагается, что все три фазы являются сжимаемыми, имеют постоянную температуру, а любые изменения температуры не учитываются.

Для частично насыщенного пористого материала пористость обозначается как

V

void

9 = "Vf, (1)

где Vvoid - объем взаимосвязанных пор в образце объема V . Закрытые поры рассматриваются как часть упругого скелета. Степень насыщения материала жидкостью Sw и газом Sa определяется как отношение занимаемого наполнителем объема Vw и Va соответственно к общему объему пор, т. е.

V V

Sw , Sa , Sa + Sw = 1. (2)

w jt ' a jt ' a w \ /

void void

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

w a

ния на границе раздела сред, возникает разность давлений в жидкости p и газе p , которая называется капиллярным давлением pc. Капиллярное давление может быть представлено как функция степени насыщения

pc = pa _ pw = pdS?. (3)

Величина pd определяет давление газа, необходимое для вытеснения жидкости из пор. Коэффициент 0, характеризующий распределение размера зерен скелета, обычно принимает значения от 0,2 до 3. Se - эффективная водонасыщенность определяется следующим образом [29]:

Se =

0, Sw < S

Sw _ Sw

S _ S '

ra rw

1, Sw > S

Srw < S < Sra , (4)

где - остаточная водонасыщенность; £га - остаточная газонасыщенность.

Определяющие соотношения для суммарного напряжения вслед за [30] могут быть записаны в виде

ау = а' -8уа(+ 8ара). (5)

Здесь а' обозначено эффективное напряжение и а = 1 - К / К - коэффициент эффективных напряжений, определяющий зависимость сжимаемости упругого скелета от его объ-

rw

емного модуля К и объемного модуля его зерен К,. Скелет считается изотропным и однородным с двумя константами упругого материала: объемный модуль К и модуль сдвига О . Исходя из этих предположений получена формулировка для напряжения:

а„ =

К - - ок^ + О (и,,; + щ,,)

(6)

Кроме того, в допущении о малых градиентах деформации полагаем, что компоненты тензора деформаций твердого тела г ^ и смещения и, связаны соотношением

8„ =2 {ии + ии)

Уравнения баланса массы для твердого тела и наполнителей имеют вид

д[(1 -ф)р, ]

дг

^ (ф^Р г). дг

д №Ра)

дг

+

(1 -Ф) Р,

ди.

Ф^Р,

Ф0а Ра

дг

д (щ + <)

дг

д (и,- + <) дг

= 0,

= РГ1 ,

~ та

= Ра ,

(7)

(8) (9)

(10)

„ и/

где р,, рг, ра - плотности твердотельного скелета, жидкости и газа соответственно, ui и и - компоненты смещений наполнителей, относительных по отношению к скелету. Формулировка предполагает, что растворения газообразного наполнителя в жидкости не происходит, а градиенты плотностей, насыщенности и пористости пренебрежимо малы. Для трехфазного материала можно записать

Ои,,, +[К + О) и]л -«(С + « )+ Ъ =Р^ + Ф^Рг ^ + Ф^аРа ^

дг2

дг2

дг2

Ф0г

диг

дг

(

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

= -К,.

о ди, Ф^^-^ = -Ка

Р. +Рг

Р" +Ра

д^и, ~дё

2л, г

• + Рг

д2и

д 2и,.

- + Ра

дг2

д 2и,а'

(11) (12) (13)

дг дг2 а дг2 )

Уравнение (11) записано для скелета, уравнения (12), (13) - для жидкости и газа соответственно. Введена объемная сила ¥{ и усредненная плотность

р = (1 - ф)р, + фОгрг + фОара. Перенос наполнителей определяется законом Дарси. Соответствующие фазовые проницаемости определяются как кг = Кгк / цг и ка = Кгак / ца, где Кг и Кга - относительные фазовые проницаемости; к - внутренняя проницаемость материала; цг и ца - вязкости наполнителей. Относительная фазовая проницаемость определяется как отношение фазовой проницаемости к абсолютной. Используем вслед за [31] следующие экспериментально установленные зависимости для относительных про-ницаемостей в случае газоводяной смеси в порах:

^ишпоу Ь.Л., РегтоуЛ.Ы. /РШИРПМвекап^ВиНеПп 3 (2016) 47-61

К* _ 8в2+30)/Э , Кт = (1 - 8е )2 [1 - ^)/0] .

Система уравнений (8)-(13) полностью описывает поведение частично насыщенного пороупругого континуума. Выберем в качестве базовых функций давление в порах р*,

а

р и перемещение скелета и1.

Для дальнейшего рассмотрения применяется интегральное преобразование Лапласа. Пусть /(г) - функция по г, /(г)= 0 для г < 0. Функция преобразования по Лапласу и ее обращение определено следующим образом:

/(5) = | f (г)е , 5 = а + 1ш,

0

1 а+г'ад

/(г) =- [ /((формула Римана-Меллина),

2л/ :

(15)

(16)

где 5 - комплексный параметр преобразования; а - вещественное число большее, чем вещественные части всех особенностей /(5"). Применяя преобразование Лапласа к уравнениям (12), (13) и исключая относительные смещения и* и и**, получаем

Р / I ~ ~ \ ла У

дг _ щ =

(2иг), %_5иа

( + Ра5\ )

ФР*5 ФРа5

где ц*, Ц - потоки наполнителей. В (17) введены следующие обозначения:

к фр 5 к фр 5

Р__^г*_ у__атГа

(17)

Ф8* + К*Р*5 Ф8а + КаРа5 Подстановка (17) в уравнения (8)—(11) позволяет получить формулировку

(18)

^ +

'К + |)- (Р - - у8аРа) 5и - (а - рад - (а - у)Бар' _-Рг, (19)

(

Ф

Л

С8 + я* - 8'и

V У

5р * + ^ р* - (С8аа8* + ^ ) 5ра _- Г

-(а-у)За5й1л -+ 8)5р* -

Р*5

Ф

С8аа8а +~ 8а - 8и

К „

5ра +У8ара11 _- Г, Ра5 '

(20) (21)

где

а-Ф

С _—[Г~ , _ -0(8* - ) , 8аа _ 8а + 0(8* - 8 т*) , 8и _ -

К5 Г

Уравнения (19)-(21) могут быть записаны в операторной форме:

0(8„ - 8т*) Г - 8т е+1)/0

*_т*

8 - V

V та гл> У

В _ р * _ I *

[ ра _ - 1а

X)

с несамосопряженным оператором

В

А18ч + Ад г 5А3д, А4д А5д] А6 А7

А9 А10

А: = СУ2-(р-р^-УБара), А2 = К + Аз =-(а-Р)^, А4 = -(а-У)^,

А5 = -(а-Р)^, Аб = -ГСЯ^т + Кг 8* -]* + —У2, А7 = -(А + ^).,

I К- ) Р-5

А =-(а-у)^, А = -(Д, +8и), Аю = -

(

Ф

Л

V Ка )

5 + ^ У 2.

Ра5

Применение интегрального преобразования позволяет строить метод граничных интегральных уравнений.

2. Граничные интегральные уравнения и фундаментальные решения

и™ 1} и Г и Т

и = Р ш р р а.

р ™ р га р аа

Фундаментальное решение находится из следующего операторного соотношения:

Ви + 15(х - у) = 0, (23)

где I - единичная матрица; 8 - дельта-функция Дирака; х, у е Я3. Итоговая матрица фундаментальных решений

(24)

содержит девять компонент - компоненты фундаментальных решений перемещений 17™

и давлений Р, Раа, а также шесть спаренных компонент. Все компоненты матрицы приведены в [29].

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

Интегральное уравнение динамической пороупругости, построенное на основе метода взвешенных невязок, имеет вид

7

У Ч

Т р №№

и,

и

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Ц7 ш р тт

рш

р га р аа

Г 7 "

- д Г =

а г

д _

П7

У

Т Т

г-7 аз

Г и 1 " 0 0" Г и 1 "

07 —— 77 та р " ^г + 0 с 0 р "

07 тТ 7аа _ Р 0 _ 0 0 с _ Р 0 _

(25)

гр ££ г д ...

где Т д жа

Т а£ д aw даа

< = Р Р^ (( + рУ

- матрица сингулярных решений; р =а^П] - поверхностная сила;

П - поток жидкости; Ца =——((" +р £2и. )пг - поток газа; п -

Ра£

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

3. Гранично-элементная (ГЭ) дискретизация

Чтобы ввести ГЭ-дискретизацию, рассмотрим регуляризованное уравнение [32]:

«оик (Х £)+\{тгк (Х У £)Ц (К £)—ТЦ (Х У £)Ц (Х ¿)—Цк (x, у ^ £)) & = 0,

(Xе5), г=[р,г2,гъ,ца]Т, и = (,и2,щ,рК,ра),

где ил, Тк - матрицы фундаментальных и сингулярных решений; Т- матрица компонент, содержащих особенности; г - вектор обобщенных усилий; и - вектор обобщенных перемещений.

Базовый процесс ГЭ-дискретизации состоит в разбиении поверхности на ЫЕ граничных элементов Ее (1 < е < ЫЕ) совокупностью четырехугольных и треугольных восьмиуз-

ловых биквадратичных элементов. При этом треугольные элементы рассматриваются как вырожденные четырехугольные элементы, каждый из которых отображается на контрольный элемент Ае (каждый Ае - это либо квадрат £ = (£^ £ 2) е [—1,1]2, либо треугольник 0 <£ +£2 <1, £1 >0, £2 >0).

Неизвестные граничные поля (и, г) определяются через узловые значения и = и(2 ) и гк = г(гк) в интерполяционных узлах. При этом для расчетного значения параметра £ будем иметь следующие выражения граничных перемещений и поверхностных сил внутри элемента 5к:

и,(у) = УЯ1(£)иГ,г), г = 1,2,3; у е ,

1=1

гг(у) = г?к 1г = 1,2,3; у е ^.

Здесь Я1 (£) - функции формы для линейного четырехугольного элемента.

Для получения дискретного аналога ГИУ применим метод коллокации. В качестве узлов коллокации будем выбирать узлы аппроксимации исходных граничных функций. В итоге формируется система линейных алгебраических уравнений

1 N 4 N

1 из ш + УУ Лш,к,' и?( к,') = У Бш,кгк

2 г V ] V ]'

^ к=1 1=1 к=1

1Ы 4 N

1 ~ т + ^^ дт,к,1 (к,1) = ^ £>т

8

к=1 1=1 к=1

Уравнения записаны соответственно в узлах аппроксимации обобщенных граничных перемещений и обобщенных поверхностных сил.

дт^1=11 [я1 сад х, / (& р)-5х(к, , ук (^))>к в)—и,

-1 -1

Врк = 11 и (хт, Ук (& р.

-1 -1

4. Шаговая схема численного обращения преобразования Лапласа

Рассмотрим метод, опирающийся на теорему об интегрировании оригинала, - шаговый метод численного обращения преобразования Лапласа. Рассмотрим следующий интеграл:

t

у(0 = 1 / (Т)Л. (26)

о

Интеграл (26) порождает задачу Коши для обыкновенного дифференциального уравнения

—х(г) = ях^)+С, X (0) = о. —

Интеграл (26) заменяется квадратурной суммой, весовые множители которой определяются с помощью изображения по Лапласу / и линейного многошагового метода [33-36]. Дальнейшее изложение идет с учетом результатов этих работ. Традиционный шаговый метод интегрирования оригинала состоит в том, что интеграл (26) считается по следующему соотношению:

у (0) = 0, у (п ^) = £ Шк (^), п = 1,...М, (27)

к=1

где

Шп(А0 = -—1 /I е гпф—Ф.

2- 0 I Аt )

Единственная аппроксимация, которая используется при выводе формулы (27), заключается в применении линейного многошагового метода для решения задачи Коши дифференциального уравнения первого порядка. Все остальные вычисления проводятся непосредственно. Многошаговый метод должен быть порядка точности р > 1, являясь

строго ноль-устойчивым или А-устойчивым, а функция /(я) ограничена в правой полуплоскости относительно прямой (с-/да,с + /да), то есть

/(я) < КЫ ц при К <да, ц > 0.

Если функция f(s) аналитична и ограничена в области

|arg(s-c)| <л-ф,

где ф<|, критерий устойчивости может быть ослаблен до А(а)-устойчивости. Выражаясь в терминах характеристической функции у(z):

• функция у (z) не должна содержать нулей и полюсов внутри замкнутого единичного круга |z| < 1) за исключением однократного нуля при |z| = 1;

• |argу(z)| < л-а, где |z| < 1 при а >ф;

• At-1у(e-At) = 1 + 0(Atp ), где At ^ 0 при p > 1.

К соответствующим примерам многошаговых методов относятся методы дифференцирования назад порядка p < 6: для А-устойчивого метода дифференцирования назад

второго порядка (а = 90°) можем записать следующую характеристическую функцию:

у (z) = 3/2 - 2 z + z 72.

Следующая аппроксимация, которая используется в методе, - численное интегрирование, при получении значений шп (At) аппроксимируем с L равными шагами 2л/L :

R-n l-1 ®n (At) = — 2 f

1=0

у(Re ~ ) At

e"nl2L. (28)

При условии, что функция /(.) в уравнении (28) вычисляется с некоторой погрешностью в выбор Ь = N и Яп = -у/в допускает погрешность вычисления шп порядка 0(>/в) .

5. Численные эксперименты

Рассматривается задача о действии торцевой силы на призматическое пороупругое тело длиной 3 м, высотой и шириной 1 м, с жестко закрепленным концом. Краевая задача представлена на рис. 1. Отклики перемещения иг и поверхностные усилия , вызванные

силой Е (г) = 1 Н/м2, наблюдаются соответственно на нагруженном и на закрепленном торцах. Параметры частично насыщенного пористого материала соответствуют параметрам песчаника: К = 1,02-109 Н/м2, О = 1,44-109Н/м2, ф = 0,23, р£ = 2650 кг/м3,

р. =997 кг/м3, ра = 1,10 кг/м3, К = 3,55-1010 Н/м2, К. = 2,25-109 Н/м2,

Ка =1,10-105Н/м2, к =2,5-10—12 м2, л. = 1,0-10—3 Нс/м2, ла = 1,8-10—5Нс/м2, = 0,9,

= 0, 5та = 1, 0 = 1,5. На временном интервале от 0 с до 0,035 с шаг по времени Аг

выбирается из условия Ь = N, параметр Я принят равным 0,997. Задача решается в трёхмерной постановке. Заданные граничные условия позволяют сравнивать решение исходной задачи с решением соответствующей одномерной задачи, представленным в [21].

При решении использовались три гранично-элементные сетки с различной степенью пространственной дискретизации (рис. 2). Наличие двух плоскостей симметрии дает возможность рассматривать только часть сетки, при этом четверть сетки а содержит 14 элементов, сетки б - 56 элементов, сетки в - 126 элементов.

Рис. 1. Геометрическая модель и граничные условия задачи о частично насыщенной

пороупругой консоли Fig. 1. Geometry and boundary conditions of a partially saturated poroelastic column

а б в

Рис. 2. Гранично-элементные сетки: a - 56 элементов; б - 224 элемента; в - 504 элемента Fig. 2. Boundary-element meshes: a - 56 elements; b - 224 elements; c - 504 elements

Рис. 3. Гранично-элементные решения, полученные на трех различных сетках:

a - перемещения; б - усилия Fig. 3. Boundary-element solutions on the three different meshes for: a - displacements; b - tractions

Рис. 4. Гранично-элементные решения, полученные при трех различных значениях параметра шагового метода: a - перемещения; б - усилия

Fig. 4. Boundary-element solutions for the three different values of the time-step method parameter for: a - displacements; b - tractions

Результаты исследования сходимости по расчетной сетке показывают, что с ростом числа элементов гранично-элементные решения для перемещения (рис. 3, а) и напряжения (рис. 3, б) принимают значения, близкие к рассчитанным по аналитическим формулам. Варьирование параметром шаговой схемы At приводит к тому, что при снижении его до значения 0,00002185c в расчетах на сетке c наблюдается незначительное уточнение амплитуды решения (рис. 4, а, б), однако при дальнейшем его уменьшении решение расходится.

Заключение

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

Работа выполнена при финансовой поддержке Российского научного фонда (проект № 15-19-10056).

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Библиографический список

1. Frenkel J. On the theory of seismic and seismoelectric phenomena in a moist soil // Journal of Physics. - 1944. - Vol. 8. - P. 230-241.

2. Косачевский Л.Я. О распространении упругих волн в двухкомпонентных средах // Прикладная математика и механика. - 1959. - Т. 23, № 6. - С. 1115-1123.

3. Biot M.A. Theory of propagation of elastic waves in a fluid-saturatedporous solid. I. Low-frequency range // Journal of the Acoustical Society of America. - 1956. - Vol. 28. - P. 168-178. DOI: 10.1121/1.1908239

4. Biot M.A. Theory of propagation of elastic waves in a fluid-saturated porous solid. II. Higher-frequency range // Journal of the Acoustical Society of America. - 1956. - Vol. 28. - P. 179-191. DOI: 10.1121/1.1908241

5. Brutsaert W. The propagation of elastic waves in unconsolidated unsaturated granular mediums // Journal of Geophysical Research. - 1964. - Vol. 69. - P. 243-257. DOI: 10.1029/JZ069i002p00243

6. Nikolaevskij V.N. Mechanics of porous and fractured media. Singapore, World Scientific, 1990. - 472 p.

7. Garg S.K., Nayfeh A.H. Compressional wave propagation in liquid and/or gas saturated elastic porous media // Journal of Applied Physics. - 1986. - Vol. 60. - P. 3045-3055. DOI: 10.1063/1.337760

8. Santos J.E., Corbero J.M., Douglas J.Jr. Static and dynamic behavior of a porous solid saturated by a two-phase fluid // Journal of the Acoustical Society of America. - 1990. - Vol. 87. - P. 1428-1438. DOI: 10.1121/1.399439

9. A model for wave propagation in a porous medium saturated by a two-phase fluid / J.E. Santos, J.Jr. Douglas, J. Corbero, O.A. Lovera // Journal of the Acoustical Society of America. - 1990. - Vol. 87. -P. 1439-1448. DOI: 10.1121/1.399440

10. Tuncay K., Corapcioglu M.Y. Body waves in poroelastic media saturated by two immiscible fluids // Journal of Geophysical Research. - 1996. - Vol. 101. - P. 25149-25159. DOI: 10.1029/96JB02297

11. Tuncay K., Corapcioglu M.Y. Wave propagation in poroelastic media saturated by two fluids // Journal of Applied Mechanics. - 1997. - Vol. 64. - P. 313-320. DOI: 10.1115/1.2787309

12. Wei C., Muraleetharan K.K. A continuum theory of porous media saturated by multiple immiscible fluids: I. Linear poroelasticity // International Journal of Engineering Science. - 2002. - Vol. 40. -P. 1807-1833. DOI: 10.1016/S0020-7225(02)00068-X

13. Lo W.-C., Sposito G., Majer E. Wave propagation through elastic porous media containing two immiscible fluids // Water Resources Research. - 2005. - Vol. 41. - P. 140-160. DOI: 10.1029/2004WR003162

14. Lu J.-F., Hanyga A., Jeng D.-S. A mixture-theory-based dynamic model for a porous medium saturated by two immiscible fluids // Journal of Applied Geophysics. - 2007. - Vol. 62. - P. 89-106. DOI: 10.1016/j.jappgeo.2006.08.002

15. Berryman J.G. Waves in partially saturated porous media / Eds. In W.E. Fitzgibbon, M.F. Wheeler // Wave Propagation and Inversion (Proceedings of SIAM Conference on Mathematical and Computational Issues in Geophysical Fluid and Solid Mechanics) SIAM, Philadelphia, 1992. - P. 1-25.

16. Smeulders D., De La Rosette J., Van Dongen M. Waves in partially saturated porous media // Transp. Porous. Med. - 1992. - Vol. 9. - P. 25-37.

17. Thomas S.D. A finite element model for the analysis of wave induced stresses, displacements and pore pressures in an unsaturated seabed I // Theory. Comput. Geotech. - 1989. - Vol. 8. - P. 1-38.

18. Thomas S.D. A finite element model for the analysis of wave induced stresses, displacements and pore pressures in an unsaturated seabed II // Model verification. Comput. Geotech. - 1995. -Vol. 17. - P. 107-132.

19. Meissner H., Becker A. Dynamic behaviour of partially saturated sand // T. Built Environ. -1995. - Vol. 14. - P. 45-55.

20. Schanz M. Poroelastodynamics: Linear models, analytical solutions, and numerical methods // Appl. Mech. Re. - 2009. - Vol. 62 (3). - 15 p.

21. Li P., Schanz M. Wave propagation in a 1-d partially saturated poroelastic column // Geophys. J. Int. - 2011. - Vol. 184. - P. 1341-1353.

22. Gatmiri B., Jabbari E. Time-domain green's function for unsaturated soils. Part I: two-dimensional solution // International Journal of Solids and Structures. - 2005. - Vol. 42. - P. 5971-5990.

23. Gatmiri B., Jabbari E. Time-domain Green's function for unsaturated soils. Part II: three-dimensional solution // International Journal of Solids and Structures. - 2005. - Vol. 42. -P. 5991-6002.

24. Jabbari E., Gatmiri B. Thermo-Poro-Elastostatic Green Functions for Unsaturated Soils // Comput. Model. Eng. Sci. - 2008. - Vol. 18 (1). - P. 31-43.

25. Gatmiri B., Maghoul P., Duhamel, D. Two-dimensional transient thermo-hydro-mechanical fundamental solutions of multiphase porous media in frequency and time domains // Int. J. Sol. Struct. -2010. - Vol. 47. - P. 595-610. DOI: 10.1016/j.ijsolstr.2009.10.022

26. Maghoul P., Gatmiri B., Duhamel D. Three-dimensional transient thermo-hydro-mechanical of unsaturated soils // Int. J. Numer. Anal. Methods Geomech. - 2009. - Vol. 34. - P. 297-329. DOI: 10.1002/nag.820

27. Maghoul P., Gatmiri B., Duhamel D. Boundary integral formulation and two-dimensional fundamental solutions for dynamic behaviour analysis of unsaturated soils // International Journal of Soil Dynamics and Earthquake Engineering. - 2011. - Vol. 31 (11). - P. 1480-1495.

28. Analytical 3d transient elastodynamic fundamental solution of unsaturated soils / I. Ashayeri, M. Kamalian, M.K. Jafari, B. Gatmiri // Int. J. Numer. Anal. Methods Geomech. - 2010. - Vol. 35. -P.1801-1829.

29. Li P. Boundary element method for wave propagation in partially saturated poroelastic continua. Computation in engineering and science. Graz: Verlag der Technischen Universität, 2012.

30. Bishop A.W., Blight G.E. Some aspects of effective stress in saturated and partly saturated soils // Geotechnique. - 1963. - Vol. 13. - P. 177-197.

31. Brooks R.H., Corey A.T. Hydraulic properties of porous media. - In Hydraulic Papers. Colorado State Univ., 1964.

32. Баженов В.Г., Игумнов Л.А. Методы граничных интегральных уравнений и граничных элементов в решении задач трехмерной динамической теории упругости с сопряженными полями. -М.: Физматлит, 2008. - 352 с.

33. Lubich C. Convolution quadrature and discretized operational calculus // Numer. Math. I. -1988. - Vol. 52 (2). - P. 129-145.

34. Lubich C. Convolution quadrature and discretized operational calculus // Numer. Math. II. -1988. - Vol. 52 (4). - P. 413-425.

35. Schanz M. Wave Propagation in Viscoelastic and Poroelastic Continua. - Berlin Springer, 2001. - 170 p.

36. Литвинчук С.Ю., Петров А.Н., Сабаева Т.А., Решение задачи о действии давления внутри полости в пороупругом полупространстве // Проблемы прочности и пластичности. - 2014. -№ 76 (2). - C. 97-105.

References

1. Frenkel J. On the theory of seismic and seismoelectric phenomena in a moist soil. Journal of Physics, 1944, vol. 8, pp. 230-241.

2. Kosachevskii L.Y. O rasprostranenii uprugikh voln v dvukhkomponentnykh sredakh [On the propagation of elastic waves in two-component media]. PMM, 1959, vol. 23 (6), pp. 1115-1123.

3. Biot M.A. Theory of propagation of elastic waves in a fluid-saturatedporous solid. I. Low-frequency range. Journal of the Acoustical Society of America, 1956, vol. 28, pp. 168-178. DOI: 10.1121/1.1908239

4. Biot M.A. Theory of propagation of elastic waves in a fluid-saturated porous solid. II. Higher-frequency range. Journal of the Acoustical Society of America, 1956, vol. 28, pp. 179-191. DOI: 10.1121/1.1908241

5. Brutsaert W. The propagation of elastic waves in unconsolidated unsaturated granular mediums. Journal of Geophysical Research, 1964, vol. 69, pp. 243-257. DOI: 10.1029/JZ069i002p00243

6. Nikolaevskij V.N. Mechanics of porous and fractured media. Singapore, World Scientific, 1990, 472 p.

7. Garg S.K., Nayfeh A.H. Compressional wave propagation in liquid and/or gas saturated elastic porous media. Journal of Applied Physics, 1986, vol. 60, pp. 3045-3055. DOI: 10.1063/1.337760

8. Santos J.E., Corbero J.M., Douglas J.Jr. Static and dynamic behavior of a porous solid saturated by a two-phase fluid. Journal of the Acoustical Society of America, 1990, vol. 87, pp. 1428-1438. DOI: 10.1121/1.399439

9. Santos J.E., Douglas J.Jr., Corbero J., Lovera O.A. A model for wave propagation in a porous medium saturated by a two-phase fluid. Journal of the Acoustical Society of America, 1990, vol. 87, pp. 1439-1448. DOI: 10.1121/1.399440

10. Tuncay K., Corapcioglu M.Y. Body waves in poroelastic media saturated by two immiscible fluids. Journal of Geophysical Research, 1996, vol. 101, pp. 25149-25159. DOI: 10.1029/96JB02297

11. Tuncay K., Corapcioglu M.Y. Wave propagation in poroelastic media saturated by two fluids. Journal of Applied Mechanics, 1997, vol. 64, pp. 313-320. DOI: 10.1115/1.2787309

12. Wei C., Muraleetharan K.K. A continuum theory of porous media saturated by multiple immiscible fluids: I. Linear poroelasticity. International Journal of Engineering Science, 2002, vol. 40, pp. 1807-1833. DOI: 10.1016/S0020-7225(02)00068-X

13. Lo W.-C., Sposito G., Majer E. Wave propagation through elastic porous media containing two immiscible fluids. Water Resources Research, 2005, vol. 41, pp. 140-160. DOI: 10.1029/2004WR003162

14. Lu J.-F., Hanyga A., Jeng D.-S. A mixture-theory-based dynamic model for a porous medium saturated by two immiscible fluids. Journal of Applied Geophysics, 2007, vol. 62, pp. 89-106. DOI: 10.1016/j.jappgeo.2006.08.002

15. Berryman J.G. Waves in partially saturated porous media. In W.E. Fitzgibbon and M.F. Wheeler, editors, Wave Propagation and Inversion. ( Proceedings of SIAM Conference on Mathematical and Computational Issues in Geophysical Fluid and Solid Mechanics) SIAM, Philadelphia, 1992, pp. 1-25.

16. Smeulders D., De La Rosette J., Van Dongen M. Waves in partially saturated porous media. Transp. Porous. Med., 1992, vol. 9, pp. 25-37.

17. Thomas S.D. A finite element model for the analysis of wave induced stresses, displacements and pore pressures in an unsaturated seabed I. Theory. Comput. Geotech, 1989, vol. 8, pp. 1-38.

18. Thomas S.D. A finite element model for the analysis of wave induced stresses, displacements and pore pressures in an unsaturated seabed II. Model verification. Comput. Geotech, 1995, vol. 17, pp. 107-132.

19. Meissner H., Becker A. Dynamic behaviour of partially saturated sand. T. Built Environ., 1995, vol. 14, pp. 45-55.

20. Schanz M. Poroelastodynamics: Linear models, analytical solutions, and numerical methods. Appl. Mech. Rev, 2009, vol. 62 (3), 15 p.

21. Li P., Schanz M. Wave propagation in a 1-d partially saturated poroelastic column. Geophys. J. Int., 2011, vol. 184, pp. 1341-1353.

22. Gatmiri B., Jabbari E. Time-domain green's function for unsaturated soils. Part I: two-dimensional solution. International Journal of Solids and Structures, 2005, vol. 42, pp. 5971-5990.

23. Gatmiri B., Jabbari E. Time-domain Green's function for unsaturated soils. Part II: three-dimensional solution. International Journal of Solids and Structures, 2005, vol. 42, pp. 5991-6002.

24. Jabbari E., Gatmiri B. Thermo-Poro-Elastostatic Green Functions for Unsaturated Soils. Comput. Model. Eng. Sci., 2008, vol. 18, no. 1, pp. 31-43.

25. Gatmiri B., Maghoul P., Duhamel D. Two-dimensional transient thermo-hydro-mechanical fundamental solutions of multiphase porous media in frequency and time domains. Int. J. Sol. Struct., 2010, vol. 47, pp. 595-610. D0I:10.1016/j.ijsolstr.2009.10.022

26. Maghoul P., Gatmiri B., Duhamel D. Three-dimensional transient thermo-hydro-mechanical of unsaturated soils. Int. J. Numer. Anal. Methods Geomech., 2009, vol. 34, pp. 297-329. DOI: 10.1002/nag.820

27. Maghoul P., Gatmiri B., Duhamel D. Boundary integral formulation and two-dimensional fundamental solutions for dynamic behaviour analysis of unsaturated soils. International Journal of Soil Dynamics and Earthquake Engineering, 2011, vol. 31 (11), pp. 1480-1495.

28. Ashayeri I., Kamalian M., Jafari M.K, Gatmiri, B. Analytical 3d transient elastodynamic fundamental solution of unsaturated soils. Int. J. Numer. Anal. Methods Geomech., 2010, vol. 35, pp. 1801-1829.

29. Li P. Boundary element method for wave propagation in partially saturated poroelastic continua. Computation in engineering and science, Verlag der Technischen Universität Graz, 2012.

30. Bishop A.W., Blight G.E. Some aspects of effective stress in saturated and partly saturated soils. Geotechnique, 1963, vol. 13, pp. 177-197.

31. Brooks R.H., Corey A.T. Hydraulic properties of porous media. In Hydraulic Papers. Colorado State Univ., 1964.

32. Bazhenov V.G., Igumnov L.A. Boundary Integral Equations & Boundary Element Methods in treating the problems of 3D elastodynamics with coupled fields [Metody granichnyh integralnyh uravneniy i granichnyh elementov v reshenii zadach tryohmernoy dinamicheskoy teorii uprugosti s sopryazhyonnymi polyami], PhysMathLit, Moscow, 2008, 352 p.

33. Lubich C. Convolution quadrature and discretized operational calculus. Numer. Math. I. 1988, vol. 52 (2), pp. 129-145.

34. Lubich C. Convolution quadrature and discretized operational calculus. Numer. Math. II. 1988, vol. 52 (4), pp.413-425.

35. Schanz M. Wave Propagation in Viscoelastic and Poroelastic Continua. Berlin Springer, 2001, 170 p.

36. Litvinchuk S.Yu., Petrov A.N., Sabaeva T.A. Reshenie zadachi o deistvii davleniia vnutri polosti v porouprugom poluprostranstve [Analyzing the problem of pressure acting inside a cavity in a poroelastic half-space]. Problems of strength and plasticity, 2014, vol. 76 (2), pp. 97-105.

i Надоели баннеры? Вы всегда можете отключить рекламу.