Torsion of a growing shaft

Full Text

Введение. Современные методы изготовления деталей сложной формы подразумевают разнообразные технологические процессы обработки, в том числе связанные с удалением части материала - резка, фрезеровка, сверление и т. д. Другой класс технологических процессов основан на синтезе деталей путем последовательного добавления материала на подложку или поверхность произвольной формы. Использование новых технологий позволяет создать пространственную деталь весьма сложной формы из широкого класса материала. Однако в процессе роста в теле неминуемо возникают остаточные деформации. С этими деформациями связаны остаточные напряжения и их уровень может вызывать такие явления, как потеря устойчивости и разрушение, причем они могут произойти уже в процессе такого изготовления. Таким образом, разработка и развитие математических моделей и теоретических методов расчета внутренних технологических напряжений и деформаций в изделиях сложной формы является актуальной научной проблемой как с позиций фундаментальных исследований, так и с точки зрения многочисленных приложений. Фундаментальными вопросами поверхностного роста занимается механика растущих тел (см., например, [1-5]). В работе [1] развита теория поверхностного роста для исследования двумерных задач наращивания деформируемых тел в случае, когда скоростью деформирования поверхности тел за счет нагрузок и натяга приращиваемых элементов можно пренебречь по сравнению со скоростью притока нового материала к этой поверхности. Предложены методы решения новых краевых задач, основанные на приведении неклассических задач наращивания вязкоупругих стареющих тел к задачам теории упругости с некоторым параметром, использовании теории аналитических функций для решения последних и восстановлении истинных параметров напряженно-деформированного состояния тел (см., например, [6-10]). Решение каждой конкретной задачи наращивания деформируемых тел представляет собой самостоятельную и трудоемкую проблему (см., например, [11-13].) Однако уже сейчас по виду полученных математических соотношений можно предсказать такие присущие растущим телам явления, как возникновение остаточных напряжений после снятия нагрузок, появление в наращиваемом теле поверхностей разрыва напряжений, зависимость напряженно-деформированного состояния вязкоупругих тел от скорости и способа их изготовления. Установлено, что в двумерных задачах концентрации напряжений возле отверстий и кручения в готовом теле без учета процесса наращивания максимум интенсивности касательных напряжений достигается на границе тела. В процессе роста максимум интенсивности касательных напряжений может достигаться на границе раздела основного тела и наращиваемой части, на границе готового тела и в произвольной точке наращиваемой части тела. Полученные результаты могут служить основой при решении прикладных задач расчета деталей и элементов конструкций, изготавливаемых при помощи методик современной технологической практики. 685 М а н ж и р о в А. В., М и х и н М. Н., М у р а ш к и н Е. В. 1. Постановка краевой задачи кручения. Рассмотрим (риc. 1) достаточно длинный вал (круговой цилиндр) длиной 2l и радиусом b0 (отношение l к b0 довольно велико), изготовленный из вязкоупругого стареющего материала в начальный момент времени. К обоим торцам вала приложены плоские круглые диски радиуса a < b0 . В момент времени τ0 на диски начинают действовать крутящие моменты M (t), поворачивая их на угол γ = 2α(t) друг относительно друга. Поверхность вала свободна от напряжений (риc. 1). В момент времени τ1 к поверхности вала начинается приток вещества. Вновь присоединяемые элементы также не напряжены, а время их изготовления совпадает со временем изготовления первоначального тела. Закон роста вала полностью определяется функцией b(t), которая определяет изменение радиуса вала со временем. При этом b(τ1 ) = b0 . Рост прекращается в момент времени τ2 . В этот момент радиус вала равен b1 (b(τ2 ) = b1 ), а его боковая поверхность остается свободной от напряжений при t τ2 . Для решения контактной задачи растущего тела воспользуемся квазистатическим приближением в отсутствие объемных сил. Будем считать вал достаточно длинным как во время процесса роста, так и после его прекращения (в этом случае отношения l/b(t) и l/b1 должны быть довольно большими). Принимая во внимание симметричность постановки краевой задачи, рассмотрим только половину вала с одной торцевой поверхностью, зажатой на жестком основании, а c другой торцевой поверхностью, связанной с диском. Воспользуемся общим подходом для решения подобных задач механики растущих тел, разработанным и описанным в работах [1-10]. M (t) b1 b(t) 2l b0 M (t) z Рис. 1. Растущий вал в условиях кручения [Figure 1. Growing Shaft under Torsion] 686 Кручение растущего вала Рассмотрим систему фундаментальных соотношений краевой задачи в интервале времени t ∈ [τ0 , τ1 ], тогда для начального вязкоупругого стареющего вала получим 2σrϕ ∂σrϕ ∂σϕz + + = 0 (∇ · T = 0), ∂r ∂z r z = 0, 0 r = b0 , 0 εrϕ = r z a : uϕ = α(t)r, l : σrϕ = 0, 1 ∂uϕ uϕ , - 2 ∂r r εϕz = T = 2G(t)(I + L(τ0 , t))E, (1) z = 0, a r b0 : σϕz = 0, z = l, 0 r b0 : uϕ = 0; 1 ∂uϕ 2 ∂z 1 E = [∇u + (∇u) ] ; 2 (I - L(τ0 , t)) = (I + N(τ0 , t))-1 ; t f (τ )K1 (t, τ ) dτ, L(τ0 , t)f (t) = K1 (t, τ ) = G(τ ) τ0 ∂ 1 + ω(t, τ ) , ∂τ G(τ ) где T и E - тензоры напряжений и деформаций с ненулевыми компонентами σrϕ , σϕz и εrϕ , εϕz соответственно; u - вектор перемещений с единственной ненулевой компонентой uϕ ; K1 (t, τ ), ω(t, τ ), G(t) - ядро ползучести, мера ползучести и модуль упругой деформации при чистом сдвиге. Положим T◦ = (I - L(τ0 , t))TG-1 (2) и подействуем на выражение (1), содержащее T и его компоненты, оператором (I - L(τ0 , t)). Тогда, учитывая (2), получим следующую краевую задачу: ◦ ◦ ◦ ∂σϕz 2σrϕ ∂σrϕ + + = 0 (∇ · T◦ = 0), ∂r ∂z r z = 0, 0 r = b0 , 0 r z a : uϕ = α(t)r, z = 0, a ◦ l : σrϕ = 0, 1 E = [∇u + (∇u) ], 2 z = l, 0 (3) r r ◦ b0 : σϕz = 0, b0 : uϕ = 0; T◦ = 2E. Исходя из (3) можно получить, что перемещение uϕ удовлетворяет уравнению ∂ 2 uϕ ∂ 2 uϕ 1 ∂uϕ uϕ Vuϕ ≡ + + - 2 = 0. (4) ∂r2 ∂z 2 r ∂r r Следуя [11], возьмем решение (4) в форме (см. [6]) uϕ (r, z, t) = ld0 (t)r z 1- + b0 l ∞ n=1 dn (t) sh[δn (l - z)] J1 (rδn ) , δn sh(δn l) (5) где dk (t) (k = 0, 1, 2, . . . ) - известные функции времени, δn - неопределенные константы, а Jν (x) - функции Бесселя порядка ν. 687 М а н ж и р о в А. В., М и х и н М. Н., М у р а ш к и н Е. В. Заметим, что выражение (5) для перемещения uϕ удовлетворяет граничному условию в (3) на зажатой торцевой поверхности вала для z = l и позволяет выписать компоненты тензора преобразованных напряжений T◦ в виде (см. (1) и (3)) ◦ σϕz (r, z, t) = - d0 (t)r - b0 ∞ ◦ σrϕ (r, z, t) =- n=1 ∞ dn (t)J1 (rδn ) n=1 ch[δn (l-z)] , sh(δn l) (6) sh[δn (l - z)] . dn (t)J2 (rδn ) sh(δn l) Используя граничное условие из системы (3) на боковой поверхности вала (r = b0 ) и условия (6), найдем неизвестные постоянные δn . Действительно, ◦ к нулю при r = b , получим, что δ = λ /b , приравнивая выражение для σrϕ 0 n n 0 где λn - корни уравнения J2 (λn ) = 0. Наконец, удовлетворяя граничным условиям при z = 0, для поиска последовательности функций dk (t) получим следующие парные сумматорные уравнения: ld0 (t)r + b0 d0 (t)r + b0 ∞ n=1 ∞ b0 dn (t) λn r J1 λn b0 = α(t)r, 0 r a, (7) dn (t)J1 n=1 λn l λn r cth = 0. a b0 b0 r b, Поскольку λn λ1 ≈ 3.8317 и l/b0 = κ0 1, выражение для cth(λn l/b0 ) может быть задано равным единице с высокой степенью точности, и уравнения (7) можно исследовать в виде (см. [11]) ∞ uϕ (r, 0, t) = κ0 d0 (t)r + d0 (t)r ◦ + σϕz (r, 0, t) = b0 n=1 ∞ λn r b0 dn (t) J1 λn b0 dn (t)J1 n=1 = α(t)r, 0 r a, (8) λn r b0 = 0, a r b0 . Ряды в уравнениях (8) описывают сформулированную контактную задачу в интервале t ∈ [τ0 , τ1 ], в который время входит параметрически. Построим теперь решение (8), получив сначала разрешающие уравнения задачи при непрерывном росте и после прекращения роста. Отметим лишь, что истинные напряжения можно восстановить по формуле t T(r, z, t) = G(t) T◦ (r, z, t) + T◦ (r, z, τ )R1 (t, τ ) dτ , (9) τ0 где R1 (t, τ ) - резольвента ядра K1 (t, τ ). Пусть t ∈ [τ1 , τ2 ]. Тогда краевая задача для растущего вала, скрученного диском, примет вид (см. [7-10]) ∂srϕ ∂sϕz 2srϕ + + = 0 (∇ · S = 0), ∂r ∂z r 688 (10) Кручение растущего вала z = 0, 0 r a : vϕ = α(t)r, ˙ r = b(t), 0 z l : srϕ = 0, t = τ ∗ (r), 1 D = [∇v + (∇v) ], 2 z = 0, a z = l, 0 r r b(t) : sϕz = 0, b(t) : vϕ = 0; ∂T◦ = S. ∂t S = 2D, Видно, что скорость перемещения vϕ удовлетворяет уравнению Vvϕ = 0 (см. (4)), в то время как выражение для vϕ и скорость преобразованных напряжений srϕ и sϕz могут быть записаны в виде vϕ (r, z, t) = ld◦0 (t)r z 1- + b(t) l sϕz (r, z, t) = - d◦0 (t)r b(t) ∞ n=1 ∞ d◦n J1 [rηn (t)] - n=1 ∞ d◦n J2 [rηn (t)] srϕ (r, z, t) = - d◦n (t) sh[ηn (t)(l - z)] J1 [rηn (t)] , ηn (t) sh[ηn (t)l] n=1 ch[ηn (t)(l - z)] , sh[ηn (t)l] (11) sh[ηn (t)(l - z)] . sh[ηn (t)l] Здесь d◦k (t) и ηn (t) - функции, подлежащие определению. Удовлетворяя граничным условиям (10) и учитывая, что l/b(t) 1, мы приходим к уравнениям с бесконечными суммами для нахождения d◦k (t): ∞ vϕ (r, 0, t) = κ(t)d◦0 (t)r + n=1 sϕz (r, 0, t) = d◦0 (t)r + b(t) λn ηn (t) = , b(t) κ(t) = b(t)d◦n (t) λn r = α(t)r, ˙ J1 λn b(t) ∞ d◦n J1 n=1 l , b(t) τ1 rλn = 0, b(t) t a r 0 r b(t), a, (12) τ2 . Если d◦k (t) найдены, что означает, что и S и v также известны, то тензор напряжений T и вектор перемещений u можно восстановить по формулам T(r, z, t) = G(t) T(r, z, τ0 (r)) 1+ G(τ0 (r)) t + t R1 (t, τ ) dτ + τ0 (r) τ S(r, z, τ ) + τ0 (r) S(r, z, ζ) dζR1 (t, τ ) dτ , (13) τ0 (r) t u(r, z, t) = u(r, z, τ0 (r)) + v(r, z, τ ) dτ. τ0 (r) Краевая задача для растущего вала, когда прекращается его рост t τ2 = τ ∗ (b1 ) (b(t) = b1 ) и выполняются краевые условия τrϕ = 0, заданные на поверхности вала, имеет вид (10). Как и раньше, ее можно свести к краевой задаче в скоростях перемещения и преобразованных напряжений с решением в виде (11) при условии b(t) = b1 . Результирующие сумматорные уравнения сохраняют вид (12), где b(t) = b1 , κ(t) = κ1 = l/b1 , ηn (t) = ηn = λn /b1 , 689 М а н ж и р о в А. В., М и х и н М. Н., М у р а ш к и н Е. В. t τ2 . После их решения напряжения σrϕ , σϕz и перемещения uϕ определяются согласно (13). Следует отметить, что зависимость S и v от времени t параметрическая. К полученным уравнениям должно быть быть добавлено условие равновесия диска a σ(ρ, t)ρ2 dρ, M (t) = -2π σ(ρ, t) = σϕz (ρ, 0, t), (14) 0 которое выполняется во всем рассматриваемом временном промежутке. На основании условий (14) могут быть также получены следующие соотношения: a M (t) M ◦ (t) = (I - L(τ0 , t)) = -2π σ ◦ (ρ, t)ρ2 dρ, τ0 t τ1 , (15) G(t) 0 t ∂M ◦ (t) 1 ∂M (t) ∂M (τ ) ∂ω(t, τ ) = + dτ + ∂t G(t) ∂t ∂τ ∂t τ0 a ∂ω(t, τ0 ) + M (τ0 ) = -2π s(ρ, t)ρ2 dρ, t τ1 , (16) ∂t 0 которые в ряде случаев более удобны для построения решения контактной задачи. 2. Решение задачи кручения. Разрешающие парные сумматорные уравнения краевой задачи в трех основных временных интервалах могут быть представлены общими соотношениями ∞ ζϕ0 x + n=1 ϕn J1 (λn x) = ψx, λn 0 x c, (17) ∞ - p(x) = ϕ0 x + ϕn J1 (λn x) = 0, c x 1. n=1 Здесь для t ∈ [τ0 , τ1 ] мы положили ζ = κ0 , ϕk = dk (t), p(x) = σ ◦ (xb0 , t), ψ = α(t), c = a/b0 , x = r/b0 ; для t ∈ [τ1 , τ2 ] мы имеем ζ = κ(t), ϕk = d◦k (t), ψ = α(t), ˙ p(x) = s(xb(t), t), c = a/b(t), x = r/b(t), и для t τ2 в отличие от предыдущего случая - ζ = κ1 , b(t) = b1 . Построим решение (17), следуя описанному в [12] методу. Пусть p(x) = c ∂ ∂x x g(ξ) ξ 2 - x2 dξ h(c - x), (18) где h(c - x) - функция Хевисайда. Ряд во втором уравнении (17) является разложением Дини [13] функции (-p(x)), коэффициенты ϕk которого заданы формулами 1 c x2 p(x) dx = 8 ϕ0 = -4 0 2 2 J1 (λn ) 2 = 2 J1 (λn ) 1 ϕn = - 690 ξg(ξ) dξ, 0 xp(x)J1 (λn x) dx = 0 c g(ξ) sin(λn ξ) dξ, 0 n = 1, 2, . . . , (19) Кручение растущего вала в случае, когда учтено условие (18). Подставляя соотношения (19) в первое уравнение (17) и используя технику описанную в работах [12, 14-16], получим интегральное уравнение Фредгольма второго рода для определения функции g(x): c 3ψx , 1 x c, (20) π 0 ∞ 16 4 K2 (y) k(x, ξ) = (1 - 2ζ)xξ + 2 8xξI2 (y) - sh(xy) sh(ξy) dy, π π 0 I2 (y) g(x) + g(ξ)k(x, ξ) dξ = где Kν (y), Iν (y) - функции Бесселя мнимого аргумента порядка ν. Решение (17), очевидно, также дает полное решение рассматриваемой контактной задачи. Его можно найти с помощью методов, рассмотренных в работах [17, 18]. Здесь же для построения приближенного решения (20) воспользуемся методом, предложенным в [11]. Заметим, что для ζ 10 отклонение приближенного решения от численного решения не превышает 8.5 % при c = 0.7, 7 % при c = 0.6 и 1 % при c 0.5. Воспользуемся тем фактом, что величина ζ достаточно велика, и ограничимся первым членом в выражении для ядра k(x, ξ) (см. (20)): g(x) + 16 (1 - 2ζ)x π c g(ξ)ξ dξ = 0 4ψx π (1 x c). (21) Тогда, подставляя g(x) = Ax в уравнения (21) и определяя константу A, в силу (18) будем иметь p(x) = - 4ψ x √ π + 16(2ζ - 1)c3 /3 c2 - x2 (1 x c). Зависимости преобразованных контактных напряжений σ ◦ (r, t) и их скоростей s◦ (r, t) от угла поворота диска α(t) могут быть вычислены по формулам σ ◦ (r, t) = α(t)W (r, b0 ), τ0 t τ1 , ◦ s (r, t) = α(t)W ˙ (r, b(t)), τ1 t τ2 , ◦ s (r, t) = α(t)W ˙ (r, b1 ), t τ2 , r 4 √ . W (r, ξ) = - 3 3 π + 16(2l/ξ - 1)a /(3ξ ) c2 - r2 (22) (23) (24) Для заданного угла поворота диска контактные напряжения σ ◦ (r, t) и их скорости s◦ (r, t) тут же могут быть найдены из соотношений (22)-(24), а контактные напряжения σ(r, t) могут быть восстановлены по описанным ранее зависимостям. Момент, действующий на диск, вычисляется согласно (14). Заметим, что при α(t) = const взаимное влияние начального вала и его вновь образующейся ненапряженной части не проявляется. На основании соотношений (9), (15) и (22) для заданного крутящего момента M (t) получим r M (t) 3M (t) √ , τ0 , α(t) = B(b0 )(I - L(τ0 , t)) 3 2 2 4πa G(t) a -r 3 2l - ξ B(ξ) = + . 3 16a πξ 4 σ(r, t) = t τ1 , (25) 691 М а н ж и р о в А. В., М и х и н М. Н., М у р а ш к и н Е. В. Воспользовавшись зависимостями (16), (19), (23) и (24), мы наконец получаем соотношение (25) для определения контактных напряжений в промежутке времени t τ1 и следующие выражения для угла поворота диска: α(t) ˙ = α(t) ˙ = ∂M ◦ (t) B(b(t)), ∂t ∂M ◦ (t) B(b1 ), ∂t t α(τ ˙ ) dτ, α(t) = α(τ1 ) + τ1 t τ2 , τ1 t α(t) = α(τ2 ) + α(τ ˙ ) dτ, t τ2 . τ1 Оказывается, что рост вала при кручении диска моментом сил оказывает незначительное влияние на распределение контактных напряжений, если радиусы диска и вала не очень велики (конкретные отношения приведены выше). Однако в этом же случае появляется существенная зависимость угла поворота диска от момента и от скорости роста, когда вал начинает расти. Для расчета параметров напряженно-деформированного состояния вала на некотором расстоянии от его краевых поверхностей (в условиях принципа Сен-Венана) можно воспользоваться результатами полученными в работе [19]. 3. Численные расчеты. Рассмотрим контактную краевую задачу для вала, изготовленного из бетона с модулем упруго-мгновенной деформации сдвига G(t) = G = const и мерой ползучести при сдвиге в виде [1] ω(t, τ ) = (D0 + F e-βτ )(1 - e-γ(t-τ ) ). Произведем замену переменных согласно формулам r , a r∗ = ρ , a t∗ = M ∗ (t∗ ) = M (t) , Ga3 ρ∗ = t , τ0 σ ∗ (r∗ , t∗ ) = α∗ (t∗ ) = α(t), σ(r, t) , G β ∗ = βτ0 , τ1∗ = τ1 , τ0 τ2∗ = τ2 , τ0 γ ∗ = γτ0 , b0 b1 b(t) l , b∗1 = , b∗ (t∗ ) = , l∗ = , D0∗ = D0 G, F ∗ = F G a a a a и, опуская звездочку в обозначениях, зададим следующие значения функций и параметров: b∗0 = 1 290 b0 (t + τ2 - 2τ1 ) , l= , b(t) = , b1 = 2b0 , M (t) = 1, 0.7 0.7 τ2 - τ1 D0 = 0.251, F = 1.818, β = 0.31, γ = 0.6, τ0 = 10 дней. b0 = Видно, что во время роста вала его радиус удваивается. Скорость роста постоянна и определяется только временем начала и прекращения роста. Крутящий момент, действующий на диск, не изменяется со временем. Кроме того, отношение длины вала к его радиусу не менее 20 в течение всего процесса роста, в то время как отношение радиуса диска и вала не превосходит 0.7, т. е. можно воспользоваться формулами для поиска приближенного решения. Что касается распределения контактных напряжений, то достаточно вспомнить соотношения (25), чтобы отметить, что они (распределения) практически не зависят от свойств материала и процесса роста в рассматриваемом случае. 692 Кручение растущего вала Поведение угла вращения диска в зависимости от контактных напряжений требует более детального анализа. Кривые на рис. 2 изображают изменение угла поворота γ = 2α в момент времени t во время первого процесса наращивания. В этом случае рост вала начинается одновременно с наличием крутящего момента (τ1 = 1) для ˙ ˙ различных скоростей роста b(t): b(t) = b0 /9 (τ2 = 10) - сплошная линия ˙b(t) = b0 /3 (τ2 = 4) - штрих-пунктирная линия, и b(t) ˙ = b0 (τ2 = 2) - пунктирная линия. Моменты остановки процесса роста отмечены вертикальными сплошными линиями. Кривые на рис. 3 соответствуют зависимостям угла вращения диска от времени для случая повторного процесса наращивания. В этом случае для вала, загруженного в момент времени t = 1, и момента начала роста - τ1 = 2, ˙ ˙ нами выбраны различные скорости роста b(t): b(t) = b0 /8 (τ2 = 10) - штрих˙ пунктирная линия, b(t) = b0 /2 (τ2 = 4) - пунктирная линия. Для сравнения, изменение угла поворота диска, скручивающего вал фиксированного радиуса b0 , показано сплошной линией. Дуги кривых, расположенные между вертикальными сплошными линиями, характеризуют поведение угла вращения диска в промежутках непрерывного роста валов. Выводы. 1. Контактные напряжения, действующие на вал из-за взаимодействия с жесткими дисками, практически не зависят от свойств материала и процесса роста, в случае, когда вал достаточно длинный (отношение длины вала к его радиусу больше или равно 20) и в то же время диски не слишком велики (отношение радиуса диска и вала не превышает 0.7) в течение всего процесса наращивания. 2. Угол поворота γ(t) существенно зависит от скорости роста. Таким образом, предельное значение приращения угла поворота диска ∆(∞) (∆(t) = γ(t) - γ(τ0 )) во время медленного роста валов может превышать такое же значение для быстрого роста в 5 раз и более. 3. Предельное значение угла поворота диска существенно зависит от временного интервала между началом загрузки и началом роста. 4. Для постоянного крутящего момента существует характерное время, начиная с которого влиянием процесса наращивания на напряженнодеформированное состояние вязкоупругого вала можно пренебречь. 693 М а н ж и р о в А. В., М и х и н М. Н., М у р а ш к и н Е. В. γ 14 10 6 1 2 4 6 8 10 12 t Рис. 2. Угол поворота диска в зависимости от различных скоростей роста в случае первого процесса наращивания [Figure 2. The rotation angle of the disk depending on the different growth rates in the case of the first growth process] γ 20 16 12 8 1 2 4 6 8 10 12 t Рис. 3. Угол поворота диска в зависимости от различных скоростей роста в случае процесса повторного наращивания [Figure 3. The rotation angle of the disk depending on the different growth rates in the case of the second growth process] Конкурирующие интересы. Заявляем, что в отношении авторства и публикации этой статьи конфликта интересов не имеем. Авторский вклад и ответственность. Все авторы принимали участие в разработке концепции статьи и в написании рукописи. Авторы несут полную ответственность за предоставление окончательной рукописи в печать. Окончательная версия рукописи была одобрена всеми авторами. Финансирование. Работа выполнена при частичной финансовой поддержке Российского фонда фундаментальных исследований (проект № 17-01-00712_a).

About the authors

Alexander V Manzhirov

A. Ishlinsky Institite for Problems in Mechanics, Russian Academy of Sciences; National Engineering Physics Institute “MEPhI”; N. E. Bauman Moscow State Technical University

Email: manzh@inbox.ru
http://orcid.org/0000-0002-7578-6031 Dr. Phys. & Math. Sci., Professor; Deputy Director; Lab. of Modelling in Mechanics of Solids; Professor; Dept. of High Mathematics; Professor; Dept. of Applied Mathematics 101, pr. Vernadskogo, Moscow, 119526, Russian Federation; 31, Kashirskoe shosse, Moscow, 115409 Russian Federation; 5/1, 2-ya Baumanskaya st., Moscow, 105005, Russian Federation

Mikhail N Mikhin

Russian State University for the Humanities

Email: mmikhin@inbox.ru
Cand. Phys. & Math. Sci., Associate Professor; Associate Professor; Dept. of High Mathematics 9, Miusskaya pl., Moscow, 125993, Russian Federation

Evgenii V Murashkin

A. Ishlinsky Institite for Problems in Mechanics, Russian Academy of Sciences; National Engineering Physics Institute “MEPhI”; N. E. Bauman Moscow State Technical University

Email: murashkin@ipmnet.ru
http://orcid.org/0000-0002-3267-4742 Cand. Phys. & Math. Sci.; Senior Researcher; Lab. of Modelling in Mechanics of Solids; Associate Professor; Dept. of High Mathematics; Associate Professor; Dept. of Applied Mathematics 101, pr. Vernadskogo, Moscow, 119526, Russian Federation; 31, Kashirskoe shosse, Moscow, 115409 Russian Federation; 5/1, 2-ya Baumanskaya st., Moscow, 105005, Russian Federation

References

Supplementary files