Вычислительные технологии Том 4, № 3, 1999
ОБ ЭФФЕКТЕ НАПРАВЛЕННОГО ПЕРЕНОСА МАСС ПРИ СЛОЖНОМ НАГРУЖЕНИИ *
О. П. БУШМАНОВА, А. Ф. РЕВУЖЕНКО Институт горного дела СО РАН , Новосибирск, Россия
Numerical modeling of the complex Earth loading process is considered acted upon by tide forces. The case is studied when the viscous properties of the medium are changing with the change in the depth. The properties of the internal nucleus, which is considered rigid, are considered. It has been shown numerically that in this case we also observe the effect of the differential rotation of the internal mass and the nucleus. The specific calculation results are presented and the comparison with the approximate calculation formula of the rotational velocity of the internal nucleus is performed.
В настоящее время имеется достаточное количество фактов, которые показывают, что в некоторых случаях динамические явления в шахтах связаны с деформированием Земли под действием приливных сил. Сами по себе приливные деформации ничтожны, но они действуют длительное время (время существования Земли) и могут приводить к необратимым изменениям. Например, если напряженно-деформированное состояние горного массива близко к критическому, то приливные деформации могут послужить пусковым механизмом для землетрясений и динамических проявлений горного давления. Благодаря эффекту накапливания этот результат может быть значительным. Гипотеза о таком накапливании высказывалась в ряде работ [1-3]. В [4] были предложены один из возможных механизмов необратимого накапливания таких изменений и способ моделирования аналогичных процессов в лабораторных условиях.
Реальный процесс приливного деформирования чрезвычайно сложен. Здесь действует множество факторов, большинство из которых можно оценить только с определенной степенью достоверности. Поэтому большое значение в данном случае приобретает роль математических моделей. Первая проблема, которая возникает на этом пути, связана с выбором модели среды. Прежде всего ясно, что модель идеально упругого тела в данном случае является слишком грубой и никаких необратимых эффектов не даст (упругий материал “не помнит” историю своего нагружения). Поэтому для исследуемого эффекта можно выбрать модель линейной вязкой жидкости [5]. При достаточно большой вязкости инерционными силами можно пренебречь. При отсутствии массовых сил уравнения движения и уравнения состояния вязкой несжимаемой жидкости в плоском случае имеют вид
d(Jxx + д(гХу _ о д(гХу + dtjyy _ Q ди + ди _ Q
дх ду дх ду дх ду
jxx _ -p +2^, jyy _ -p +2^, jxy _ „( g + g). (1)
*Работа выполнена при финансовой поддержке фондов Интас — РФФИ, проект №95-0742. © О. П. Бушманова, А. Ф. Ревуженко, 1999.
Здесь, как обычно,ахх, ауу,аху — компоненты тензора напряжений, п,ь — компоненты вектора скорости, р — давление, ц — коэффициент вязкости.
В качестве первого шага по усложнению модели в данной работе предлагается учесть неоднородность распределения вязких свойств (ц = ц(х,у)). Это связанно с тем, что реальные процессы происходят в условиях ярко выраженной слоистости материала, которая присутствует как на глобальном уровне, так и на уровнях меньших масштабов. После преобразования к безразмерным величинам запишем уравнения (1) в переменных Ф — П (функция тока — вихрь):
АП + А\Пх + А^Пу + 2АзФху + А^Фуу — Фхх) = 0,
АФ + П = 0. (2)
Здесь используются общепринятые обозначения производных по х и у, А — оператор Лапласа, функция тока Ф определяется из уравнений
Фу = и, Фх = —V,
а коэффициенты имеют вид
А 2цх А 2цу А 2цху А цхх — цуу
ц ’ ц ’ ц ’ ц
Единицами масштаба (эталонными величинами) длин, напряжений и давления, скоростей, вязкости являются соответствующие данной задаче характерные величины.
Будем считать, что процесс деформирования осуществляется в эллиптической области, на границе которой вектор скорости направлен по касательной и его величина выбрана в качестве характерной скорости. В качестве характерного линейного размера выбираем полуось эллипса, в качестве характерной вязкости — вязкость на границе области.
Не ограничивая общности, можно считать, что на границе заданы следующие условия:
Ф„ = —1, Фт = 0. (3)
Здесь Фп, Фт — производные вдоль нормали и касательной к границе. Реальная деформация Земли, как уже отмечалось, ничтожна. Однако ясно, что остаточное перемещение за один цикл при очень малой деформации будет пропорционально деформации. Таким образом, если в расчетах значительно увеличить приливную деформацию, то увеличится и остаточное смещение, то есть при численном моделировании остаточное смещение накопится за меньшее число циклов. Поэтому для того, чтобы проверить качественное влияние переменной вязкости, в дальнейшем выбираем большой эксцентриситет эллипса.
Для получения конечномерного аналога системы (2) используется метод конечных элементов [6], а исследуемая эллиптическая область разбивается на треугольные конечные элементы с линейными функциями формы (рис. 1). Картина смещения материальных частиц строится при помощи численного решения дифференциальных уравнений:
йх . . йу . .
— = и(х,у), — = v(x,y),
йг у ’ йг у ’
где г — время.
В исследуемом случае задание функции тока и вихря скорости возможно лишь в табличном виде (в точках сетки), поэтому распределение скоростей определяется в точках
Рис. 1.
Рис. 2.
сетки, одно из семейств линий которой совпадает с линиями тока. Построение сетки конечных элементов, обладающей указанным свойством, осуществляется на основе решения задачи на произвольно выбранной сетке, а для аппроксимации линий сетки используются сплайн-функции третьего порядка.
Расчеты проводились по предложенной схеме для различных законов изменения вязкости и для различных эксцентриситетов эллипса. Как и следовало ожидать, увеличение вязкости материала при приближении к центру приводит к уменьшению скоростей частиц, то есть к большему запаздыванию. На рис. 2 показано положение частиц, первоначально находившихся на большой оси эллипса, при различных циклах нагружения в случае изменения вязкости по закону
ц(г) = ехр(а(г — 1)),
/ У
здесь г = \/ х2 + —, Ь = 0, 54 — малая полуось эллипса, а = —1.
V Ь2
Хорошо заметно, что при большом числе циклов срединный участок почти прямолинеен, то есть вращение внутренней области становится все более близким к вращению абсолютно твердого тела.
В работе [7] была предложена схема моделирования области с жестким ядром и дана приближенная формула для расчета скорости его вращения. В данной работе решается краевая задача для эллиптической области с жестким ядром радиуса Я в центре. На внешней границе эллипса также задан единичный вектор скорости, направленный по касательной к границе (условия (3)). На границе с ядром вектор скорости направлен по
касательной к границе раздела и имеет постоянное значение:
V = шЯ, ш = const.
(4)
Необходимо найти такую угловую скорость ядра ш , чтобы момент всех сил, действующих на ядро, при выполнении граничных условий (3), (4) был близок к нулю. Для этого на основе решения краевой задачи (2) - (4) с определенными конкретными значениями ш, строится итерационный процесс путем минимизации момента сил:
Здесь ат — касательное напряжение на границе с ядром 7 . Исследуемая область разбивается на конечные элементы, как показано на рис. 3. Предполагается, что в деформируемой области справедливы уравнения (2). Расчеты проводились для однородного материала = 1) при различных эксцентриситетах эллипса и размерах жесткого ядра. Напряжения вычислялись при помощи метода сопряженных аппроксимаций [8]. Момент сил, действующих на ядро, считался близким к нулю при |М| < 0.02. Полученные результаты сравнивались с решением, найденным в [7] на основе приближенной модели:
Ь = 0.9, 0.8, 0.7, 0.6, 0.5. Точками обозначены значения, полученные при решении задачи
где Н = Ь — Д. На рис. 4 кривые 1 - 5 являются графиками функции (5) при значениях
1.1 -
1.0
П
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 h
2
1
Рис. З.
Рис. 4.
на основе описанного алгоритма. Точность совпадения является высокой. Относительная ошибка не превышает 0.2 % при Ь = 0.9, 0.8; 0.5 % при Ь = 0.7; 1.5 % при Ь = 0.6, 0.5. Сравнение результатов проводилось при 0.1 < Н < Ь — 0.1, так как сетка конечных элементов при Н = 0 вырождается.
Таким образом, рассмотренные выше модели с переменной вязкостью и жестким ядром также указывают на возможность необратимого накапливания малых приливных деформаций. Количественный анализ и учет трехмерности требуют более сложных моделей и постановок.
Список литературы
[1] Вегенер А. Происхождение континентов и океанов. Наука, М., 1984.
[2] Штауб Р. Механизм движений земной коры в приложении к строению земных горных систем. Гл. ред. геологоразвед. и геофиз. лит., Л. — М., 1938.
[3] Надаи А. Пластичность и разрушение твердых тел. Т. 2, Мир, М., 1969.
[4] Бобряков А. П., Ревуженко А. Ф., Шемякин Е. И. О возможном механизме переноса масс Земли. Докл. АН СССР, 272, №5, 1983, 1097-1099.
[5] Седов Л. И. Механика сплошной среды. Т. 1, Наука, М., 1973.
[6] Сегерлинд Л. Применение метода конечных элементов. Мир, М., 1979.
[7] Ревуженко А. Ф. О приливном механизме переноса масс. Известия АН СССР. Физика Земли, №6, 1991, 13-20.
[8] Оден Дж. Конечные элементы в нелинейной механике. Мир, М., 1976.
Поступила в редакцию 5 октября 1998 г., в переработанном виде 5 ноября 1998 г.