CC BY
76
ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
2010 Математика и механика № 2(10)
УДК 539.3+539.4+622.833.5+622.834
А.В. Кузнецова, И.Ю. Смолин
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ МЕХАНИЧЕСКОГО ПОВЕДЕНИЯ ГОРНЫХ ПОРОД ВОКРУГ ВЫРАБОТКИ ПРИ РАЗНЫХ СКОРОСТЯХ
ПОДВИГАНИЯ ЗАБОЯ1
Для описания различия в шагах обрушения кровли очистных горных выработок при разных скоростях подвигания забоя предложено использовать вязкопластическую модель среды и условие обрушения кровли по максимальным накопленным неупругим деформациям. Приведено описание использованного варианта вязкопластической модели и численная постановка задачи, которая была решена методом конечных элементов. Представленные результаты выполненных модельных расчетов свидетельствуют о работоспособности и перспективности предложенного подхода для описания скоростной чувствительности в явлениях разрушения геосред.
Ключевые слова: горные породы, горная выработка, напряженно-деформированное состояние, вязкопластическая среда, скоростная чувствительность, обрушение кровли
Практика работы угольных шахт показывает, что в связи со сменой горнодобывающего оборудования, отличающегося большей производительностью, возросли скорости подвигания забоя. При этом увеличились также и шаги обрушения кровли, значения которых не соответствуют показателям, рассчитанным по действующим нормативным документам. Поэтому актуальной для обеспечения безопасности проведения подземных работ стала проблема пересмотра основ расчета напряженно-деформированного состояния в окружающем выработку полезного ископаемого горном массиве при высокопроизводительной выемке угольных пластов.
Для объяснения и описания этого явления в работе [1] предлагается использовать экспериментально установленную зависимость прочности горных пород от скорости их нагружения [2, 3]. Как известно, прочность горных пород растет с увеличением скорости нагружения. В работе [1] продемонстрировано, что принципиально такой подход позволяет при численном моделировании предсказать насколько увеличится шаг первичного обрушения кровли при возрастании скорости подвигания забоя. Однако это не единственный способ. Например, в работах [4, 5] используются эволюционные модели накопления повреждений, которые при разных скоростях подвигания забоя обеспечивают разные картины поврежденно-сти в кровле при численном моделировании и, как следствие, позволяют получить различные шаги обрушения кровли. Третий подход связан с учетом вязких свойств геосреды. Как показано в работах [6, 7], для предела текучести пластичных материалов, природа увеличения предела текучести с ростом скорости нагружения кроется в вязких свойствах материала. За счет того, что вязкие напряжения уменьшаются на разную величину за разное время, прочностные свойства материалов при разных скоростях нагружения будут отличаться. Аналогичная ситуация имеет место и для прочностных свойств горных пород.
1 Работа выполнена при поддержке РФФИ, грант № 10-05-00509.
Цель настоящей работы состоит в демонстрации того, что применение вязкопластических моделей позволяет описать и смоделировать численно увеличение шагов обрушения кровли очистных горных выработок при увеличении скорости подвигания очистного забоя.
Для этого проведем сравнение распределений напряжений и деформаций при изменении размеров горной выработки, полученных при разных скоростях подви-гания забоя. Решение такой задачи будет выполнено методом конечных элементов в двухмерной постановке в условиях плоского деформированного состояния.
Метод расчета и определяющие соотношения вязкопластической среды
В данной работе для расчета использовалась программа Tochnog [8], которая позволяет решать широкий класс задач, в том числе и геомеханики, методом конечных элементов. Для описания реакции геоматериалов в ней имеются такие модели материалов, как упругопластические модели Друккера - Прагера, Кулона -Мора, вязкоупругие и вязкопластические модели. Рассмотрим подробнее используемую вязкопластическую модель.
Для описания процессов неупругой деформации геоматериалов в некоторых случаях удобно воспользоваться математическим аппаратом теории упругопластического течения [9]. При этом, не вдаваясь в микроскопические механизмы неупругого поведения геоматериалов, под пластическими деформациями понимают любые необратимые деформации независимо от их природы. Условие в напряжениях для перехода от упругого поведения к пластическому задается функцией текучести foj), которой в пространстве напряжений соответствует некоторая поверхность. Если это условие выполнено, то тензор скорости пластических деформаций (или их приращений в инкрементальном подходе) определяется формулой, выражающей закон пластического течения (или закон градиентальности), как величина, пропорциональная градиенту от пластического потенциала g(aj), также являющегося функцией напряжений. Если функция текучести и пластический потенциал совпадают, то закон пластического течения называется ассоциированным, а в противном случае - неассоциированным.
В работе Д. Друккера и В. Прагера [10] предложена модель на основе ассоциированного закона течения со следующей функцией текучести (пластическим потенциалом):
f = g = aJi + J-Y . (1)
Здесь а - коэффициент внутреннего трения; Y - сдвиговая прочность материала, или сцепление, J1 и J2 - соответственно первый и второй инварианты тензора напряжений.
При использовании такой модели компоненты тензора скоростей пластических деформаций имеют следующий вид:
^ sij ^
adi] +~1=
1 2J.
(2)
Здесь \ - пластический множитель (неопределенный множитель Лагранжа), который вычисляется из условия, что напряжения находятся на поверхности текучести, 8,у - символ Кронекера, а ^ означают компоненты девиатора тензора напряжений. Таким образом, пластическая деформация имеет не только сдвиговой, но и объёмный характер, а между их характеристиками имеется следующая связь:
стей пластических деформаций. Физически обусловленной причиной появления объёмных изменений в ходе пластической деформации является возникновение микротрещин и микропор. Коэффициент дилатансии Л, связывающий объёмную
в этом случае 3 а и в силу используемого закона пластичности жестко привязан к коэффициенту внутреннего трения а.
Однако в таком случае величина объёмной пластической деформации оказывается существенно завышенной по сравнению с экспериментальными данными. Чтобы получить согласование с экспериментами по дилатансии, можно ввести второй независимый параметр модели - коэффициент дилатансии, например, с использованием неассоциированного закона течения и введением пластического потенциала в виде
В этом случае скорости пластических деформаций определяются выражением
а коэффициент дилатансии Л = 3р, т.е. определяется вновь введенным параметром модели.
Для случая плоской деформации можно ввести понятия угла внутреннего трения как соответствующую величину в законе Кулона. Формула, связывающая коэффициент а и угол ф внутреннего трения, имеет вид
Аналогичным образом можно ввести и угол дилатансии.
В случае применения вышеуказанных теорий пластического течения величина пластических деформаций определяется уровнем напряжений и пластическим множителем. Таким образом, скорость пластических деформаций не зависит от времени, материал является не чувствительным к скорости нагружения.
Для того чтобы учесть скоростную чувствительность, используются вязкопластические модели - это модели пластичности с учётом влияния скорости деформации. Они часто применяются для высокоскоростных условий деформирования и нестационарных процессов. В отличие от вязкоупругих моделей, влияние на напряжение в вязкопластических моделях оказывается не через упругие, а через пластические деформации. Одним из вариантов построения вязкопластических моделей является обобщение закона текучести путем замены пластического множителя на некоторую безразмерную функцию от функции текучести, умноженную на скоростной параметр.
В качестве примера приведем формулу для степенной модели, которая имеется в программе Tochnog и использовалась в расчетах:
Здесь п - параметр текучести (величина, обратно пропорциональная вязкости),
(3)
(4)
2бш ф
>/3(3 - БШ ф)
(6)
/гЄ - исходное значение функции текучести, р - параметр модели. Физической основой этой феноменологической модели является широко применяемый степенной закон ползучести Нортона для второй стадии ползучести (установившаяся ползучесть), который устанавливает пропорциональность скорости деформации напряжению в некоторой степени. В качестве пластического потенциала можно использовать функции для различных моделей пластичности, в частности и Друк-кера-Прагера.
В отличие от чисто пластических моделей в этом случае напряжения могут превышать на некоторое время предел текучести. Возникающие при этом вязкие напряжения постепенно уменьшаются со временем по мере развития пластических деформаций. Скорость релаксации вязких напряжений определяется параметрами модели.
Описание геометрических характеристик и физических свойств расчётной модели
Численное моделирование было проведено для модельного объекта, характерного для шахт Кузбасского региона. Размер расчетной области и краткая характеристика неоднородности горного массива по вертикали (глубине )представлены на рис. 1.
Давление вышележащих пород
/ Основная кровля \ 90 м /
Непосредственная кровля \ 6 м /
Пласт угля ^ ^4 м
/ Непосредственная почва ) ^10 м (
Основная почва \ \ 90 м !
Рис. 1. Схема расчетной области по глубине и характеристика вмещающих пород
Расчетная схема представляет собой геометрическую модель в плоской двухмерной постановке, включающую пласт угля и слои вмещающих пород, окружающие выработку полезного ископаемого. Разделение на слои обусловлено разным составом и физико-механическими свойствами горных пород. Сверху полезного ископаемого выделяют непосредственную и основную кровли, а снизу - непосредственную и основные почвы.
Физико-механические свойства угля и вмещающих пород, использованные в расчетах в данной работе, представлены в табл. 1. Они были выбраны на основе имеющихся в научной литературе данных [11 - 13]. Поскольку значений для параметров выбранной вязкопластической модели найти не удалось, то величиина
параметра текучести была выбрана из условия соответствия результатов расчётов по порядку величин с данными по шагам обрушения кровли из практики работы угольных шахт, а параметр р был принят равным единице.
Т аблица 1
Физико-механические характеристики пород
Название Плотность, г/см3 Модуль Юнга, ГПа Коэффициент Пуассона Сцепление, МПа Угол внутреннего трения, град Угол дилатан-сии, град Начальное значение функции текучести, МПа Параметр текучести, с-1
Основная кровля 2,52 20 0,14 0,6 38 13 6 10-12
Непосредственная кровля 2,41 20 0,14 0,2 38 13 2 10-12
Уголь 1,39 5 0,3 0,1 30 10 1 10-12
Непосредственная почва 2,42 20 0,14 0,6 38 13 6 10-12
Основная почва 2,5 20 0,14 0,6 38 13 6 10-12
Граничные условия, определяющие особенности нагружения горного массива, были выбраны следующие. На верхней грани расчетной области прилагается давление вышележащих слоев пород, соответствующее заданной глубине залегания угольного пласта. Нагрузка задается также действием силы тяжести во всей области. Нижняя грань закреплена в вертикальном направлении. Боковые грани -вертикальные оси симметрия (запрещено смещение в горизонтальном направлении). Кроме того, в процессе расчета в пласте угля появляется и растет в одном направлении горная выработка, что обуславливает изменение геомеханической обстановки в исследуемом объеме горного массива (рис. 2). В начальный момент времени очистная выработка отсутствует. В процессе расчёта часть конечных элементов убирается и тем самым задается образование очистной горной выработки и изменение ее размеров за счет продвижения забоя (рис. 3).
Рис. 2. Геометрическая модель расчётной области в разные моменты времени
Размер расчетной области составил 400 х 200 м. Расчётная сетка при отсутствии горной выработки - 134x200 ячеек. Использовались прямоугольные четырехугольные элементы. Примеры конечно-элементных аппроксимаций при разных размерах выработки показаны на рис. 3. Каждый конечный элемент вдоль горизонтальной оси X имеет размер 3 м, поэтому минимальное расстояние, на которое может измениться размер горной выработки, также составляет 3 м.
Рис. 3. Сетка конечных элементов в нетронутом горном массиве и при продвижении горной выработки
Результаты расчётов и их обсуждение
В программе Tochnog расчет ведется с явным указанием шагов по времени. То есть, хотя задача решается как квазистатическая и на каждом шаге определяется равновесное напряженно-деформированное состояние системы, но время отслеживается и используется при интегрировании определяющих соотношений. В расчете указывается, за какой период времени убираются конечные элементы. Таким образом, изменяя это время можно варьировать скорость подвигания забоя. Расчёты были проведены для двух таких скоростей: 30 и 60 м/мес.
При такой постановке задачи имеется возможность оценить влияние скорости подвигания забоя на шаг первичного обрушения кровли. Первичное обрушение кровли является более опасным по своим проявлениям, чем последующие обрушения, поэтому его исследование является и более важным с практической точки зрения. Для определения условия обрушения нужно выбрать соответствующий
критерий. В работе [14] при проведении подобного типа расчетов в качестве критерия обрушения пород кровли предложено условие, когда выбранный критерий прочности в напряжениях достигается в большей части конечных элементов, покрывающих кровлю. Однако в указанной работе не учитывались вязкие свойства среды. В нашем случае более подходящими следует признать критерии, основанные на развитии неупругих деформаций со временем или работе напряжений на таких деформациях.
Результаты проведенных расчетов представлены на рис. 4. Анализ распределения интенсивности пластических деформаций, рассчитанных по формуле єр = ^О.бвуву , где ву - компоненты девиатора тензора пластических деформаций,
показал, что в результате учета вязких свойств геосреды, при продвижении горной выработки на одно и то же расстояние с разной скоростью, в непосредственной кровле накапливаются разные неупругие деформации.
Рис. 4. Распределение интенсивности пластических деформаций в непосредственной кровле над выработкой (границы выработки показаны пунктирными вертикальными линиями) длиной 36 м (а) и 66 м (б) при разных скоростях подвигания забоя: 1 - 30 м/мес.; 2 - 60 м/мес.; 3 - 90 м/мес.
Если сопоставить данные расчетов и практики работы угольных шахт, то за критическое значение накопленных неупругих деформаций можно принять значение порядка 1,5-Ш-6.
Полученное в расчетах соотношение между шагами начального обрушения кровли при скоростях подвигания забоя 30 и б0 м/мес. хорошо согласуется с данными, приведенными для различных шахт Донбасса [15]. К сожалению, подобных сведений для Кузбасса авторам найти в научной литературе не удалось, хотя согласно неопубликованным данным оно близко к данным для Донбасса.
Заключение
Проведенные расчеты с применением выбранной модели вязкопластической среды и модельных параметров подтвердили возможность описать увеличение шага обрушения кровли при увеличении скорости подвигания забоя в угольных шахтах. Предложенный алгоритм оценки начального шага обрушения непосредственной кровли после привязки параметров вязкопластической модели к данным экспериментального исследования вязких свойств вмещающих пород и уточнения критерия обрушения кровли может быть использован в реальных условиях шахт Кузбасского региона.
ЛИТЕРАТУРА
1. Иванов А.С., Сдвижкова Е.А., Бабец Д.В. Численное моделирование влияния скорости обнажения горных нород на механические процессы вблизи сопряжения очистной и подготовительной выработки // Матеріали Міжнародної конференції «Форум гірників -2009». Т. 3. Донецк: Національний гірничий університет, 2009. С. 37 - 44.
2. Мансуров В.А., Протосеня А.Г. Поведение горных нород нри различных скоростях нагружения. Фрунзе: Илим, 1982. 88 с.
3. Красько Н.И. Исследование зависимости прочности нород на растяжение от времени их нагружения // Наукові праці ДонНТУ. Серія: «Гірничо-геологічна». 2002. Bm. 54. С. 154 - 159.
4. Макаров П.В., Смолин И.Ю., Евтушенко Е.П. и др. Моделирование обрушения кровли над выработанным пространством // Физ. мезомех. 2008. Т. 11. № 1. С. 44 - 50.
5. Макаров П.В., Смолин И.Ю., Евтушенко Е.П. и др. Сценарии эволюции горного массива над выработкой // Физ. мезомех. 2009. Т. 12. № 1. С. 75 - 82.
6. Макаров П.В., Жукова Т.В., Платова Т.М. и др. Исследование вязких и релаксационных свойств металлов в ударных волнах методами математического моделирования // ФГО. 1987. № 1. С. 29 - 34.
7. Макаров П.В. Сдвиговая прочность и вязкость металлов в ударных волнах // Ударные волны и экстремальные состояния вещества / нод ред. B.E. Фортова, ЛВ. Альтшулера, Р.Ф. Трунина, А.И. Фунтикова. М.: Наука, 2000. С. 219 - 255.
8. Руководство пользователя TOCHNOG [Электронный ресурс], 2004. - URL: http:// tochnog.sourceforge.net, свободный.
9. Стефанов Ю.П. Некоторые особенности численного моделирования поведения упругохрупкопластичных материалов // Физ. мезомех. 2005. Т. 8. № 3. С. 129 - 142.
10. Друккер Д., Прагер В. Механика грунтов и пластический анализ или предельное проектирование // Механика. Новое в зарубежной науке. Bbm. 2. Определяющие законы механики грунтов. М.: Мир, 1975. С. 1бб - 177.
11. Прочность и деформируемость горных нород / Ю.М. Карташов, БВ. Матвеев, r.B. Михеев, А.Б. Фадеев. М.: Недра, 1979. 2б9 с.
12. Каталог механических свойств горных нород нри широкой вариации видов напряженного состояния и скорости деформирования / А.Н. Ставрогин, Е.Ю. Семенова, B^. Авксентьева и др. - Л.: BE^^IH, 197б. 171 с.
13. Механические и абразивные свойства горных пород / под общ. ред. Л.А. Шрейнера. М.: Гостоптехиздат,1958. 202 с.
14. Хозяйкина Н.В. Закономерности изменения предельного напряженного состояния в сложноструктурной кровле лав пологопадающих угольных пластов: дис. ... канд. техн. наук. Днепропетровск, 2004. 144 с.
15. Іванов О.С. Аналіз факторів впливу на крок обвалення норід покрівлі лави в умовах високого ступеню метаморфізму норід // Наукові нраці ДонНТУ. Серія «Гірничо-геологічна». 2009. Bm.10(151). С. 148 - 151.
С B EДEНИЯ ОБ АBТ ОРАХ:
КУЗНЕЦОВА Анастасия Викторовна - аспирантка кафедры прочности и проектирования физико-технического факультета Томского государственного университета. E-mail: KuznetsovaA@sibmail.com
СМОЛИН Игорь Юрьевич - доктор физико-математических наук, доцент, профессор кафедры прочности и проектирования физико-технического факультета Томского государственного университета, старший научный сотрудник Института физики прочности и материаловедения Сибирского отделения РАН. E-mail: smolin@ispms.ru
Статья принята в печать 07.05.2010 г.
