Проблемы идентификации неоднородных свойств твердых тел

Вестник СамГУ — Естественнонаучная серия. 2007. №4(54).

УДК 539.3

93

ПРОБЛЕМЫ ИДЕНТИФИКАЦИИ НЕОДНОРОДНЫХ СВОЙСТВ ТВЕРДЫХ ТЕЛ

© 2007 А.О. Ватульян1

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

Введение

Математическая модель однородной изотропной теории упругости, сформированная в XIX веке, была основана на положении, что свойства упругости характеризуются двумя упругими постоянными (модуль Юнга и коэффициент Пуассона, которые определяются из простых макроэкспериментов). Эта модель сыграла огромную роль в становлении на научную основу расчетов на прочность, устойчивость и колебания различных конструкций и сооружений [1]. Учет же затухания в реальных конструкциях обычно основывался на простых моделях вязкоупругости (типа Максвелла, Кельвина-Фойхта, Зинера), которые характеризуются также некоторым конечным (3-4) набором постоянных. В то же время ряд новых задач естествознания и техники требует отказа от гипотезы однородности для адекватного описания изучаемого объекта. При этом физические свойства (модули Ляме или компоненты тензора упругих постоянных, а также плотность среды) уже задаются при помощи функциональных зависимостей, которые должны быть предварительно определены из некоторых экспериментов. Наиболее распространенным способом определения таких зависимостей является анализ колебаний исследуемого объекта при варьировании способа нагружения и частоты колебаний, приводящий к исследованию довольно сложных нелинейных некорректных проблем.

1 Ватульян Александр Ованесович (vatulyan@math.rsu.ru), кафедра теории упругости Южного федерального университета, 344090, Россия, г. Ростов-на-Дону, ул. Мильчакова, 8а.

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

Отметим, что к коэффициентным обратным задачам теории упругости, в которых требуется определить коэффициенты дифференциальных операторов по информации о граничных полях смещений, приводят три типа задач.

Первый тип—собственно коэффициентные задачи, в которых требуется найти модули Ляме (или один из них) и плотность как функции координат по измерению поля смещений на границе. Этот тип задач, инициированный в первую очередь проблемами геофизики, исследовался достаточно давно, начиная с классической работы Герглотца [2], в которой построено решение одномерной обратной задачи для шара. Далее это направление развивалось на базе уравнений нестационарной теории упругости, и в рамках которой решен ряд обратных одномерных задач для полупространства [3, 4].

Второй тип — это задачи, к которым приводятся геометрические обратные задачи об определении форм полостей или включений малого характерного размера для стержней и пластин; для таких задач коэффициенты дифференциальных операторов отличны от постоянных значений лишь в некоторой малой подобласти исследуемого тела и их требуется восстановить по информации о граничных полях смещений. Отметим, что наиболее часто такие постановки встречаются для стержней и пластин с переменной жесткостью [5].

Третий тип коэффициентных обратных задач — об определении структуры существенно неоднородного предварительного напряженного состояния — практически неисследован, хотя в литературе имеется достаточно публикаций о расчете скоростей в предварительно напряженных телах при однородном поле предварительных напряжений в канонических областях типа слоя, бесконечного цилиндра [6, 7] и о влиянии предварительного напряженного состояния на резонансные частоты конструкций [8]. Отметим, что многие из исследуемых искусственных и природных конструкций находятся в условиях неоднородного предварительного напряженного состояния, особенно важен учет этого обстоятельства для задач идентификации биологических тканей и органов в рамках моделей линейной теории упругости и вязкоупругости. Кроме того, задача контроля предварительного напряженного состояния сооружений является одной из важнейших технических проблем, имеющей приложения в совершенствовании методик неразрушающего контроля, в горной механике и геофизике, биомеханике, а также в определении уровня предварительных напряжений в аварийных конструкциях [9]. Отметим, что акустические методы исследования, наряду с интерференционными, являются главными источниками получения информации о внутренних напряжениях и служат инструментом идентификации их структуры и уровня. Колебания упругих тел в условиях неоднородного

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

Область приложения задач об определении коэффициентов дифференциальных операторов по амплитудно-частотным зависимостям граничных полей постоянно расширяется. Если первые коэффициентные обратные задачи сформулированы в нестационарной постановке и посвящены в первую очередь проблемам геофизики и сейсморазведки [3, 4, 10-12], то в последние годы область приложения этих задач постоянно расширяется. Это проблемы неразрушающего контроля элементов конструкций, идентификации новых композиционных материалов сложной структуры, задачи эластогра-фии в медицинской диагностике мягких тканей, мониторинг имплантантов и контроль скорости восстановления костной ткани в месте перелома, проблемы акустического контроля при создании функционально-градиентных материалов, заменивших во многих технических приложениях слоистые композиты. Отметим, что значительная часть работ, выполненных в этой области, опирается на моделирование в рамках кусочно-однородной или слоистой структуры исследуемого объекта ввиду особой важности моделей такого типа в приложениях, и приводят решение исходной некорректной задачи к определению конечного числа параметров в некоторой ограниченной области и-мерного пространства. Такой поиск в последние годы осуществляется на основе метода регуляризации на компактных множествах [13], а среди конечномерных вариантов отметим как традиционные градиентные методы нахождения минимумов функционалов невязки [14], так и нейросети и генетические алгоритмы.

Главная трудность при исследовании коэффициентных обратных задач, которые являются нелинейными, состоит в достаточно сложной процедуре построения операторных соотношений, связывающих искомые и измеряемые в эксперименте функции. Это обусловлено переменностью коэффициентов дифференциальных операторов и невозможностью построения в явном виде общих представлений решений для соответствующих операторов, как это имеет место для операторов с постоянными коэффициентами при помощи фундаментальных и сингулярных решений и формул Сомильяны. Особую сложность представляют случаи, когда модули резко меняются (на один или два порядка), имеют разрывы первого рода в достаточно малой области. В ситуации, когда коэффициенты дифференциальных операторов переменны, методы решения прямых задач опираются либо на аппарат интегральных уравнений Фредгольма второго рода (для стержней и пластин), либо на конечноэлементные технологии. Что же касается построения операторных соотношений в обратных задачах, то эту сложность можно обойти, формулируя итерационные процессы отыскания неизвестных функций, основанные на обобщении теоремы взаимности [16], и сводя исходную задачу к поэтапному решению интегральных уравнений Фредгольма первого рода

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

1. Постановка обратных коэффициентных задач различных типов

Рассмотрим установившиеся колебания с частотой ш ограниченной области V с границей Ь = ЬииЬа, а Пу — компоненты единичного вектора внешней нормали к Я. Сформулируем несколько типов коэффициентных обратных задач.

Задача 1 типа. Уравнения колебаний имеют вид:

агу у + рш2иг = 0, г = 1,2,3,

= СЩик,1, к, I = 12,3,

иг |5 и = ° аг]П] 1^ и = Рг

Здесь Сгу1к — компоненты тензора упругих модулей, являющиеся кусочнонепрерывными функциями координат и удовлетворяющие обычным условиям симметрии и положительной определенности. Сформулируем задачу определения коэффициентов дифференциального оператора теории упругости по информации

ЩЬ0 = Жх, Ш), Ш € [Ш1, Ш2],

соответствующей измерению поля перемещений на части границы Ь а, на которой осуществляется нагружение, и моделирующей эхо-режим в акустических методах неразрушающего контроля. Вместе с тем отметим, что если носитель нагрузки занимает часть Ь а, то на этой части необходимо измерять смещения в соответствии с (1.4).Отметим, что часто в постановке различных типов обратных задач считается известным поле смещений на части границы Ь а° с Ь а, свободной от нагрузок

иг|ьа° = Жх, ш), ш € [Ш1, Ш2].

Так, например, обратная задача об определении модуля Юнга в задаче для стержня при продольных колебаниях описывается следующей краевой задачей

V Е°

и

д х

Задача 2 типа. В рамках постановки задачи 1 будем считать, что тело V содержит внутреннюю полость У°, свободную от нагрузок, и задача

ее определения (фактически геометрическая обратная задача) может быть сформулирована как коэффициентная, в частности, в нахождении области Vo, внутри которой обращаются в нуль модули упругости и плотность. Так, например, таковой является задача для стержня с полостью малого размера, площадь поперечного сечения которого F(х) есть функция продольной координаты, при этом соответствующая постановка имеет вид

№)и')' + к2 F

к

дх

где по функции /(к) требуется восстановить функцию F(х) , которая характеризует положение полости и ее размер. Решение простейших задач в рамках этого подхода можно найти в [5].

Замечание 1. Аналогичным образом может быть проанализирована задача идентификации включения малого размера; при этом неизвестные коэффициенты являются кусочно-постоянными функциями.

Задача 3 типа. Будем считать, что в теле имеется неоднородное предварительное напряженное состояние, характеризующееся компонентами тензора напряжений аУ и компонентами вектора перемещений и(0). Будем считать, что колебания вызываются нагрузкой, приложенной на части Ьа, а часть Ьи закреплена. Линеаризованные относительно дополнительных смещений иг уравнения колебаний после отделения времени имеют следующий вид [17]:

Ту у + рш2иг = 0,

Ту = ат]

ату = Ст]'к!ик,1,

иг\ьи = ^ ТгИ]\ь а = Рь

где Тгу — компоненты несимметричного тензора, определяемого согласно (1.11).

Замечание 2. Отметим что при отсутствии предварительного напряженного состояния задача (1.10)—(1.13) переходит в постановку (1.1)—(1.3).

Наличие предварительного напряженного состояния приводит к изменению амплитудно-частотных характеристик точек тела V, которое может быть положено в основу процедуры идентификации начальных напряжений. Отметим, что акустические методы их идентификации в случае однородного начального напряженного состояния подробно обсуждались в литературе [6, 7], причем для определения его уровня достаточно измерять либо изменения скоростей упругих волн, либо изменения нескольких первых резонансных частот. К сожалению, такая методика оказывается неприемлемой в случае неоднородного предварительного напряженного состояния, поскольку даже процедура решения прямой задачи (1.10)—(1.13) требует исследования линейного дифференциального оператора в частных производных с переменными коэффициентами. Отметим, что сильно неоднородные поля возникают в окрестности концентраторов, таких как полости и

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

поля смещений в зависимости от частоты колебаний вида (1.4). Требует-

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

2. Формулировка обобщенного соотношения взаимности

Для получения аналитических соотношений, связывающих заданные и искомые функции, в обратных коэффициентных задачах, которые являются нелинейными, наиболее часто используется метод линеаризации [13, 14]. При этом для нахождения искомой операторной связи необходимо обращать дифференциальные операторы с переменными коэффициентами, что практически невозможно для произвольных видов зависимостей. Как правило, большинство авторов, исследуя различные задачи неоднородной теории упругости, осуществляют замену исходного дифференциального оператора на дифференциальный оператор, который может быть обратим аналитически (как правило, это операторы с постоянными коэффициентами). Отметим, что процедура линеаризации требует вычисления производных операторов по Фреше. В то же время возможно избежать этой громоздкой процедуры, используя обобщенное соотношение взаимности. Поскольку постановка задачи третьего типа наиболее обща, сформулируем это соотношение для нее. Пусть имеются два решения задачи (1.10)—(1.13), отвечающие различным предварительным напряженным состояниям, различным модулям упругости и плотностям. Основные характеристики задачи снаб-

о (1) (01) ^(2) (2) (02) ™(1) тд

дим соответственно индексами 1 и 2: и. , а). , 1:. и и. , а). , 1.. . Используя

г гу гу г гу гу

обычную технику для уравнения (1.10), имеем

или, используя преобразование первого объемного интеграла в (2.1) в поверхностный в соответствии с теоремой Гаусса-Остроградского, окончатель-

V

V

но получим следующее соотношение:

0 = Г(о<02)-о<01)) иУи™ёУ + Г(о<2У02)Л1) -о^У01^)ёУ+

J \ г] г] / г,] г,т ^ V т] г,т г,] т] г,т г,] )

У /Км-сть) икЦифУ + “2 /

У

£ рг(и<2 - иг1^ ёБ = 0, ш е [Ш1, Ш2].

У +Г р (и<2)_ и <^

Я о

Из этого соотношения при отсутствии предварительного напряженного состояния вытекает соотношение взаимности для задач первого типа, полученное в [16],

5 С - ] и^и^.ёУ + / рг (и!? - иР) ёБ +

+ш2/

У

а при постоянных модулях упругости и плотности (в задаче третьего типа) — следующее соотношение

(о<02) - о<01)) иУи™ёУ + [(и(.02) - и(.01)) С0 .И^-<2)'

‘ г,т 3 V г,т г,т ) тты к,1

Г(о(02) - о(01)) и^У2ёУ + Г (и(°2) - и<01А С0 .И^У)ёУ+

J \ т] т] / г,] г,т Л \ г,т г,т / т.к к,1 г,]

f рг

^ Рг I и) - и)

Б о

Отметим, что если в модели, описывающей предварительное напряженное состояние, пренебречь величинами и™ по сравнению с единицей [17] (что часто осуществляется в постановках [6]), то соотношение (2.4) упрощается и принимает вид

Г

I \ т] т] / г, ] г,т ! \г г/ ’ ь 1? ^ у /

У Б о

Построенные соотношения взаимности (2.2)—(2.5) позволяют сразу формулировать итерационные процессы во всех типах обратных задач, минуя процедуру решения задач первого приближения, которая получается традиционным способом, описанным, например в [5].

3. Построение итерационных процессов

На основе соотношения (2.2) и его упрощенных вариантов можно построить итерационные процессы, сформулировать последовательность линейных интегральных уравнений Фредгольма 1-го рода с суммируемым ядром.

з.1. Задача первого типа

В этом случае будем исходить из соотношения (2.3). Положим в нем

„(1) = „<«-!) „(2) = „<п-1) + ,(п) с<1) = с(га-1) с<2) = С<п-1) + Сп) и сохраняя в

и. — и. , и. — и. + и. , с.7 — с..»* , с.7 — с.7 + с.7 и, сохраняя в

г г 1 г г г 1 г]к1 г]к1 1 г]Ы г]к1 г]к1 1 Г

нем линейные по величинам п-го приближения слагаемые, с учетом дополнительного условия (1.4) получим следующее уравнение относительно С^:

/йУ ^Р1

V Б0

п = 1,2, ••• , N. При этом ядра соответствующих интегральных операторов выражаются через деформации, найденные на предыдущем этапе из реше-

Г'

ния прямой задачи с известным тензором упругости Сг^ .

Такой подход требует знания некоторого начального приближения, которое может быть определено в классе некоторых простых функций, например линейных, из условия минимума функционала невязки. Пример решения конкретных задач о восстановлении модулей упругости как функций координат можно найти в [16].

3.2. Задача второго типа

В этом случае будем также исходить из соотношения (2.3). Положим в нем

и. = и. , и. = и. + и. ,

г г ’ г г г ’

С] = с

г]к1 г]к1 г]к1 х]к1 х]к1 у 7

р(1) = р(п-1), р(2) = р(п-1) + р(п),

аналогично предыдущему получим

ГС(пк;икпГ1>и(п-1>йУ - ш2Гр(п>и(п-1)и(п-1>йУ+

3 г]к1 к,1 г,] 3 г г г

У У + Г Рг

/ Рг / - и(п 1^ йБ = 0, ш е [Ш1, Ш2].

Соотношение (3.3) можно трактовать как интегральное уравнение относительно компонент С(п> и р(п>(х), если предварительно решена прямая задача о нахождении полей смещений и деформаций внутри области У и на ее границе Б с упругими характеристиками С]Г)(х) и р(п-1>(х). Отметим, что подынтегральные выражения в объемных интегралах в (3.3) представляют собой по форме аналог удвоенной удельной потенциальной энергии деформаций и удельной кинетической энергии, в которой перемещения и деформации соответствуют п - 1 итерации(предыдущей), а модули— п итерации (последующей). Мерой выхода из итерационной процедуры является функционал

Ш2

(/ - и(п 1>) йБ йш,

Ш1 Б0

и если его значение становится меньше погрешности входной информации, то процесс необходимо остановить.

Поскольку в изотропном случае тензор упругих постоянных выражается через две функции Ляме Х(х) и ц(х), то уравнение (3.3) при условии постоянной плотности примет вид

dV+

^Рх\1; - и(п dS — 0, ш е [шь Ш2].

+

3.3. Задача третьего типа.

В этом случае будем исходить из соотношения (2.4). Положим в

(01) /га-1) (02) ,(п-1) , /га) (1) (га-1) (2) (га-1) , (га) т-т

нем О.. — Г.. , О.. — Г.. + V. , и. — и. , и. — и. + и. Положим

X] X] X] X] X] X X X X X

(01) (п-1) (02) (п-1) 1 (п) гг. /о ,ч

и) — V; , и; — V; + V: . Тогда, сохраняя в (2.4) линейные относитель-

но возмущений слагаемые, получим

Г г(пУй-1) u(ц-r)dV + О Cmiklu(;l:1)u(n71)dv+

I т] х, ] х,т I х,т к,1 х,]

V V Г

+ I Р;(/х- иХ0))dS — 0, ш е [ш1,Ш2].

So

Это уравнение можно использовать для определения поправок компо-

нент тензора предварительных напряжений г • по отношению к некоторому

т]

выбранному начальному состоянию, которое также может быть выбрано простейшим образом, например однородным.

Замечание 3. В общем случае отыскания предварительного напряженного состояния требуется определить 6 компонент симметричного тензора т]х) и одного интегрального уравнения типа (3.5) этого недостаточно, чтобы осуществить процедуру реконструкции; для формулировки дополнительных интегральных уравнений такого же вида необходимо либо изменить структуру функций нагружения р;, либо область приложения нагрузки S О, при этом получатся аналогичные уравнения с другими ядрами и правыми частями.

Уравнения (3.1), (3.3) и (3.5) порождают операторные уравнения Фред-гольма 1-го рода с вполне непрерывными операторами, при обращении которых необходимо использовать регуляризующую процедуру в той или иной форме при наличии априорной информации о том, что искомые функции неотрицательны. Множество поиска можно сузить также при наличии дополнительной информации о носителе неизвестных функций (полость, включение) или при наличии ограничения, что они зависят лишь от одной координаты, что часто встречается в приложениях. Как правило, интегральное уравнение или система сводятся к решению системы линейных алгебраических уравнений относительно узловых переменных на основе одной из проекционных схем, в частности, и конечноэлементной аппроксима-

ции [16]. При численной реализации процедуры обращения был использован метод регуляризации А.Н. Тихонова [18, 19]. Отметим также, что вопрос о выборе отрезка изменения частоты колебаний является весьма важным как с точки зрения обоснования единственности поставленных обратных задач, так и с точки зрения построения эффективных численных схем. Как правило, этот отрезок необходимо выбирать в нерезонансном диапазоне для обеспечения единственности [20].

Работа выполнена при финансовой поддержке РФФИ (код проекта 05-01-00734).

Литература

[1] Ляв, А. Математическая теория упругости / А. Ляв. - М.; Л., 1935. -674 с.

[2] Herglotz, G. Uber die Elastizitat der Erde bei Boruckshtiqung unterVariablen Dichte / G. Herglotz // Zeitschr. fur Math. und Phys. -1905. - Bd. 52. - No. 3. - P. 275-299.

[3] Яхно, В.Г. Обратные задачи для дифференциальных уравнений упругости / В.Г. Яхно. - Новосибирск: Наука, 1990. - 304 с.

[4] Кабанихин, С.И. Проекционно-разностные методы определения коэффициентов гиперболических уравнений / С.И. Кабанихин. - Новосибирск: Наука, 1988. - 168 с.

[5] Бочарова, О.В. Реконструкция малых полостей в упругих стержнях / О.В. Бочарова, А.О. Ватульян, Р.С. Жарков // Известия вузов, Сев. кавк. рег. Сер. естеств науки. - 2006. - №2. - С. 28-32

[6] Гузь, А.Н. Волны в слое с начальными напряжениями / А.Н. Гузь,

А.П. Жук, Ф.Г. Махорт. - Киев: Наукова думка, 1976. - 104 с.

[7] Гузь, А.Н. Упругие волны в сжимаемых материалах с начальными напряжениями и неразрушающий ультразвуковой метод определения двухслойных остаточных напряжений / А.Н. Гузь // Прикладная механика. - 1994. - Т. 30. - №1. - C. 3-17.

[8] Филиппов, А.П. Колебания деформируемых систем / А.П. Филиппов. -М.: Машиностроение, 1970. - 734 с.

[9] Остаточные напряжения в деформируемых твердых телах / Г.Н. Чернышев [и др.]. - М.: Наука, 1996. - 310 с.

[10] Романов, В.Г. Обратные задачи математической физики / В.Г. Романов. - М.: Наука, 1984. - 261 с.

[11] Bui, H.D. Inverse Problems in the Mechanic of Materials: An

Introduction / H.D. Bui. - CRC Press, Boca Raton, FL, 1994. -224 p.

[12] Isakov, V. Inverse problems for PDE / V. Isakov. - Springer-Verlag, 2005. -284 p.

[13] Денисов, А.М. Введение в теорию обратных задач j А.М. Денисов. -М.: Изд-во МГУ, 1994. - 206 с.

[14] Oberai, A.A. Solution of inverse problems in elasticity imaging using the adjoint method j A.A. Oberai, N.H. Gokhale, G.R.Feijoo jj Inverse problems. - 2003. - №19. - P. 297-313.

[15] Ватульян, А.О. Интегральные уравнения в обратных задачах определения коэффициентов дифференциальных операторов теории упругости j А.О. Ватульян jj Доклады Российской АН. - 2005. - Т. 405(3). -С. 343-345.

[16] Ватульян, А.О. Об итерационном подходе в обратных задачах теории упругости / А.О. Ватульян, А.Н. Соловьев // Экологический вестник научных центров ЧЭС. - 2006. - №1. - С. 23-29.

[17] Ватульян, А.О. К формулировке интегральных уравнений в проблеме идентификации предварительного напряженного состояния /

A.О. Ватульян jj Экологический вестник научных центров ЧЭС 2006. - №2. - С. 23-25.

[18] Самарский, А.А. Численные методы решения обратных задач математической физики j А.А. Самарский, П.Н Вабищевич. - М.: Едиториал УРСС, 2004. - 480 с.

[19] Тихонов, А.Н. Методы решения некорректных задач / А.Н. Тихонов,

B.Я. Арсенин. - М.: Наука, 1986. - 287 с.

[20] Ватульян, А.О. Обратные задачи в теории упругости: материалы IX Всероссийского съезда по теоретической и прикладной механике / А.О. Ватульян. - Нижний Новгород, 2006. - Т. 3. - C. 54-55.

Поступила в редакцию 15jVj2007;

в окончательном варианте — 15j Vj2007.

PROBLEMS OF SOLIDS HETEROGENEOUS PROPERTIES IDENTIFICATION

© 2007 A.O. Vatulyan2

In the paper various settings of inverse problems of solids heterogeneous properties identification (the Lame module, prestressed state) are analyzed by the aid of coefficient inverse problems of elasticity theory solution. On the basis of generalized reciprocity relation iterative processes (the Fredholm integral equations of first type are solved on each step) are formulated.

Paper received 15/V/2007; Paper accepted 15/V/2007.

2Vatulyan Alexander Ovanesovich (vatulyan@math.rsu.ru), Dept. of Theory of Elasticity, South Federal University, Rostov-on-Don, 344090, Russia.