ВОЗДЕЙСТВИЕ НЕСТАЦИОНАРНОЙ РАСПРЕДЕЛЕННОЙ НАГРУЗКИ НА ПОВЕРХНОСТЬ УПРУГОГО СЛОЯ Текст научной статьи по специальности «Математика»

Электронный журнал «Труды МАИ». Выпуск № 71

www.mai.ru/science/trudy/

УДК 539.3

Воздействие нестационарной распределенной нагрузки на

поверхность упругого слоя

Кузнецова Е.Л.*, Тарлаковский Д.В.**, Федотенков Г.В.,***Медведский А.Л.

Московский Авиационный Институт (национальный исследовательский университет), МАИ, Волоколамское шоссе, 4, Москва, А-80, ГСП-3, 125993, Россия

*e-mail: vida_ku@mail.ru **e-mail: tdvhome@mail.ru ***e-mail: greghome@mail.ru

Аннотация

Рассматривается плоская нестационарная задача о воздействии произвольно

распределенной поверхностной нагрузки f (S, t) = [ f (x, т) f3 (x, x)J на

тонкослойную обшивку летательного аппарата, моделируемую линейно упругим однородным изотропным слоем постоянной толщины. Движение слоя рассматривается в прямоугольной декартовой системе координат Oxyz. Ось Oz направлена вглубь упругого слоя, а Ox - вдоль его свободной поверхности z = 0. Предполагается, что внешняя нагрузка не зависит от координаты y, что приводит к плоской постановке задачи. Для решения используется принцип суперпозиции фундаментальных решений, согласно которому перемещения в упругом слое представляются так [1, 2]

т b(t)

u(x,т,z) = J J G(xт-1,z)f(S,t)d^dt,

0 a (t)

Здесь G (x, т, z) - матрица фундаментальных решений (функций влияния).

Функции влияния являются решением начально-краевых задач о воздействии на границу упругого слоя сосредоточенной нагрузки распределенной по закону 5(х)5(т), где 5 - дельта-функция Дирака. Для построения функции влияния

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

о /о о о

д + я и ^ +п s на экспоненту с показателем - линейной комбинации этих корней (д, я - параметры преобразований Фурье и Лапласа соответственно; п -безразмерный параметр, характеризующий свойства материала). Это позволяет использовать метод совместного обращения преобразований Лапласа и Фурье, основанный на построении аналитических представлений изображений [4]. В случае наличия в степени экспоненты одного корня этот алгоритм позволяет получить оригиналы в явном виде. Если степень экспоненты содержит более одного корня, то для получения оригиналов используется модифицированный алгоритм совместного обращения изображений Фурье-Лапласа [3, 5].

Для получения решения разработан и реализован численно-аналитический алгоритм, основанный на методе квадратур [6] и формулах Симпсона. Приведены примеры расчетов.

Ключевые слова: тонкослойная обшивка летательных аппаратов, нестационарные воздействия, упругие волны, интегральные преобразования,

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

Введение

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

Постановка задачи

В начальный момент времени к поверхности г = 0 предварительно невозмущенного линейно упругого слоя, моделирующего наружную обшивку ЛА, толщины к прикладывается нестационарная распределенная нагрузка. Движение слоя рассматривается в прямоугольной декартовой системе координат Охуг. Ось Ог направлена вглубь упругого слоя, а Ох - вдоль его свободной поверхности г = 0 (рис. 1). Предполагается, что внешняя нагрузка не зависит от координаты у , что приводит к плоской постановке задачи.

^ (х, т) = / (х, т) Я (т) [ Н (х - а (т))-Н (Ь (т)-х)]

(х, т) = /з (х, т) Н (т) [ Н (х - а (т))-Н (Ь (т)-х)]

х

ХХХХКХХХХХХХ] а(т)

кххххххххххххх

Ь(т)

Рис.1.

Границы носителя внешней нагрузки по пространственной переменной в общем случае изменяются во времени по заданным законам: а (т) - закон изменения левой границы, Ь (т) - закон изменения правой границы (рис. 1).

В работе используем следующую систему безразмерных величин (при одинаковом начертании они обозначены штрихом, который в дальнейшем опускается):

, х , г сл , ф , ш ,, к , и , w , х = —, г = —, т = -Ц ф =-пГ, ш = -¡2, к = —, и = —,= —,а. . = -—1—, Ь Ь Ь Ь2 Ь Ь Ь Ь 1 Х + 2|

/;=/ (к=1,3), с2, с2=|,п=к:

X

Х + 2|

р

р

Х + 2|

2 п2

где t - размерное время, Ь - некоторый линейный размер, ф и ш - скалярный и ненулевая компонента векторного потенциала упругих смещений и вдоль оси Ох и

w вдоль оси Ог; и /к

компоненты тензора напряжений и вектора

поверхностной нагрузки, с1 и с2 - скорости распространения волн растяжения-

сжатия и сдвига; X, | и р - параметры Ламе и плотность среды.

г

2

Уравнения движения среды в потенциалах упругих смещений и с учетом плоской постановки задачи имеют вид [1]

ф = Дф, у = п2Ду, А

д2 д2

(1)

дх2 дг2'

Связь перемещений и напряжений (указаны лишь напряжения, которые могут входить в граничные условия) с потенциалами определяется соотношениями Коши и законом Гука

дф ду дф ду

дх дг ' дг дх'

а33 =

дw ди + К

1

дг

дх' °13 п2

ди дw — + —

дг дх

(2)

В начальный момент времени среда находится в невозмущенном состоянии, что соответствует однородным начальным условиям

ф|т=0 =ф1 т=0 = 4=0 = У1 т=0 = 0 (3)

Для определенности предположим, что нижняя граница упругого слоя жестко закреплена, тогда приходим к следующим граничным условиям

а

33 1г=0

= -д (х, т) , а1з |г=0 =-Р (х, Т) , и|г=, = Нг=А = 0

г=к \г=к

Метод решения

В основу метода решения положим принцип суперпозиции, согласно которому перемещения в упругом слое представляются так [1, 2]

т Ь(г)

и (х, т, г ) = 11 О (хт — г, г)£(£, Ы^г,

0 а (г)

(4)

и = Г и ^ , О = ^13 Л , £ = ' /1 ^

V ^31 ^33 у V /3 У

Здесь Оц - функции влияния упругого слоя, которые представляют собой перемещения как решения задачи (1) - (3) с граничными условиями (задача 1)

а331г=0 = 0, а13 1г=0 = 5(х)5(т) , Нг=к = Ч=к = 0,

или (задача 2)

&331г=0 =§(х)§(т), а13 Iг=0 = 0, Нг=к = Ч=к = 0. Здесь 5(^) - дельта-функция Дирака.

В обозначениях функций влияния упругого слоя первый индекс равный 1 соответствует перемещению по направлению оси Ох , равный 3 - по направлению оси Ог; второй индекс равный 1 соответствует задаче 1, а равный 3 - задаче 2.

Для построения функции влияния, как решения задач 1-2, используются интегральные преобразования Фурье по координате х и Лапласа по времени т. В пространстве преобразований Фурье-Лапласа приходим к следующим краевым задачам [3]

А2,г\ ЬГ Л 2,, ГЬР

а ф /2/2 2\,ЛЬГ п а ¥ /2/2 2\„ ЬР А

-—2--к1 (Я , * )ф = а—-2--к2 (Я , ^ )ш = 0 _

аг аг (5)

к1 (я, я) = я + я, к2 (я, я ) = у1 я + п2я;

Г^ЬГ • „ЬГ аш ^ЬГ аф •

О1к = -Яф--1—, °3к =—,--1ЯШ ,

аг аг (6)

ьг л ( г^ЬГ \ (6)

аО,,

„ЬГ = О- - цкОЬГ, < =-12

аг п

1к ;„гЬр --1ЯО3к

V аг у

ЬГ „33

к3 13

г=0

31

г=к

5к3, „13 0 =5к1, О11 к = О31 к =0. (7)

г=0

г=к

Здесь 5г1 - символ Кронекера.

Функция влияния О33 найдена и исследована в работе [3]. Здесь ограничимся описанием методики построения функции влияния О13 .

Решения уравнений (5) имеют вид (А1, В1, А2 и В2 - постоянные интегрирования):

фьр (д, s, г) = Де^)+ Ве^^)

у1Р (д, г) = А/*2^)+ В2е гк2^)

Подставляя их в равенства (6), получаем следующие выражения для функций

влияния и напряжений:

^ (д, г) = —гд

— к2 (д2, *2)

С33 (д,г) = к1 (д2, я2)

А/1^2)+ Vгк1( ^)'

Аге

к (д2,^2)— В — гк2 (д2,*2)' В2е

А/1^2)— В1е"^)

гд

\е* ^)+ В2е" гк2(^)'

П2а3^ (д,г ) = к3 (д2,*2)

— 2гдк2 (д2, я2) Л2а^ (д,г ) = —2гдк1 (д2, я2)

г^д2'*2К « д2,я2)'

Де ' ; + В1е

, гк2 (д2,*2) —гк2(д2,.2)

Де — В2е

А1егк1(д2,*2)— В1е" )'

к3 (д \ *2)

Де^2 (д2У)+ В2е" гк2(д2у)'

к3 (д, *) = к22 (д, *) + д = 2д + п2 *

(8)

Используя теперь граничные условия (7), получаем систему линейных алгебраических уравнений относительно постоянных интегрирования, решение которой записываем так:

П2 Р„ (д2, *2; §,, §2) 1 = й (д2, §„ ^2)

, А2 = 2гдк1 (д2,*2) 2

й (д2, *2; §1, §2)

В1 =

Л Р (д , * ; §р §2) „ „,. , / 2 х3

2 ^ с /, В2 = 2гдк1 (д2,*) „Л 2

п2 Р3 (д2, *2; §1, §2)

(9)

й (д2, *2; §1, §2)

й (д2, *2; §1, §2)

где

-кк1(д2,$'2) ^ - Нк2 (я2,я2)

$1 = е , $2 = е ,

В(Я, я;$1,^2 ) = (1 + $2$2) Do (Я, я) + ($? + $2) В, (Я, я) + $&В2 (я, я),

В0 (Я, я) = -^1- (Я, я) ^2- (Я, я), В1 (Я, я) = ^1+ (Я, я) ^2+ (Я, я), В2 ( я, я) = -16як1 ( я, я) к2 ( я, я) к3 ( я, я),

(Я, я) = Я ± к1 ( Я, я) к2 ( Я, я), ^2±( Я, я) = к32 ( Я, я) ± 4як1 ( Я, я) к2 ( Я, я), Р (Я, я; $1, $2 ) = Р (Я, я )$2 + Р (Я, я )$& + Рм (Я, я )$2$2, Р (Я, я; $1, $2 ) = Р (Я, я) + Р (Я, я )$2 + Р (Я, я )$1$2, Р, (Я, я; $1, $2 ) = Р-2 ( Я, я )$2 + Р23 ( Я, я )$1$2 + Р24 ( Я, я )$12$2,

Р (Я, я; $2 ) = р30 (Я, я) + р31 (Я, я )$2 + р33 (Я, я )$l$2, р01 (Я, я ) = р12 (Я, я ) = к3 (Я, я) Р+( Я, я),

р03 (Я, я ) = Р3 (Я, я) = -4Як1 (Я, я) к2 (Я, я),

Р04 (Я, я ) = р10 (Я, я ) = -к3 (Я, я) Р-( Я, я), р22 (Я, я )=-Р31 (Я, я ) = Р+( Я, я), р33 (Я,я) = -р23(Я,я) = к3(Я,я), р24(Я,я) = -р30(Я,я) = Р-(Я,я).

Подставляя равенства (9) в (8) определяем трансформанту искомой функции

влияния О13

ЬГ

О.Т(Я,я,г) = - В( 2% $ )I"Р(Я2,я2;$1,$2)^) +

В ( Я , я ; $1, $2 ) -

+р(Я2,$1,$2)4^) + 2*1 (Я2,я2)к2(Я2,^)Р(Я2,я2;$1,$2)е*() - (10)

2к1 (Я2, *2) к2 (Я2, <2) Р (Я2, <2; $1, $2) е" * Я')

2 2 \ 1 I 2 2 \ п ( 2 2 й й \ -гк2(Я2,я2

Вычислить оригинал функции влияния О13 непосредственно по его изображению (10) не представляется возможным. Поэтому, учитывая, что в некоторой полуплоскости Яе я >а0 имеют место неравенства < 1 и |$2| < 1, используем следующее разложение в степенной ряд:

й0 (g, *)

й (д, §1, §2) 1 + (§? + §2) К (д, * ) + §1§2К2 (д, * ) + §2§2

= I (—1)" Г(§2+§2) К (д. * )+§,§2*-2 (д, * )+§2§2 Т = (")

п=0

= I(—1)" I (п;"2, "3, "4 )К""3. (д,*)§2(п'+п4>+"3§2("+"4>+"3

2' 3' 4) n,n2n3\^, )

+"3 +"4 ="

ГЦГ^

"=0

где

К, (д, * ) = МЦ, К2 (д, * ) = •

й0 (g, *) й0 (g, *)

КЩП2Щ(g,*) = КГ+"2(^*)К"3(^*), ("2,nз,"3) =

"1 !"2 !"3 !"4!

Подставляя ряд (11) в формулу (10), получаем следующее представление для изображения функции влияния:

(g, ^ г ) = I ^ (g, ^ г ) (g, ^ г ) = g, ^ г) + (g, ^ г),

"=0

<^1^3^"(g,^г) = (—1)"(g,^г) I ("рn2,nз,"4)Kn1n2n3(g2,)х (12)

^+"2 + +"3 +"4 ="

—[2(nJ +"4)+"3]Лк1 (д2,*2 ) —[2("2 +"4)+"]й*2(д22)

х е е

Здесь

(g, ^ г) = <3о,1 (g, ^ г) + <£,2 (g, ^ г),

1 ( ) гдк3 (д2,*2 ) „ ( 2 2 ) — к:(д2,*2)г <130,1 (g,^г)= 2 / 2 ^^2+(д ,* )е ,

* 2б3 (д2, *2)

^ (д,,г)=)(13)

* 63 (д , * )

п2 ¿а (д, * ) = к2—( д, *) к2+( д, *),

б3 (д, *) = 8 (1 + к) д3 + 8п2 (2 + к) д2 * + 8п V2 + Л6

1

оо

(ь, г Ь-г-^Ч I ЗД2, * V ^ ")(2к-г) +

,* ) I-

+р03(д2, * 2)в~к1(ь 2Лк-г е^2*2» + р04(ь * 2)е-к(г^)(2н-г )е-к2(г^)2ь +

+р12(ь2,*2)е"к2( )2ке')г + р13(ь2,* Vк1(** )(к+г)е"к2(« ' )к +

22

+2к1 (ь ,^)к2(^,^)Р22(^,^+ (14)

+2к1 (ь2, * 2)к2(ь2, *2) Р23(ь2, * 2)е" ^^е" к-г) + +2к1 (ь2,*2)к2(ь2,*2)Р24(ь2,*г^-^^е"^2'*2*2' ' -2к1 (ь2, * 2)к2(ь2, *2) Р31(ь2, * 2)е" ^^е" г

2

22

Как видно из представления (12), каждый член ряда представляет собой сумму слагаемых в виде произведений однородных степени (-1) рациональных функций от

параметров преобразований и квадратных корней к1 (ь2,*2) и к2 (д2,*2) на

экспоненту с показателем - линейной комбинации этих корней. Это позволяет с помощью замены ь = X* представить каждый из членов ряда в виде:

/1Р (ь,г) = ^(*)Ио(\)е-*&(Х\ gL(*) = *-1 (15)

Здесь к0(Х) - рациональная функция от X и квадратных корней к1 (X2,1) и к2(X2,1), а

ш(Х) - линейная комбинация этих корней с неотрицательными коэффициентами.

Для таких функций используется метод совместного обращения преобразований Лапласа и Фурье, основанный на построении аналитических представлений изображений [4]. Оригинал функции (15) имеет вид [3,4]:

/ (x, г^ Z2, т) = - [Д+( х гР Z2, т)] Д +(х Z2, т)-

-ко [ Д- (х, г1, г2, т) ] Д - (х, г1, г2, т)} Н (т - ш), С = х + ¡у, (16)

Х(С, г1, г2, т-шо) = Д(£, г1, г2, т), Д±( х, г1, г2, т)= Иш д(£, г1, г2, т),

у^±0

где Я(т) - функция Хевисайда, а г1, 12,т) - функция, неявно задаваемая уравнением

В случае наличия в функции г1, г2) одного корня этот алгоритм позволяет получить оригиналы в явном виде. Если г1, г2) содержит более одного корня, то для получения оригиналов используется модифицированный алгоритм совместного обращения изображений Фурье-Лапласа [3, 5]. Он основан на сведении уравнения (17) к задаче Коши для обыкновенного дифференциального уравнения в комплексной плоскости р. Интегрирование проводится вдоль прямых с начальной и конечной точками р = у и р = х (рис. 2):

о>(р, г1,12) + ¿рС = т.

(17)

1т р

Яе р

Рис. 2.

После построения функций влияния упругого слоя процесс получения решения сводится к вычислению интегралов в представлении (4). Для этого удобно использовать численный алгоритм.

Алгоритм вычисления перемещений упругого слоя Для вычисления интегралов входящих в (4) используем численный алгоритм, основанный на методе Симпсона [6]. На пространственно-временную область Я наносим ортогональную сетку с равномерным шагом А (рис. 3): г, = /А, 5. = -А, Я\ = ииК ((= 0,1,2,..., - е г) К- = {(г,5) | г_1 <г<г,\-_1 1}

Функциям и (х,т, г), w(x, т, г), / (х,т) и /3 (х,т) ставим в соответствие сеточные функции итп (г ) = итп (хт, т„, г), мтп (г ) = мтп (хт, тй, г), = !х (хт, \),

/зшп = /3 (хт, тп ) ( хт = тА , \ = пА ).

Границу области дБ(т) = {а(г)и Ь(г) |г е [0,т]} приближенно заменяем кусочно-постоянными функциями а(г) и в(г), определяемыми по формулам а(г) = А[а(г^А], в(г) = А[Ь(г)/А], где [•] означает операцию округления до ближайшего целого числа. В этом случае в момент времени т = гп пространственно -временная область интегрирования Бп = {(г,£)|0 < г < гп, а (г )<5< Ь (г)} аппроксимируется выпуклым многоугольником Вп, состоящим из объединения

п 121 -1

элементарных вертикальных полос Вп = и и К, А1 = а/А, 121 = в ¡/А, а 1 = а(г 1),

/=1 1ц

Рг=Р( г).

Для интегралов в (4), используем квадратурные формулы, основанные на методе Симпсона

г Ь г) Л 2 ^ X

I | 01к(х-5,т-г,г)(5,ЫУЬ ££,

0 а (г) 36 ¡=1 }=1и

= 16^/к (хт -57 - 0.5Л, тя - г, - 0.5Л, г) 11 (5j + 0.5Л, г + 0.5Л) + +4[о*(хт -5j -0.5Л,тя -г,,г) 11 (5j + 0.5Л,гг) + +0,к (хт - 5; - 0.5Л, тя - г, - Л, г)/1 (5^ + 0.5Л, гг + Л) + +0/к (хт -5;, тя - г1 - 0.5Л, г) (5,, Ь + 0.5Л) + +01к (хт - 5; - Л, тя - г - 0.5Л, г)/1 (5^ + Л, гг + 0.5Л)] + +о/к(хт-5 j, тя- г, г) /1 (5 j, г) + +о/к(хт - 5 j - Л, тя- г, г) /1 (5 j + Л, г) + +о/к(хт - 5 j, тя- г - Л, г) /1 (5 j, ^ + Л) + +Ок(хт - 5 j - Л, тя- г - Л, г) /1 (5 j + Л, г + Л) (I, к = 1,3).

dD: Ç=b(t)

dD: Ç=a(t)

Рис. З.

Примеры расчетов

В качестве примера рассмотрим различные варианты задачи о воздействии на границу упругого слоя нестационарной нормальной нагрузки

-д(х,т) = Н [х - а(т)] - Н [Ь(т)- х].

Материал слоя обладает следующими безразмерными параметрами: к = 1, П = 1.7 (значение параметра п = 1.7 соответствует механическим характеристикам низкоуглеродистой стали). Значение безразмерного линейного параметра Ь принято равным 1.

На рис. 4 представлено распределение нормальных перемещений по координате х при постоянном носителе внешней нагрузки: -а(т) = Ь(т) = 1.

Распределения построены при значении г = 0.1. Видно, что в начальные моменты времени перемещения распределяются по близкому к равномерному закону, повторяя профиль нагрузки. Затем отчетливей начинают проявляться волновые эффекты.

На рис. 5 представлены распределения нормальных перемещений, соответствующие симметричному распределению нагрузки по координате х , границы носителя которой расширяются по законам: -а (т) = Ь (т) = т2 (г = 0.1). Из

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

Аналогичные зависимости представлены на рис. 6 при симметричном расширении границ носителя нагрузки по законам -а (т) = Ь (т) = л/т на глубине г = 0.1.

Распределения перемещений на рис. 7 построены при значении г = 0.2 и при несимметричном характере расширения границ носителя внешней нагрузки:

а(т) = -т2, Ь(т) = >/т. Как и следовало ожидать, это приводит к несимметричному

характеру распределений перемещений.

w

0.2

0.1

0

0.5

х

Рис. 4.

w

0.2

' т = 0.3 \

Л т = 0.5

0

0.5

х

Рис. 7.

Выводы

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

2. С помощью модифицированного алгоритма совместного обращения интегральных преобразований Фурье-Лапласа построены оригиналы функций влияния упругого изотропного слоя.

3. С использованием принципа суперпозиции фундаментальных решений

построены интегральные представления нестационарных решений задач о

воздействии на границу упругого слоя внешней нагрузки, распределенной по произвольному пространственно-временному закону.

4. Разработан и реализован на ЭВМ численный алгоритм расчета нестационарных процессов распространения упругих волн перемещений в слое под влиянием внешних нагрузок.

5. Приведен ряд расчетных примеров и проанализированы результаты.

Работа выполнена в Московском авиационном институте (национальном исследовательском университете) при финансовой поддержке работ по проекту Минобрнауки «Инновационный спускаемый с орбиты аппарат-демонстратор внедрения аэроупругих развертываемых при полете в космосе и в атмосфере элементов конструкций в космической технике», а также при поддержке грантов РФФИ № 12-08-31298-мол_а, № 14-08-01169-а, № 12-0100273, № 12-01-31220-мол-а, и гранта ФЦП (соглашение № 14.В37.21.0381 от 02.08.2012). Библиографический список

1.Горшков А.Г., Медведский А.Л., Рабинский Л.Н., Тарлаковский Д.В. Волны в сплошных средах: Учебное пособие для вузов. М.: Физматлит, 2004. 632 с.

2.Горшков А.Г., Тарлаковский Д.В., Федотенков Г.В. Плоская задача о вертикальном ударе цилиндрической оболочки по упругому полупространству // Известия Российской академии наук. Механика твердого тела. 2000. № 5. С. 151158.

3.Кузнецова Е.Л., Тарлаковский Д.В., Федотенков Г.В. Распространение нестационарных волн в упругом слое // Известия Российской академии наук. Механика твердого тела. 2011. № 5. С. 144-152.

4.Суворов Е.М., Тарлаковский Д.В., Федотенков Г.В. Плоская задача об ударе твердого тела по полупространству, моделируемому средой Коссера // Прикладная математика и механика. 2012. Т. 76. № 5. С. 850-859.

5.Вестяк В.А., Садков А.С., Тарлаковский Д.В. Распространение нестационарных объемных возмущений в упругой полуплоскости // Известия Российской академии наук. Механика твердого тела. 2011. № 2. С. 130-140.

6.Горшков А.Г., Амар Абдул Карим Салман, Тарлаковский Д.В., Федотенков Г.В. Удар деформируемым цилиндрическим телом по упругому полупространству // Известия Российской академии наук. Механика твердого тела. 2004. № 3. С. 82 - 90.