И.А.Бригаднов
Прямые методы решения вариационной задачи.
УДК 539.3
ПРЯМЫЕ МЕТОДЫ РЕШЕНИЯ ВАРИАЦИОННОЙ ЗАДАЧИ ДЛЯ МНОГОКРИТЕРИАЛЬНОЙ ОЦЕНКИ НЕСУЩЕЙ СПОСОБНОСТИ ГЕОМАТЕРИАЛОВ
И. А. БРИГАДНОВ
Санкт-Петербургский горный университет, Санкт-Петербург, Россия
В статье рассматриваются прямые методы решения вариационной задачи в напряжениях для многокритериальной оценки несущей способности образца из геоматериала в текущей конфигурации, которая может быть как отсчетной (недеформированной), так и актуальной (деформированной). Поставленная задача состоит в минимизации интегрального квадратичного функционала от различных компонент напряжений в выбранной контрольной подобласти на множестве полей напряжений, статически уравновешенных с внешними воздействиями. Для простейших конфигураций образца предлагается использовать метод обобщенных рядов Фурье в гильбертовых пространствах. Для сложных конфигураций образца с концентраторами напряжений предлагается использовать конечно-элементную аппроксимацию с последующей минимизацией конечномерной квадратичной функции с линейными ограничениями равенств. Приводится содержательный численный пример по оценке несущей способности образца из геоматериала при чистом сжатии.
Ключевые слова: геоматериал; вариационная задача в напряжениях; несущая способность; многокритериальная оценка; обобщенные ряды Фурье; конечно-элементная аппроксимация
Как цитировать эту статью: Бригаднов И.А. Прямые методы решения вариационной задачи для многокритериальной оценки несущей способности геоматериалов // Записки Горного института. 2018. Т. 232. С. 368-374. DOI: 10.31897/РМ1.2018.4.368
Введение. Горные породы и бетон являются одними из основных конструкционных материалов, поэтому оценка их несущей способности является весьма актуальной научно-технической проблемой [3, 13]. В работе [1] автором предложен оригинальный подход, основанный на многокритериальной оценке несущей способности конечного образца из геоматериала в отсчетной (недеформированной) или актуальной (деформированной) конфигурации. В рамках этого подхода поставлена вариационная задача для напряжений в заданной подобласти, где в зависимости от инженерных соображений, оцениваются среднеквадратичные значения различных компонент напряжений и по их совокупности делается вывод о несущей способности текущей конфигурации тела по отношению к заданным внешним воздействиям. Предложенный подход можно отнести к нелокальным методам, в которых разрушение геоматериала рассматривается как физический процесс в некоторой окрестности концентратора напряжений [7], а не в материальной точке, где достигается максимальное значение напряжения. При помощи нового метода можно построить глобальные оценки снизу для внешних силовых воздействий, гарантированно разрушающих конечный образец из геоматериала, что крайне важно для процесса дробления горных пород [13].
Постановка задачи. Рассмотрим конечный образец из геоматериала как деформируемое твердое тело, занимающее в текущей конфигурации ограниченную липшицеву область QeR3. Тело находится в равновесном напряженно-деформированном состоянии под действием внешних стационарных воздействий: на непустой части границы Г1 с лебеговой мерой | Г1! > 0 задано фиксированное перемещение и, в О приложена объемная сила с плотностью / а на части границы Г = дО\Г приложена поверхностная сила с плотностью g. Тогда, согласно принципу виртуальной мощности, справедливо основное равенство статики [2, 5]:
П(П,Т0, V) = Л(П, /, g, V), Уу еГ (П), (1)
где П(П, Т0, V) = | Т0--(V® v)TdQ и Л(П,/, g, V) = | / ■ vdQ +| g ■ V dy - мощности внутренних
П П Г2
напряжений Т0 и внешних сил (/ g) на допустимых скоростях перемещений V е Г(П) = {у : П ^ Я3; х) = 0, х е Г1} соответственно. Здесь и далее точка обозначает скалярное произведение векторов, а две последовательные точки - двойное скалярное произведение (свертку); V = д/дх - символический оператор Гамильтона; символ ® - тензорное (диадное) произведение; верхний индекс Т обозначает операцию транспонирования [4]. В отсчетной конфигурации Т0 - первый (несимметричный) тензор номинальных напряжений Пиола - Кирхгоффа, а в акту-
Ж И.А.Бригаднов 001: 10.318971РМ1.2018.4.368
:,и Прямые методы решения вариационной задачи...
альной конфигурации Т0 - симметричный тензор истинных напряжений Коши [2, 5, 9]. Заметим, что актуальная конфигурация образца может быть получена в результате решения, например, общей задачи упруго-вязкопластического деформирования [5].
При оценке несущей способности деформируемого твердого тела по отношению к заданным внешним воздействиям необходимо учитывать особенности топологии текущей конфигурации, поскольку наиболее критичными являются зоны концентрации напряжений (отверстия, выточки, разрезы и т.д.). Поэтому выделим некоторую ограниченную липшицеву подобласть ю с О, содержащую концентратор напряжений. Воздействие оставшейся части тела на выделенную (далее контрольную) подобласть заменим следующими краевыми условиями: часть границы подобласти у1 с Г1 закреплена, а на оставшейся части границы у2 = бю/у1 задана внешняя сила пТ0, где п -внешняя по отношению к контрольной подобласти нормаль к границе у2. Заметим, что закрепленная часть границы у1 может отсутствовать.
Введем в рассмотрение множество допустимых скоростей перемещений в ю:
V(ю)= {V е W 1,2(ю,R3): v(х) = 0, х е у1}
и квадратичный функционал на гильбертовом пространстве напряжений L2(ю, М 3) вида
К(ю,Т) = -"^| ГТ-В■■TdQ, (2)
2 ю
I I ю
где В - тензор 4-го ранга, зависящий только от координаты х. Здесь и далее М3- пространство вещественных матриц 3 х 3.
В зависимости от выбора тензора В функционал К(ю,Т) имеет различный физический смысл:
• если В = Е - первый единичный тензор 4-го ранга [4], тогда К12 - среднеквадратичная интенсивность напряжений в ю;
• если в = Е -11 ® I, где I - единичный тензор 2-го ранга [4], тогда К12 - среднеквадратич-
3
ная интенсивность касательных напряжений в ю ;
1 1/2
• если В = — I ® I, тогда К - среднеквадратичная интенсивность гидростатического давле-
ния в ю;
если в = — | Е —— I ® I I - тензор жесткости материала, обратный тензору упругости, где
2ц ^ 1 + у
ц - модуль сдвига, а V - коэффициент Пуассона [2, 3, 5, 9], тогда К - средняя удельная внутренняя энергия деформированного твердого тела в ю .
В общем случае выбор контрольной подобласти ю и тензора В определяется инженерно-техническими соображениями.
Как и ранее, в контрольной подобласти справедливо основное равенство статики П(ю,Т= А(ю, f, пТ0, V) для любого допустимого поля скоростей перемещений vеV(ю), где функционал П определен в задаче (1), а функционал А имеет вид
А(ю, f, п Т0, V) = Г f ■ V dQ + Г пТ0 ■ vdу .
Поставим следующую задачу: среди всех статически равновесных полей напряжений в контрольной подобласти ю найти такое, которое минимизирует функционал К(ю,Т), а именно,
Т* = а^тДК(ю,Т): Т е G}, (3)
О = {Т е Ь2(ю,М3): П(ю,Т0,V) = А(ю, f,пТ0,V), Vv е V(ю)} .
И.А.Бригаднов
Прямые методы решения вариационной задачи.
С математической точки зрения в рассматриваемой задаче ищется минимум квадратичного функционала на линейном аффинном многообразии в гильбертовом пространстве напряжений L (ю, M ). Далее будем рассматривать только тензоры ©, удовлетворяющие сильному условию
I 12
Коулмана-Нолла, для которых существует постоянная а > 0 такая, что Q-- © • • Q > a Q для любого Q G M3 и почти всех x ею [2, 9]. В этом случае функционал К(ю, T) является строго выпуклым и коэрцитивным на L2(raM3) и поэтому гарантируется существование и единственность решения задачи (3) [6, 9].
Итак, в контрольной подобласти ю рассматриваются самые слабые поля напряжений, глобально уравновешенные в текущей конфигурации тела Q, т.е. предполагается, что в ю твердое тело сопротивляется внешним воздействиям самым слабым образом. В действительности тело может уравновешивать существенно более интенсивные внешние воздействия. Если величина К1/2 превышает заданный уровень напряжений, тогда топология подобласти ю не может обеспечить несущую способность текущей конфигурации тела по отношению к заданным внешним воздействиям, поскольку любой материал может выдержать только конечный уровень напряжений (касательных, обжимных, разрывных, комбинированных и т.д.). Таким образом, ищется оценка снизу для несущей способности текущей конфигурации твердого тела по отношению к заданным внешним воздействиям и контрольной подобласти. Предлагаемый подход можно отнести к нелокальным методам оценки разрушения геоматералов [7].
В работе [1] аналитически показана содержательность предлагаемого подхода на примере задачи о сферически симметричном нагружении твердой сферы внутренним давлением в отсчет-ной конфигурации.
Методы решения задачи. Особенностью вариационной задачи (3) является то, что допустимое множество G содержит первые производные искомых функций, а интегральный квадратичный функционал К(ю,Т) нет. Поэтому для ее решения удобно использовать прямые методы. Наиболее общим из них является метод обобщенных рядов Фурье в гильбертовых пространствах [10]. Однако он применим только для простейших областей, для которых можно относительно легко построить счетный базис. Для сложных областей наиболее приемлемым является метод конечных элементов [5, 8], который можно трактовать как частный случай конечных рядов Фурье. В этом случае в результате конечно-элементной аппроксимации исходная вариационная задача (3) сводится к простейшей конечномерной задаче минимизации квадратичной функции с линейными ограничениями равенствами, решение которой легко реализуется в среде инженерных расчетов MATLAB.
Метод обобщенных рядов Фурье. Выберем в гильбертовом пространстве Ь2(ю) счетный базис {фк(х)}ке^ в гильбертовом пространстве W 1,2(ю) счетный базис ^m(x)}meN такой, что фт(х) = 0 на у1, и в гильбертовом пространстве L2(y2) счетный базис {yr(x)}reN. Тогда справедливы следующие разложения в ряды Фурье [10] (здесь и далее используется правило суммирования по повторяющимся индексам от 1 до œ):
а) для допустимых напряжений Т(х) = T^k(x), где Tke M3 - вещественные матричные коэффициенты;
б) для допустимых скоростей перемещений v(x) = vmфm(x), где vm GR3 - вещественные векторные коэффициенты;
в) для заданных внешних сил f (x) = fk фк (x) и g(x) = gryr (x), где fkeR3 и g еЯ - соответствующие вещественные векторные коэффициенты.
В силу произвольности допустимых скоростей перемещений, т.е. векторных коэффициентов vm (meN), основная задача (3) принимает следующий вид:
T }keN = arginf {iT • • Kl] • • Tj : T'cim = bm }, (4)
где тензоры 4-го ранга Ây = Aj , удовлетворяющие сильному условию Коулмана - Нолла, определяются по формуле
И.А.Бригаднов DOI: 10.31897/РМ1.2018.4.368
Прямые методы решения вариационной задачи...
векторы с1т по формуле
с =
{Ф^фт ¿О
а векторы Ьт по формуле
т г*
ьт =
fk {Ф1 Фт¿О + ^ {уМ
г Т т
У2
Таким образом, задача свелась к минимизации сильно выпуклой и коэрцитивной квадратичной формы на линейном аффинном многообразии или гиперплоскости в бесконечномерном евклидовом (гильбертовом) пространстве, что обеспечивает существование и единственность решения [6].
Для решения задачи (3) воспользуемся методом множителей Лагранжа [6]. Введем в рассмотрение функцию Лагранжа
L(Tk, Г) = IТ1 • А • •Т1 + Г (Т1с1т - Ьт),
необходимое условие стационарности которой имеет вид
дТ дЬ
—= А* •-Т1 + Гт ® с*т = О, — = Т*с*т - Ьт = 0, (5)
Т дГ V >
где О - нулевой тензор 2-го ранга.
У каждого тензора 4-го ранга А11, удовлетворяющего сильному условию Коулмана - Нолла, всегда существует обратный тензор 4-го ранга В11, также удовлетворяющий сильному условию Коулмана - Нолла, такой, что В11 = В11 и В11 • -А11 = А11 • -В11 = Е - первый единичный тензор 4-го ранга [2, 4]. Тогда из первого условия (5) находим
Т* = -В*1 • •Г ® с1т = -Гт • В*1 • с1т .
Подставив это выражение во второе условие (5), получаем бесконечную систему линейных алгебраических уравнений (СЛАУ) относительно множителей Лагранжа Г. Из свойств тензоров 4-го ранга В11 следует невырожденность тензора 2-го ранга с'т • В11 • с1 [2, 4], поэтому полученная СЛАУ однозначно разрешима:
гт = -Ь (с1т • В1 • с1' )-1.
Таким образом, искомые матричные коэффициенты в разложении поля напряжений в ряд Фурье по базисным функциям гильбертова пространства Ь (о) имеют вид
Т* = -Ь (с1т • В1 • с1 )-1 • В*г • сгт . (6)
Например, если взять © = Е и ортонормировать базисные функции в Ь2(о) по процедуре Шмидта [10], тогда А11 = В11 = Е , где 5ц - символ Кронекера. В этом случае решение (6) принимает более простой вид:
Т* = -(с1т • с1' )-1 Ь' ® с*т .
К сожалению, построенное аналитическое решение может быть использовано только для простейших подобластей о, поскольку выбор счетного базиса в гильбертовых пространствах для произвольных областей является самостоятельной трудноразрешимой проблемой [6, 10].
Вариационно-разностный метод. Для решения исходной вариационной задачи (3) воспользуемся стандартной конечно-элементной аппроксимацией с последующим численным решением соответствующей конечномерной проблемы.
¿К И.А.Бригаднов 001: 10.318971РМ1.2018.4.368
:,и Прямые методы решения вариационной задачи...
Как правило, текущая конфигурация является либо отсчетной недеформированной, либо определяется численно в результате решения краевой задачи для всей области О в рамках какой-либо модели (для малых или конечных деформаций, для упругой или неупругой среды, для перемещений или напряжений). При этом используется стандартная конечно-элементная аппроксимация. Для области ОеЯп (п = 1,2,3) строятся множества Ои и Ги = сОи такие, что I О \ ОИ| ^ 0 и | Г \ ГИ| ^ 0 для И ^ +0 регулярным образом, где И - наибольший шаг триангуляции, а ОИ представляет собой совокупность простейших согласованных между собой
симплексов - отрезков, треугольников или тетраэдров для пространств размерности п = 1, 2 или 3 соответственно [8].
Выбор контрольной подобласти ш И для дальнейшего анализа сводится к выделению множества конечных элементов в окрестности концентратора напряжений с переносом краевых условий на границу подобласти. Далее можно воспользоваться более мелкой триангуляцией в контрольной подобласти. Каждая триангуляция характеризуется множеством узлов {Х}ке1, где I -множество номеров узлов в ШИ. Обозначим через II и 12 множества номеров узлов на частях гра-1 2
ницы у И и у И соответственно.
Для допустимых скоростей перемещений, допустимых напряжений и внешних сил воспользуемся стандартной кусочно-линейной непрерывной аппроксимацией вида
^И (X) = Vм Ф я (х), Ти (х) = Тк Ф к (х),
(7)
fh (х) = fk Фк (х), gh (х) = gr % (х),
где нижние индексы т е I \ 11, к е I и г е 12, а величины с верхним индексом обозначают соответствующие узловые значения. Здесь используются стандартные непрерывные скалярные функции формы [8]:
• линейная на каждом симплексе функция Ф к : шИ ^ Я такая, что Фк (х') = 5кг., финитный
носитель которой состоит из соседних симплексов с общей вершиной хк ;
• линейная на каждой грани соответствующего симплекса функция % : у2h ^ R такая, что %(х') = 8ti, финитный носитель которой состоит из соседних граней на уh с общей вершиной хг.
Хорошо известно [8], что линейное пространство, натянутое на линейно независимый базис (ф к }ке1, является конечномерным подпространством гильбертова пространства L2(rah), которое
всюду плотно в L2(rah) при h ^+0 регулярным образом. Поэтому аппроксимации (7) можно рассматривать как конечные ряды Фурье соответствующих полей, что позволяет все предыдущие рассуждения полностью перенести на конечномерный случай.
Задача (3) аппроксимируется следующей конечномерной проблемой:
(Tk }ке1 = arg min ( i T • • А11 • • T1 : T'c'm = bm }, (8)
где индексы i, j e I и m e I \ Ix. Здесь тензоры 4-го ранга А'1 = А1 вычисляются по формуле
А' = 71г Г ЮФ'.Ф,dQ, I щ ®h
векторы cim по формуле
c'm = |Ф'УФ m dQ
fflh
и векторы bm по формуле
bm = fk |фкФmddQ + gr J%%mdy .
®h Ущ
J\ И.А.Бригаднов DOI: 10.31897IPMI.2018.4.368
:,и Прямые методы решения вариационной задачи...
Преобразуем задачу (8) к виду, удобному для численного решения. Для этого соберем все матричные переменные {Tk }keI в глобальный вектор y размерности Mn2, где M =| I1 - число узлов в шй. Далее по стандартной методике строим глобальную матрицу G квадратичной функции размерности Mn2 х Mn2, глобальную матрицу Q линейных ограничений размерности Kn х Mn2, где K =| I \ /j | - число узлов в юй \ у\, и глобальный вектор d правых частей линейных ограничений размерности Kn. Тогда конечномерная задача принимает простой вид
y* = arg min {\ yTGy : Qy = d}, (9)
где индекс T - транспонирование. Заметим, что число переменных в задаче (9) всегда больше числа линейных ограничений равенств, поскольку M > K и n = 1, 2, 3.
Для тензора 4-го ранга © в функционале (3), удовлетворяющего сильному условию Коул-мана - Нолла, регулярная конечно-элементная аппроксимация с последующей стандартной сборкой гарантирует сильную положительную определенность матрицы G квадратичной функции [8], поэтому всегда существует единственное решение конечномерной задачи (9). При помощи метода множителей Лагранжа решение легко записывается в явном виде, удобном для реализации в системе инженерных расчетов MATLAB,
У* = U (QU) 1 d, (10)
где U = G-Q .
Вычислительные эксперименты. Эффективность вариационно-разностного метода решения основной задачи (3) легко показать на примере оценки несущей способности длинного бруса квадратного сечения при чистом сжатии вдоль продольных боковых поверхностей [12]. В рамках модели плоского деформированного состояния сечение бруса (единичный квадрат) разбивалось на 2N2 одинаковых треугольных конечных элементов, где N - число равных дроблений сторон квадрата. В качестве текущей рассматривалась деформированная конфигурация сечения бруса, найденная в результате решения следующей модельной задачи линейной теории упругости в вариационной форме для перемещений при жестком кинематическом деформировании [12]:
1 1
u = arginf {J(u): u eV} , J(u) = JJ(I2(V® u) + XI2(V® u))dxdy, (11)
0 0
V = {ueW 12(Q,R2): ux(0,y) = 0,ux(1,y) = Ux,uy(0,y) = uy(1,y) = 0,y e(0,1)},
где и = (их, иу) - плоский вектор перемещений, /¡(V ® и) = их х + иу у - интенсивность объемной деформации /2 (V ® и) = (| (и^ х - иххиуу + и2уу)+ 1и2х у)12 - интенсивность сдвиговой деформации
(первый и второй инварианты тензора малых деформаций Коши соответственно [2,9]). Здесь безразмерный параметр Х = 3(1 - 2у)/(1 + у) зависит от коэффициента Пуассона уе (0,05) [2, 9]; V -множество кинематически допустимых перемещений - таких, что левая сторона квадрата жестко закреплена, а правая кинематически жестко перемещается на величину их без вертикального
проскальзывания. Как и в работе [12], конечномерный аппроксимационный аналог задачи (11) решался в среде инженерных расчетов МАТЬАВ при помощи стандартной процедуры $т1псоп для числа дроблений сторон квадрата N = 20 .
Известно, что разрушение образца из геоматериала при одноосном сжатии (т.е. при их < 0) начинается на свободной боковой поверхности за счет отслаивания [3, 11, 13]. Поэтому в силу симметрии задачи в качестве контрольной подобласти была выбрана окрестность верхней свободной поверхности сечения бруса из трех слоев конечных элементов, где оценивались минимальные среднеквадратичные значения вертикальных (отслаивающих) напряжений по формуле (10).
Например, для параметра X = 2, отвечающего коэффициенту Пуассона V = 0,125 (гранит, мрамор [3, 13]), и перемещения их = -0,2 минимальные среднеквадратичные вертикальные напряжения Т*, найденные в результате решения задачи (3) в контрольной подобласти, в сред-
И.А.Бригаднов
Прямые методы решения вариационной задачи..
У 1
0,8
0,6
0,4
0,2
-0,2
У У У У У / /
У У ц У У У У У У / г ш У У У ■у У У У
у Щ / У У У У У: У / / У У У У / У У щ
у У / У щ У У У У f / У У У У Ш У / У У
у У у У У У У У / / У У У У / У У У
щ. У у У У У У У у / / У У У У м У ш А У
у 1 / У % У У У. I 1 / У У 1 У У. У У у У
у & У У У У У У я / / У У / У. щ У У У У
i ш / У У У У У / / / У У / % У У У У У
/ щ У У У У У У / / / У У и У 1 У у У У
% щ У У У У У У / / / У У У У / У У У У
У у У У У У У У / / / У У У У / У У У У
У 1 У У У У У У / / / У У У У ё i У У У
У У / У / / У У
нем на 18 % меньше таких же напряжений, найденных в результате решения полной задачи (11). Таким образом, если значение К112 окажется больше предельного напряжения хрупкого разрушения для рассматриваемого материала, то образец гарантированно разрушится [11]. На рисунке красным цветом показаны области, где может происходить потеря глобальной устойчивости материала.
В заключение еще раз отметим, что на основе предлагаемого многокритериального нелокального анализа несущей способности образца из геоматериала можно построить абсолютные оценки снизу для внешних воздействий, гарантированно разрушающих тело, что крайне важно для процесса дробления горных пород [13].
0
0,2
0,4
0,6
0,8 х
Деформированная конфигурация сечения бруса. Зоны возможного разрушения образца из-за отслаивания выделены красным цветом
Выводы
1. Многокритериальная оценка несущей
способности образца из геоматериала сводится к вариационной задаче минимизации интегрального квадратичного функционала от различных компонент напряжений в выбранной подобласти на множестве полей напряжений, статически уравновешенных с внешними воздействиями.
2. Прямой метод решения вариационной задачи на основе обобщенных рядов Фурье в гильбертовом пространстве применим только для простейших областей, для которых можно относительно легко построить счетный базис.
3. Прямой метод решения вариационной задачи на основе конечно-элементной аппроксимации применим для сложных областей с концентраторами напряжений и сводится к простейшей конечномерной задаче минимизации квадратичной функции с линейными ограничениями равенствами.
ЛИТЕРАТУРА
1. Бригадное И.А. Многокритериальная оценка несущей способности геоматериалов // Записки Горного института. 2016. Т. 218. С. 289-295.
2. Лурье А.И. Нелинейная теория упругости. М.: Наука, 1980. 512 с.
3. Николаевский В.Н. Геомеханика и флюидодинамика. М.: Недра, 1996. 447 с.
4. Пальмов В.А. Элементы тензорной алгебры и тензорного анализа. СПб: Изд-во Политех. ун-та, 2008. 109 с.
5. ПоздеееА.А. Большие упруго-пластические деформации / А.А.Поздеев, П.В.Трусов, Ю.И.Няшин. М.: Наука, 1986. 232 с.
6. СеаЖ. Оптимизация. Теория и алгоритмы. М.: Мир, 1973. 244 с.
7. Сукнее С.В. Применение нелокальных и градиентных критериев для оценки разрушения геоматериалов в зонах концентрации растягивающих напряжений // Физическая мезомеханика. 2001. Т. 14(2). С. 67-75.
8. Сьярле Ф. Метод конечных элементов для эллиптических задач. М.: Мир, 1980. 512 с.
9. Сьярле Ф. Математическая теория упругости. М.: Мир, 1992. 472 с.
10. ТреногинВ.А. Функциональный анализ. М.: Наука, 1980. 496 с.
11. Черепанов Г.П. Механика хрупкого разрушения. М.: Наука, 1974. 640 с.
12. Brigadnov I.A. Regularization of non-convex strain energy function for non-monotonic stress-strain relation in the Hencky elastic-plastic model // Acta Mechanica. 2015. Vol. 226. Iss.8. P.2681-2691.
13. VerruijtA. Computational geomechanics. Dordrecht: Springer Science+Business Media, B.V., 1995. 384 p.
Автор И.А. Бригаднов, д-р физ.-мат. наук, профессор, [email protected] (Санкт-Петербургский горный университет, Санкт-Петербург, Россия).
Статья поступила в редакцию 03.05.2017. Статья принята к публикации 02.03.2018.