CC BY
10
ГИАБ. Горный информационно-аналитический бюллетень / MIAB. Mining Informational and Analytical Bulletin, 2021;(3):83-100 ОРИГИНАЛЬНАЯ СТАТЬЯ / ORIGINAL PAPER
УДК 622.02 DOI: 10.25018/0236-1493-2021-3-0-83-100
РЕШЕНИЕ НЕОДНОРОДНОЙ НЕЛИНЕЙНОЙ ОСЕСИММЕТРИЧНОЙ ЗАДАЧИ ГОРНОГО ДАВЛЕНИЯ ПРИ НАЛИЧИИ ЗАПРЕДЕЛЬНЫХ ЗОН С ИСПОЛЬЗОВАНИЕМ ЛОГАРИФМИЧЕСКИХ ДЕФОРМАЦИЙ
Н.П. Немчин1, С.В. Ветров1
1 Забайкальский государственный университет, Чита, Россия, e-mail: nemchin_np@mail.ru
Аннотация: Приведено решение нелинейных задач определения напряженно-деформированного состояния (НДС) массива вокруг протяженной горизонтальной выработки. Для этого разработан вариант численного метода, представляющий собой комбинацию нелинейного программирования и явного метода конечных разностей. НДС массива определяет величину давления на контуре выработки и радиальное перемещение этого контура, которые зависят от размера зоны запредельного деформирования. Решается задача определения горного давления путем перебора значений размеров зоны запредельного деформирования. Учтены неоднородность массива в радиальном направлении от центра выработки, история нагружения. Приводятся методы определения параметров задачи для ослабленного массива. Дан пример расчета реальной выработки. Указано, что существует область оптимальных перемещений крепи с заданной несущей способностью. Дана оценка упрощающих предположений в задачах горного давления. Подчеркивается необходимость ограничения относительной объемной деформации массива вблизи выработки при решении задач горного давления. Проведено сравнение решений с использованием логарифмических деформаций и деформаций Коши. Использование логарифмических деформаций существенно уточняет картину НДС при наличии запредельной зоны. Ключевые слова: НДС массива, логарифмические деформации, зона упругих деформаций, зона запредельных деформаций, явный метод конечных разностей, история нагру-жения, физическая нелинейность, геометрическая нелинейность, протяженная выработка, горное давление.
Для цитирования: Немчин Н.П., Ветров С.В. Решение неоднородной нелинейной осе-симметричной задачи горного давления при наличии запредельных зон с использованием логарифмических деформаций // Горный информационно-аналитический бюллетень. -2021. - № 3. - С. 83-100. DOI: 10.25018/0236-1493-2021-3-0-83-100.
Solution of nonuniform and nonlinear axially symmetric problem on confining pressure in the presence of post-limiting zones using logarithmic strains
N.P. Nemchin1, S.V. Vetrov1
1 Transbaikal State University, Chita, Russia, e-mail: nemchin_np@mail.ru
© Н.П. Немчин, С.В. Ветров. 2021.
Abstract: The article describes solution of nonlinear problems on the stress-strain behavior of enclosing rock mass around a long horizontal tunnel. The dedicated numerical method combines nonlinear programming and the explicit method of finite differences. The stress-strain behavior of rock mass governs the pressure and the radial displacement of the tunnel boundary, which in their turn depend on the size of post-limiting strain zone. The problem on the confining pressure is solved by search of values of post-limiting strain zone sizes. The nonuniformity of rock mass in the radial direction from the tunnel center, and the stress history are included. The methods to define the problem parameters for a weak rock mass are presented. A case-study of a real-life mine excavation is described. It is shown that there exists an optimal displacement zone of roof support with the preset load-bearing capacity. The simplified assumptions in the confining pressure problem are assessed. It is required to limit relative volumetric strain in the adjacent rock mass of the tunnel in the problems on confining pressure. The solutions obtained using logarithmic strains and the Cauchy strains are compared. The use of logarithmic strains considerably improves the stress-strain pattern in the presence of post-limiting zone.
Key words: rock mass stress-strain behavior, logarithmic strains, elastic strain zone, post-limiting strain zone, explicit method of finite differences, stress history, physical nonlinearity, geometrical nonlinearity, long tunnel, confining pressure.
For citation: Nemchin N. P., Vetrov S. V. Solution of nonuniform and nonlinear axially symmetric problem on confining pressure in the presence of post-limiting zones using logarithmic strains. MIAB. Mining Inf. Anal. Bull. 2021;(3):83-100. [In Russ]. DOI: 10.25018/0236-14932021-3-0-83-100.
Введение
В одиночных горизонтальных выработках определение горного давления прежде всего необходимо для прогнозирования их устойчивости и выбора типа и параметров крепи. При этом в массиве может возникать неоднородность, оказывающая существенное влияние на величину горного давления и меняющая характер распределения напряжений вокруг выработки. Это видно из работ И.В. Баклашова, А.М. Алимжанова и др.
Логарифмические или истинные деформации учитывают историю нагру-жения, поэтому их использование позволяет получить корректный результат, когда конечные деформации складываются из серии промежуточных [1, 2]. Применение логарифмической зависимости для деформаций дает более точное решение по сравнению с линейными деформациями при больших перемещениях контура выработки, порядка десятка сантиметров. При этом решение
задачи усложняется в виду нелинейной зависимости между деформациями и перемещениями.
Из работ, использующих логарифмические деформации, остановимся на [3]. В этой статье используется условие прочности Мора-Кулона, в котором коэффициент сцепления и угол внутреннего трения меняются скачкообразно при переходе через радиус границы упругой и запредельной зон. Упругие постоянные при таком переходе не меняются. Учитывается геометрическая нелинейность. Учитывается начальное поле напряжений. Разрыхление материала учтено через константу, зависящую от угла внутреннего трения в соответствии с ассоциированным законом течения. Относительное объемное расширение не ограничено. Массив, вмещающий выработку, однородный. Метод решения — аналитический.
В статьях [4, 5] используется простая полуэмпирическая зависимость для опи-
сания разрыхления. Решение своеобразное, но очень упрощенное. Один из вариантов решения использует логарифмические деформации (Генки).
Представим массив вокруг выработки как две кольцеобразных зоны вокруг нее. Первая — зона запредельных деформаций — находится непосредственно вблизи выработки между радиусами га и г*. Вторая — зона упругих деформаций — расположена сразу за первой на расстоянии г' от центра выработки. Считаем, что на расстоянии гь (на контуре Ь) от центра выработки ее влиянием на массив можно пренебречь и что в зоне запредельных деформаций прочность массива на сжатие и относительное объемное расширение изменяются до предельных значений. Эти зависимости на рис. 1 |ст| = £ (|е|) и 0 = £ (|е|) взяты из [6].
Из рис. 1 видно, что область запредельного деформирования разбивается в общем случае на 3 зоны границами, соответствующими точкам А и В. При
этом следует учитывать, что последовательность этих точек при движении от контура выработки может меняться.
В литературе известны следующие упрощающие предположения:
• точки А и В совмещаются на рисунке по вертикальной линии;
• точка А совмещается по вертикальной линии с пределом прочности (скачок прочности);
• точка А совмещается по вертикальной линии с контуром выработки (обычно с предположением о прочностной неоднородности);
• точку В не учитывают (не ограничивают величину относительного объемного расширения).
Следует также обратить внимание на то, что в зависимости от конкретных условий одна или несколько зон, непосредственно примыкающих к контуру выработки, могут не образовываться.
Для решения рассматриваемой задачи будет использован вариант численного метода, представляющего собой комби-
а(8)
аппроксимация о(е)
8(E)
аппроксимация 0(e)
зоны
деформирования
Рис. 1. Зависимость напряжений и деформаций от радиуса: I — область упругого деформирования, II — IV — область запредельного деформирования; II — область уменьшения прочности, III — область предельного разрушения, IV — область предельного разрыхления
Fig. 1. Stresses and strains versus radius: I — elastic strain zone, II—IV — post-Limiting strain zone; II — strength reduction zone, III —ultimate stress Limit zone, IV — ultimate disintegration zone
Рис. 2. Расчетная схема. Дано: г , г", гь, условие прочности на г". Найти: p(r"), ua(r") и затем r"(ua). Начальное состояние массива (а); массив в равновесном состоянии после проведения выработки (6) Fig. 2. Analytical model. Given: г, г, г , strength condition at r". Find: p{r"), ua(r") and r"(ua). Initial condition of rock mass (a); equilibrium condition of rock mass after tunnel construction (b)
нацию нелинейного программирования и явного метода конечных разностей. Явный метод — это метод, который позволяет получить состояние системы в каждой последующей точке, зная ее состояние в предыдущей [7, 8]. Предложенный метод удобен для решения задач с нелинейностью и при решении задач с несколькими зонами, границы которых подлежат определению. В этом его преимущество перед МКЭ и обычным МКР. Но следует учитывать, что он разработан для симметричных задач.
Задача будет решаться в осесим-метричном приближении (см. рис. 2). Заданы: г, г, г , условие прочности на г
К =/(<*".<*сж (Отбудем искать р(г) и и (г), затем через них получим р(иа). Для определения этих функций переберем значения г и получим кривую горного давления р(«а), а также радиуса запредельной зоны г(иа), через который можно построить кривую давления при обрушении. Таким образом, предлагаемый метод расчета требует изменений в обычном порядке работы с расчетной схемой, которая отражена на рис. 2.
Зона упругих деформаций
Допредельная зона упрощенно принята авторами как зона упругих деформаций.
Под влиянием горного давления неупругая зона может: а) возникать, и тогда рассматривается упругая зона с границами г и г. б) не возникать, и тогда рассматривается упругая задача с границами га и г. Метод, использованный в этой статье, состоит в комбинации нелинейного программирования и явного метода конечных разностей. Применим его для описания НДС массива с использованием логарифмических деформаций. Сжимающие напряжения будем считать отрицательными.
Приведем выражения для вычисления дифференциалов по времени от деформаций
d, s0 =
ш
d'£r =ТТ\
ЕЛг)
(dt<Jr -Mt°e),
где
ЕЛг) = ~
Е(г)
, Mi =
1-ц '
причем Е(г) — модуль деформации, зависящий от радиуса для неоднородного массива, ц — коэффициент Пуассона.
Для получения приращений деформаций интегрируем эти выражения по времени от момента т = 0, предшествующего возникновению выработки, т.е., когда напряжения сг = а0 = -5 , и = 0 (см. рис. 2, а) до момента, когда выработка уже существует, а деформации стабилизировались во времени. В результате получаем следующие выражения: 1
Дее = —-[сте + 5 - ц(стг + 5)],
Е1(г Г 1
Дег = у^) [аг + 5 5)], (1)
где Дв0, Двг — приращения деформаций вследствие проведения выработки.
Скорости деформаций и в ла-гранжевых координатах в осесиммет-ричной задаче задаются выражениями:
е _ ^ е _ V
Ьг _ ~Т~, Ье _ _ > йр р
где V — скорость, р — координата деформированного массива.
Эти формулы можно записать в виде:
¿¿и й{и
гг = —, ^ ее = —
¿Р Р
где р = и + г — радиус после проведения выработки; г — расстояние от центра выработки в момент т = 0, предшествующей ее возникновению, т.е. когда сг = а0 = -5. Таким образом, и — приращение перемещений при переходе от первоначального напряженного состояния (характеризуемого координатой г) к конечному состоянию (характеризуемому координатой р). То есть, и — видимые перемещения, вызванные проведением выработки.
Так как г не меняется со временем, то можно внести эту величину под знак дифференциала:
dtd (и + г) dt (и + г)
d,б, = -/, dtее = -
После интегрирования по времени от момента, предшествующего возникновению выработки, до момента т = ^ получим:
Дег = 1п |> (и + г)]т=( - 1п |> (и + г
Дее= 1п (и + г)т=( - 1п (и + г)т=о откуда, учитывая, что и| = 0 ,
Де = 1п\а—
Дее = 1п—
л г - (2)
аг ] г
Этот способ учета истории нагруже-ния с помощью уравнения типа (1) принят также в [3] для логарифмических деформаций (2) и в [9] для деформаций Коши. Приращения деформаций Дв0, Двг — логарифмические функции, зависящие от г и р. Таким образом, история нагружения в нелинейной постановке учитывается иначе, чем в линейной [10, 11]. В линейном случае также может быть использован предлагаемый метод учета истории нагружения. Приращения деформаций тогда представляют собой разности конечных и начальных деформаций Коши. Нами была решена неоднородная упругая задача в этой постановке с коэффициентом Пуассона ц = = 0,5. Результаты совпали с теми, которые приведены в [11], где также ц = 0,5, что говорит о корректности метода.
Определяющие уравнения для зоны упругих деформаций: р = и + г,
Де = 1п | ^
йг
Дее = 1п
Де =-
Де =
5а,
ЕМ (г) 1
ЕЩ
а - а
(( + 5 -Й1 (ае+ 5)) ((+ 5 ( + 5))
d (и + г )
и + г
5р
(3)
г
где р — радиус после проведения выработки; г — радиус до проведения выработки (когда напряжения сг = а0 = -5); Дв0, Двг — приращения радиальных и тангенциальных деформаций соответственно, вызванные проведением выработки; Е'м(г) — величина Е^г), уменьшенная для учета трещиноватости массива и длительности деформирования; ц1 — коэффициент Пуассона, причем модуль упругости и коэффициент Пуассона скорректированы для плоского деформированного состояния; 5 — естественное давление в массиве на глубине залегания выработки. Зависимость Е'м (г) существенна для учета технологических воздействий, например, буровзрывных работ.
Граничные условия для упругой задачи:
а
стг = -5, аг = -р , для упругой зоны
<УЬГ = -5, СУ* =
/(а*.<7* (г'))!
Для вычисления целевой функции используется явный метод конечных разностей. На радиус в области г, ^ г ^ гь наложим сетку узлов / = 1...А/. Узел с номером 1 лежит на контуре выработки в упругой задаче и на радиусе г, в задаче с упругой зоной. Узел с номером N лежит на контуре Ь. Для каждой точки повторяем вычисления для уравнений в указанном ниже порядке:
2. Д8д = 1п
где р — давление на контуре выработки; /— функция, входящая в условие прочности.
Выполнение граничных условий будет достигаться за счет минимизации целевых функций (ф и ф2).
Для упругой задачи целевая функция
°в=Ем('0Аее-5 + ^(<+5)>
5.М = ехр(Дв:),
6.^Ж-1,
дг дг
-, ;+1 г ди'
7. и =и +-Дг,
8.
дг
дг
дг
да.
9. ст'+1=с'+^Дг, г дг
ф, (Е1) = (СУ* + 5)2 -> тIп (4) 10. = ссж (П ),
6)
где ог — значение радиальных напряжений на контуре Ь; х1 = из — приращение перемещений на контуре выработки. Для упругой зоны целевая функция
Ф2(х1,Х2) = (су'+5)2 +
(5)
где индексом * обозначены значения величин на радиусе г", х1 =и,, х2 =су*. Точное значение минимума обоих целевых функций равно нулю.
где Дг — расстояние между соседними узлами; а (г) — функция, описывающая прочность массива с учетом технологических воздействий, их =х1; а) =х2 для упругой зоны и а]=-р для упругой задачи. Для вычисления величин и иав точках / + 1 используются только вычисленные ранее значения в точке / (явный метод). После этого выявляются все переменные НДС в точке / + 1 с помощью определяющей системы уравнений.
После того как будут получены значения всех величин в точках 1...Ы, вы-
числяется целевая функция (4) — для упругой задачи или (5) для задачи с упругой зоной.
Изменяя величины х1 и х2 от некоторых произвольных начальных значений (каждый раз при этом идет вычисление целевой функции по последовательности уравнений (6) для всех точек), можно добиться минимизации целевой функции по какому-либо методу нелинейного программирования. При этом удовлетворяются как определяющие уравнения, так и граничные условия.
Для иллюстрации вышеизложенного приведем два примера. В этих примерах используем конкретные выражения для описания Е'м(г) и асж(г). Модуль упругости, учитывающий влияние буровзрывных работ на массив [10,11]:
.У^
Е' fr)= Е
СМ V I I
1-я
K'aJ
(7)
где Ем — модуль упругости массива, скорректированный для плоского деформированного состояния; а, п — параметры, описывающие влияние буровзрывных работ.
Для описания зависимости прочности массива на одноосное сжатие от расстояния до контура выработки используем выражение [10,11]:
1-6
г
г
K'aJ
(8)
нородности:а выработки г -
Пуассона ц = 0,5. Результаты решения изложены в табл. 1.
Пример 2. Величины s, EfM, а, п, г, ц взяты из примера 1. Радиус упругой зоны г" = 2г, условие прочности
= + 6,773а,. -1,425а,2 а^ = = 5,964 M па; параметры неоднородности: b = 0,9; к = 3,007. Результаты решения изложены в табл. 2.
Решение для случая деформаций Ко-ши получено на основе [12]. Решение для случая логарифмических деформаций получено по формулам, приведенным в этой статье. Значение ц = 0,5 взято для того, чтобы сравнить с решением [11], где оно подразумевается.
Выводы из анализа таблиц:
1. Совпадение предложенного решения с решением [11] получено с погрешностью меньше 2%. Практическое совпадение результатов говорит о корректности предложенного численного метода решения.
2. Влияние геометрической нелинейности на решение упругой задачи и задачи для упругой зоны незначительно (менее 2,5%).
Зона запредельных деформаций
Считая величину г" заданной, опишем зону запредельных деформаций уравнениями:
р = г + и,
где <зсж — предел прочности на одноосное сжатие массива горных пород; Ь, к — параметры, описывающие влияние буровзрывных работ.
Пример 1. Расчетное давление по внешнему контуру выработки 5 = = 12,06 МПа. Выработка не закреплена, так что а" = 0. Модуль упругости массива Е;м = 1963 МПа, параметры неод-
Ае. = Ln
Де8 = Ln
dp
\dr у
V
= 0,941, п = 3,007. Радиус 1900 мм, коэффициент
AQ = f1{ar, Дее), ае =/2(аг,Д9 асж(г)), Д0 = Дег +Де8, dar _ а9 - аг dp р
(9)
vo о
Таблица 1
Упругая задача Elastic problem
1 2 3 Различия
Предложенное решение. Логарифмические деформации Предложенное решение. Деформации Коши Решение [11] Относительная разница, % колонки 1 и 3 Относительная разница, % колонки 1 и 2
r/r G, МПа Г сте, МПа U, мм а, МПа Г сте, МПа U, мм а,МПа Г сте, МПа U, мм <5 Г и <5 Г и
1 0 -2,3114 -37,578 0 -2,317 -38 0 -2,278 -37,4 0 1,7 1,7 0 0,3 1,2
1,3 -2,6688 -15,8346 -28,765 -2,611 -15,919 -29,2 -2,581 -15,672 -28,8 1,2 1,6 1,7 2,2 0,5 1,6
2 -7,6509 -16,1495 -18,598 -7,584 -16,234 -19 -7,466 -15,979 -18,7 1,6 1,6 1,6 0,9 0,5 2,0
Таблица 2
Упругая зона. НДС на границе запредельной зоны Elastic zone. Stress-strain behavior at post-limiting zone boundary
Предложенное решение (деформации Коши) Предложенное решение (логарифмические деформации) Относительная разница, %
r'lr ' а МПа а'е, МПа и", мм а", МПа Г а'е, МПа и", мм а" г и"
1,8 -2,174 -20,52 -49,2 -2,1584 -20,405 -48,6 0,7 0,6 1,3
где Д0 при малых значениях (до предельного значения 0,12, используемого в этой статье) приближенно равно изменению относительного объемного расширения вследствие проведения выработки и характеризует степень разрыхления и разрушения материала в процессе деформаций; £ — функция, описывающая дилатансию горных пород; £ — функция, входящая в уравнение состояния, зависящее от степени разрыхления Д0.
Пусть Д0* — приращение относительной объемной деформации на пределе прочности (со стороны запредельной области), принятое в дальнейшем равным нулю. При этом считаем, что эта деформация состоит из необратимой части растяжения Д0*р и упругой части сжатия Д0*е. Причем:
де* = де*' + де*р = о.
Это соответствует скачку на графике аппроксимации на рис. 1.
Для непрерывности материала желательно, чтобы радиальные перемещения, а следовательно и Дв0 не менялись при переходе через границу запредельной зоны. Тогда скачок в Д0* обуславливает скачок в т.е. в Де* производной йи/йг.
Условие непрерывности напряжений и деформаций при переходе из упругой зоны деформаций в запредельную:
а = а„
а0 = а0
где индексом 1е обозначен первый узел в зоне упругих деформаций, т.е. узел, лежащий в упругой зоне на радиусе г*; т — последний узел в зоне запредельных деформаций на радиусе г* (нумерация с контура выработки). Кроме того, из Д0* = 0 следует, что Де™ = -Де^ .
Приведем алгоритм вычисления НДС для запредельной области деформирования.
Вычислить для точки с = т лежащей на радиусе г* в запредельной области:
1. гт = г ,
т '
3.
др" дг
ди"
2. ^ = ехр(Де"),
др" дг дг 4. р = и" + г ,
г т т'
да
-1,
5.
дг
ай - а.
Л
др дг
Вычислить для каждой точки радиуса с = т...2 используя явный метод конечных разностей, где НДС в точке ¿ — 1 определяется через НДС в точке ¿. Точка с индексом 1 лежит на контуре выработки.
1 '-1 ' ди' л
1. и = и--Дг,
дг
2. Дг,
г г дг
3. р.-1 = и-1 + г
4. Де;-1 = 1л -
5. ДО'-1 = /1 (а*-\ Де;-1),
6. Де'-1 = ДО'-1 - Де!-1,
8.
9.
аО-1 = / 2 К
др'-1 дг ехр(
ди -1 др'-1
дг дг
даг-1 Га
-1,
дг
- 1Л
др
дг
(10)
Перебирая значения г* с некоторым шагом и решая последовательно задачи для упругой и запредельной зоны для каждого значения г*, получим соответствующие этим значениям величины: Р (г* ) = -а;1 (г*) — давление на контуре выработки; иа [г') = -и1 [г— перемещение контура выработки. Из этих за-
-1
г
-1
-1
висимостеи можно получить решение задачи в виде р(и_).
Конкретизируем функции/,^ и/2 соответственно:
условие прочности
■■ 0,75(0,0168 + 0,7ЪЪкос +
Д0 =
дилатансионное состояние (1 - ^ - 2с27стг ) (Леэ - Ае;), если Д@<©„„
(И)
(12)
А© , если А© >0
пр ' пр
уравнение состояния в запредельной области
-асж(г) + ТАв + С1аг+с2а2г, если =°сЖ(г)-ТАв>ап
(13)
если а = а (7)-7~Д0<а
ост еж \ / П)
+ 0,119к]с)-^-2 1 — ц
где с , с2 — коэффициенты предельной поверхности и уравнения состояния; сп — предельное значение остаточной прочности массива; ©п — предельная величина относительного объемного расширения (см. рис. 1); Т= Ц(с1 — 1) — параметр, характеризующий разрушение массива вследствие разрыхления; /. — модуль спада. В дилатансионном соотношении в грубом приближении можно принять с[=с1 и с'2= с2 что соответствует статье [12].
Следует обратить внимание на то, что вычисления явным методом с зависимостями (12) и (13) учитывают задачи, где некоторые из зон (см. рис. 1) отсутствуют. Правильная последовательность точек А и В получается тоже автоматически.
Определение параметров
Величина 0,75 приблизительно учитывает длительность деформирования. Зависимость от коэффициента Пуассона ц в знаменателе учитывает плоский характер задачи. Зависимость от коэффициента ослабления кос получена на основе экспериментальных данных с использованием формул для трещиноватых массивов [10].
Определение параметров с1 и с2
1. По пределам прочности на одноосное растяжение и сжатие а и а , определенных лабораторным путем, строится огибающая наибольших кругов Мора, например огибающая М.М. Протодья-конова.
2. Определяется предел прочности в массиве^ =кдлкососж через коэффициент длительной прочности и коэффициент ослабления, строится соответствующий круг Мора.
3. Огибающая, построенная в п. 1, переносится параллельно вниз до касания с кругом Мора, построенным в п. 2.
4. Определяются с1 и а3 для ряда предельных кругов Мора (огибающей из пункта 3) с центрами между |стс| = 5 и к1 = <*1/2 • Диапазон ст7ж/2<|сс| <5 приближенно соответствует диапазону изменения напряжений в запредельной зоне.
5. По этим значениям методом наименьших квадратов определяются параметры с и с2 для формулы:
М 2
= + СЛ + С2°1
Пункты 1 — 5 реализуются с помощью авторской программы pasp-coefs.
Параметр Е'м
Параметр определяется по формуле [13]:
Определение параметра Т В опытах на стабилометре с1 = а2. Это очевидно для образцов прямоуголь-
ного сечения. Для образцов цилиндрической формы примем это равенство для запредельной области в качестве гипотезы. Тогда уравнение состояния задается двумя функциями.
|-М 2
h = -ъсж + qcj-L + с2аг + +ГА0 - ст3 = О
г- М 2
h = ~ЪсЖ + С1°2 + С2<*2 + +ТА0 — су, = О
(15)
где Д0 = 0 — 0,.
Из принципа максимума диссипации энергии Э = кЪ,рк следует:
dFi
dF
да.
да,
(16)
да2 да2 даъ даъ
где X , X — множители Лагранжа.
Задача математически симметрична относительно функций Р и (напомним, что ст = ст ).
Поэтому множители Лагранжа
^ = = ^
Из третьего уравнения системы (16) определим
к = -
Тогда
dF. dF2
—L +—-
даъ даъ
1 +
dF
v да1 да2
' dF. 3F2
—L +—-
v5g3 daZJJ
После подстановки F и F2 из (15) следует
dtQp =(l-q -2c2al)dt
pp fb3
Проинтегрируем это уравнение по времени от предела прочности до текущего момента. При этом учтем, что опыты в стабилометрах осуществляются при постоянном боковом напряжении а, = const.
Тогда
= (е3'-Е;')(1 -c.-lc.a,).
Перейдем к полным деформациям:
е* = < + < ■
Тогда
е - es - е* + e*s =
= (е3 - е3 - б; + e;s)(l - Cl - 2с2а1)
Опыты [6, 14] показывают, что в запредельной области модуль упругости уменьшается вместе с остаточной прочностью. При этом упругую часть деформации можно считать почти постоянной еек « const в запредельной области.
Тогда 0s = 0*s, е'„=ек' и
0-0*=(е3-е;)(1-с1-2с2а1) (17)
После подстановки этого выражения в (13) можно получить as=-al+r(es-E;)(l-c1-2c2a1) +
Коэффициент хрупкости L = после дифференцирования:
да.
де.
L = Т 1-Cj -1с2аЛ
В случае одноосного сжатия ст = О и тогда
~Г = 1-г-
Практика расчетов показывает, что С1> 1 и тогда
Т = — (18)
При одноосном сжатии массива авторы предлагают следующую полуэмпирическую формулу (на основании аналогичной формулы из [13])
L = 2,5
0,5-ц сх-1
0,5-0,2 Cl+1
Для массива можно вычислить параметр 7 после подстановки этого значения /. в формулу (18).
Параметры с[ и с2 дилатансионного соотношения Вариант 1. В статье [12] было показано, что для путей нагружения со^
= COnSt можно принять с[ = С1 и cj = с2.
В рассматриваемой задаче не постоянно. Тем не менее желательно сохранить формулу, введенную в [12] для плоской деформации:
Р =
= Ра;-ст 1 + sincp
AG = (Лее - Ае*е)(1 - с1 -2c2<j1 ) + AG'
(19)
1 - sirup
где ф — угол внутреннего трения.
поскольку опыты указывают на зависимость 0 от аг Тогда соотношение (19) можно принять как грубо приближенное. Возникающие при этом погрешности должны корректироваться введением в теоретический расчет эмпирических параметров, например коэффициента перегрузки [15].
Вариант 2. При выводе дилатансионного соотношения принять линейную аппроксимацию уравнения состояния.
Тогда в (19)
п 1 + зтф с2 = 0. с1=р = --
А 1 — БШф
где ф — угол внутреннего трения. Величина А0* мала, и ею можно пренебречь. Деформация А^ в уравнении (12) соответствует моменту прохождения границы запредельной зоны г* через точку с координатой г. Численные опыты показали, что е^ мало зависит от величины г, поэтому ее можно определить на границе упругой зоны также для точек в запредельной зоне.
Угол внутреннего трения
После определения напряжений на границе запредельной области (из решения для упругой зоны) угол внутреннего трения определяется из системы:
Пример расчета горного давления для реальной выработки
Решение задачи горного давления будет складываться из множества решений в упругой и запредельной областях для набора задаваемых значений г. В результате будет получено множество значений напряжений на контуре выработки <(/■ ) и перемещений контура выработки и" (г*).
Используя эти значения, получим зависимость
М=р(0>
которая представляет собой кривую горного давления.
Приведем пример расчета с параметрами выработки «Квершлаг № 1 горизонт 750» Новоширокинского рудника [16]. Ширина выработки 3,8 м; А г = = 5,7 мм; гь = Юг, число точек в методе конечных разностей N = 3000. Радиус запредельной зоны (для рис. 3 — 6) г = 2г .
Естественное вертикальное давление на глубине залегания выработки 5В = 6,7 МПа; модуль деформации Е = = 16125,1 МПа; коэффициент Пуассона ц = 0,26; лабораторные пределы прочности: на сжатие а = 56,8 МПа, на
сж ' '
растяжение а = 9,8 МПа; выработка проведена буро-взрывным способом с расстоянием между оконтуривающими шпурами 50 см; предельный уровень разрыхления 0п = 0,12 (такое значение используется в полуэмпирическим методе для рассматриваемых горно-геологических условий [15]); коэффициент структурного ослабления массива кж = = 0,14; коэффициент длительной прочности кд = 0,75; плотность пород — 2959 кг/м3.
Рис. Ъ. Напряжения Fig. 3. Stresses
Расчетным путем были получены:
• радиус выработки, равный ее полуширине ra = 1,9 м;
• длительный модуль деформации структурно ослабленного массива для плоской деформации Е'м = 1579 МПа;
• длительный предел прочности на сжатие массива = 5,964 МПа
(°сж = кдлкососж)',
• параметры, описывающие неоднородность массива (вычисленные по
методике [9]): п = 3,01; а = 0,94; к = 3,01; b = 0,9;
• зависимости и <тсж(г) из формул (6) и (7) соответственно;
• предельное значение остаточной прочности массива <зпр =
• коэффициенты предельной поверхности и уравнения состояния сх = 5,238, с2 = 0,122; коэффициенты дилатанси-
онного состояния массива как с
j - ^ .
5 150
100
50
i ¡форм, м. Коши
2 1 2 — — дефор
Г. [Га]
Рис. 4. Радиальные перемещения (по модулю) Fig. 4. Radial displacement (by modulus)
1 3
1
2 — — деформ . Коши 1 .. 2
1.0
1.5
2.0
Г. [Га]
2.5
3.0
3.5
4.0
Рис. 5. Остаточная прочность Fig. 5. Residual strength
• расчетное давление по внешнему контуру s (s = кпе sß), где коэффициент перегрузки кпе = 2, по исследованиям [15]; sB — вертикальное напряжение до проведения выработки.
• параметр, характеризующий разрушение массива вследствие разрыхления, определенный по формуле (18): T= = 472 МПа.
Метод расчета реализован в программе «Горное давление в одиночных
горизонтальных выработках с учетом физической и геометрической нелиней-ностей» на языке программирования Python. Результаты показаны на следующих графиках (рис. 3—8).
Анализ графиков
В приведенном примере (рис. 3) абсолютные значения радиальных напряжений на контуре выработки составляют 176 кПа и 92 кПа для деформаций
1.2
1.0
3 0.6
0.4
0.2 AfV, пр
0.0
I 2--деформ. Коши \ 3 — деформ. Коши; Дб не ограничено
\ 3
1
2
1.0 1.2 1.4 г, [Га] 1.6
Рис. 6. Относительная объемная деформация Fig. 6. Relative volumetric strain
2.0
\l Yi « « \\ \\ u \v \\ \\ 1- дав/ 2- дав/ 3---дав/ 4----- дав/ \ , \ \ \ 1ение р, лог юние в случ 1ение р, де( тние в случ деформ. ае обрушени )орм. Коши ае обрушени Я пород Роб-Я пород Роб- лог. деформ деформ. Кои НИ
\ \ \ \ \ \ лхзч
----
100 200 300 400 5
Перемещение контура выработки, иа [мм]
Рис. 7. Кривые горного давления и давления в случае обрушения пород Fig. 7. Curves of confining pressure and rock failure pressure
Коши и логарифмических деформаций соответственно. Различия на 84 кПа играют существенную роль при выборе несущей способности крепления выработки.
Приращения перемещений (рис. 4) на контуре выработки составляют 224 и 251 мм для деформаций Коши и логарифмических деформаций соответственно.
Массив разрушен на протяжении почти всей запредельной области (оста-
точная прочность а = а ), и только
г ост пр'
вблизи области упругих деформаций остаточная прочность резко возрастает (см. рис. 5).
График (рис. 6) относительной объемной деформации построен по формуле &у - е0 -1. Формула следует из
Л0 = Лег + Аее = 1п (1 + Ае^ ,
где Asv
AdV dV
, причем dV — объем бес-
конечно малого элемента в начальном
1-лог. деформ. 2 --■ деформ. Коши У
у f *
* / V / у /
г s у
/J
200 300 400 500 600 700 Перемещение контура выработки, иа [мм]
Рис. 8. Радиус запредельной зоны Fig. 8. Post-limiting zone radius
(перед проведением выработки) напряженном состоянии. Предельный уровень Дву для рассматриваемой выработки распространяется примерно на треть запредельной области, для деформаций Коши и логарифмических деформаций. Если не ограничивать предельный уровень объемной деформации, то Дву равно 1,26 (Д0 = 0,82) для деформаций Коши, что в условиях массива невозможно. Построить кривую объемной деформации для логарифмических деформаций Дву в этом случае не удалось. Но эта кривая идет выше кривой, построенной по деформациям Коши. Из этого следует, что без ограничения объемной деформации решение будет некорректным.
Таким образом, из упрощающих предположений, перечисленных для рис. 1, оправдалось только одно: изменение прочности при переходе через границу запредельной области близко к скачку.
Давление на крепь при устойчивом характере деформаций падает с увеличением перемещений контура выработки, а давление на крепь при потере устойчивости в форме обрушения пород в кровле — возрастает (см. рис. 7). Последнее вычислялось по формуле Роб = У(г> — О, где У — удельный вес пород в кровле выработки. Различия порядка десятков кПа присутствуют и на кривых горного давления р (см. рис. 7), причем соответствующие значения для кривой, полученной с применением логарифмических деформаций, по абсолютной величине меньше. С увеличением радиальных перемещений контура выработки для решения,
СПИСОК ЛИТЕРАТУРЫ
где деформации учтены как линейные, радиус запредельной области растет быстрее, чем для решения с логарифмическими деформациями (см. рис. 8).
Кривая горного давления р(иа) как функция позволяет прогнозировать давление на контуре выработки, что является основанием для принятия решения о способе крепления выработки. Рассмотрим, какие выводы в отношении крепи следуют из анализа кривых, построенных для логарифмических деформаций. Пусть несущая способность крепи 200 кН. Если характеристика крепи такая, что равновесие наступает в области правее перемещений 170 мм или левее перемещений 400 мм, то разрушения не будет. Пусть несущая способность 100 кН. Такое значение было в реальной выработке [16]. Давление на крепь превышает эту несущую способность. Для рассматриваемой крепи при этой несущей способности действительно были разрушения [16].
Выводы
1. Использование логарифмических деформаций существенно уточняет картину НДС при наличии запредельной зоны.
2. Необходимо учитывать как геометрическую нелинейность, так и нелинейность уравнений состояния.
3. Рассмотренный численный метод решения позволяет решать такие задачи.
4. Метод расчета давления на крепь позволяет устанавливать оптимальный диапазон перемещений контура выработки.
1. Работнов Ю. Н. Сопротивление материалов. - М.: Физматгиз, 1962. - 456 с.
2. Rees D. Basic engineering plasticity: An introduction with engineering and manufacturing applications. Butterworth-Heinemann, 2006. 528 p.
3. Qiang Zhang, Cheng Li, Ming Min, Binsong Jiang, Liyuan Yu Elastoplastic analysis of circular openings in elasto-brittle-plastic rock mass based on logarithmic strain // Mathematical Problems in Engineering. 2017. DOI: 10.1155/2017/7503912.
4. Литвинский Г. Г. Запредельное поведение пород вокруг горной выработки (порождающее решение) // Сборник научных трудов Донбасского государственного технического университета. - 2017. - № 6(49). - С. 5-14.
5. Литвинский Г. Г. Статика разрушения и деформирования пород вокруг горной выработки // Сборник научных трудов Донбасского государственного технического университета. - 2017. - № 7(50). - С. 19-30.
6. Шашенко О. М., Сдвижкова О. О., Гапеев С. М. .Деформовашсть та мщшсть масивiв прських порщ: Монографiя. - Д.: Нацюнальний прничий ушверситет, 2008. - 224 с.
7. Smith G. D. Numerical solution of partial differential equations: finite difference methods. 1985, 337 p.
8. Гулин А. В., Самарский А. А. Численные методы. - М.: Наука, 1989. - 432 с.
9. Qiang Zhang, Chuanhu Zhang, Binsong Jiang, Ning Li, Yingchao Wang Elastoplastic coupling solution of circular openings in strain-softening rock mass considering pressure-dependent effect // International Journal of Geomechanics. 2018. Vol. 18. No 1. DOI: 10.1061/ (ASCE)GM.1943-5622.0001043.
10. Баклашов И. В. Геомеханика. Т. 1. - М.: МГГУ, 2004. - 208 с.
11. Баклашов И. В., Картозия Б. А., Шашенко А. Н., Борисов В. Н. Геомеханика. Т. 2. -М.: МГГУ, 2004. - 249 c.
12. Немчин Н. П., Ветров С. В. Решение задачи горного давления с нелинейной функцией разупрочнения явным методом конечных разностей // Известия вузов. Горный журнал. - 2015. - № 2. - С. 49-58.
13. Немчин Н. П., Терентьев П. Ю. Полуэмпирические формулы определения модуля деформации и модуля спада массива горных пород // Горный информационно-аналитический бюллетень. - 2016. - № 11. - С. 305-314.
14. Алимжанов А. М. Модель неоднородного упругопластически деформируемого массива горных пород с подземной выработкой // Механика и технологии. - 2014. -№ 1. - С. 6-12.
15. Немчин Н. П., Терентьев П. Ю. Оценка давления на крепь с помощью программы для ЭВМ «Горное давление в одиночных горизонтальных выработках-3» // Вестник Забайкальского государственного университета. - 2013. - № 1(92). - С. 32-38.
16. Немчин Н. П., Жувак А. С. Геомеханическая характеристика капитальных и подготовительных выработок ОАО Ново-Широкинского рудника как основа для тестирования методов расчета давления на крепь // Аспирант. Приложение к журналу Вестник Забайкальского государственного университета. - 2011. - № 2(10). - C. 106-114. итш
REFERENCES
1. Rabotnov Yu N. Soprotivlenie materialov [Strength of materials], Moscow, Fizmatgiz, 1962, 456 p.
2. Rees D. Basic engineering plasticity: An introduction with engineering and manufacturing applications. Butterworth-Heinemann, 2006. 528 p.
3. Qiang Zhang, Cheng Li, Ming Min, Binsong Jiang, Liyuan Yu Elastoplastic analysis of circular openings in elasto-brittle-plastic rock mass based on logarithmic strain. Mathematical Problems in Engineering. 2017. DOI: 10.1155/2017/7503912.
4. Litvinskiy G. G. Post-limiting behavior of enclosing rock around an underground excavation (generating solution). Sbornik nauchnykh trudov Donbasskogo gosudarstvennogo tekhnicheskogo universiteta. 2017, no 6(49), pp. 5-14. [In Russ].
5. Litvinskiy G. G. Statics of deformation and failure of surrounding rock mass around an underground excavation. Sbornik nauchnykh trudov Donbasskogo gosudarstvennogo tekhnicheskogo universiteta. 2017, no 7(50), pp. 19-30. [In Russ].
6. Shashenko O. M., Sdvizhkova O. O., Gapeev S. M. Deformovanist ta mitsnist masiviv girs'kikh porid. Monografiya [Деформовашсть та мiцнiсть масивiв гiрських порiд: Monograph]. D., Natsional'niy girnichiy universitet, 2008, 224 p.
7. Smith G. D. Numerical solution of partial differential equations: finite difference methods. 1985, 337 p.
8. Gulin A. V., Samarskiy A. A. Chislennye metody [Numerical methods], Moscow, Nauka, 1989, 432 p.
9. Qiang Zhang, Chuanhu Zhang, Binsong Jiang, Ning Li, Yingchao Wang Elastoplastic coupling solution of circular openings in strain-softening rock mass considering pressure-dependent effect. International Journal of Geomechanics. 2018. Vol. 18. No 1. DOI: 10.1061/(ASCE) GM.1943-5622.0001043.
10. Baklashov I. V. Geomekhanika. T. 1 [Geomechanics. Vol. 1], Moscow, MGGU, 2004, 208 p.
11. Baklashov I. V., Kartoziya B. A., Shashenko A. N., Borisov V. N. Geomekhanika. T. 2 [Geomechanics. Vol. 2], Moscow, MGGU, 2004, 249 p.
12. Nemchin N. P., Vetrov S. V. Solution rock mass on confining pressure with nonlinear function of weakening by explicit method of finite differences. Izvestiya vysshikh uchebnykh zavedeniy. Gornyyzhurnal. 2015, no 2, pp. 49-58. [In Russ].
13. Nemchin N. P., Terent'ev P. Yu. Half-empirical formulas of modulus of deformation and modulus of stress relaxation in rock mass. MIAB. Mining Inf. Anal. Bull. 2016, no 11, pp. 305314. [In Russ].
14. Alimzhanov A. M. Model of nonuniform elastoplastic deformation of surrounding rock mass around an underground excavation. Mekhanika i tekhnologii. 2014, no 1, pp. 6-12. [In Russ].
15. Nemchin N. P., Terent'ev P. Yu. Estimation of pressure on mine support using computer program Confining Pressure in Single Tunnels-3. Vestnik Zabaykal'skogo gosudarstvennogo universiteta. 2013, no 1(92), pp. 32-38. [In Russ].
16. Nemchin N. P., Zhuvak A. S. Geomechanical description of permanent and temporary roadways in Novo-Shirokino Mine as a framework for testing computational techniques for pressure on mine support. Aspirant. Prilozhenie k zhurnalu Vestnik Zabaykal'skogo gosudarstvennogo universiteta. 2011, no 2(10), pp. 106-114. [In Russ].
ИНФОРМАЦИЯ ОБ АВТОРАХ
Немчин Николай Павлович1 - канд. техн. наук, доцент, e-mail: nemchin_np@mail.ru,
Ветров Сергей Владимирович1 - старший преподаватель, e-mail: eridani_alpha@hotmail.com, 1 Забайкальский государственный университет, Для контактов: Немчин Н.П., e-mail: nemchin_np@mail.ru.
INFORMATION ABOUT THE AUTHORS
N.P. Nemchin1, Cand. Sci. (Eng.), Assistant Professor,
e-mail: nemchin_np@mail.ru,
S.V. Vetrov1, Senior Lecturer,
e-mail: eridani_alpha@hotmail.com,
1 Transbaikal State University, 672039, Chita, Russia.
Corresponding author: N.P. Nemchin, e-mail: nemchin_np@mail.ru.
Получена редакцией 27.09.2019; получена после рецензии 14.04.2020; принята к печати 10.02.2021. Received by the editors 27.09.2019; received after the review 14.04.2020; accepted for printing 10.02.2021.
