Моделирование пространственного напряженно-деформированного состояния балки-стенки Текст научной статьи по специальности «Физика»

СТРОИТЕЛЬНЫЕ КОНСТРУКЦИИ, ЗДАНИЯ И СООРУЖЕНИЯ

УДК 624.072/2:539.3

БАРАШКОВ ВЛАДИМИР НИКОЛАЕВИЧ, докт. физ.-мат. наук,

ст. науч. сотр.,

[email protected]

МАТВЕЕНКО АНАСТАСИЯ АЛЕКСАНДРОВНА, студентка,

Томский государственный архитектурно-строительный университет, 634003, г. Томск, пл. Соляная, 2

МОДЕЛИРОВАНИЕ ПРОСТРАНСТВЕННОГО НАПРЯЖЕННО-ДЕФОРМИРОВАННОГО СОСТОЯНИЯ БАЛКИ-СТЕНКИ

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

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

BARASHKOV, VLADIMIR NIKOLAYEVICH, Dr.of phys.-math. sc., [email protected]

MATVEENKO, ANASTASIYA ALEKSANDROVNA, student,

Tomsk State University of Architecture and Building,

2 Solyanaya sq., Tomsk, 634003, Russia

MODELLING OF SPATIAL STRESS-DEFORMED STATE OF THE BEAM-WALL

In three-dimensional presentation the calculation of the elastic stress-deformed state of a beam-wall is carried out. Geometrical ratio is taken in the form of Cauchy equations. Numerical realization of a problem is carried out by variation-difference method. The obtained system of the linear algebraic equations is realized by a method of the top relaxation with a choice

© В.Н. Барашков, А.А. Матвеенко, 2010

of optimum factor of a relaxation. The analysis of influence of the support structure on a kind of the stress-deformed state is performed.

Keywords: beam-wall, problem of the theory of elasticity, moving, deformation, stress, Lagrange principle, functional of full potential energy, variation-difference method, digitization, system of the linear algebraic equations, iterative methods of the solution.

Работа посвящена расчёту упругого пространственного напряжённо-деформированного состояния (НДС) балки-стенки (рис. 1), нагруженной ступенчато-равномерной внешней нагрузкой интенсивностью q и qf . Собственный вес конструкции не учитывается. Материал балки-стенки считается изотропным.

Рис. 1. Балка-стенка

Вариационно-разностный метод решения задачи теории упругости

Определение трёхмерного НДС осуществляется приближённым вариационно-разностным методом (ВРМ) [1, 2], реализующим с использованием конечно-разностных представлений вариационный принцип Лагранжа:

5Э = 0. (1)

В выражении (1) под знаком вариации находится функционал полной потенциальной энергии системы «тело-нагрузка», который определяется классическим соотношением

Э = и - А1 - А2,

где и - потенциальная энергия деформации тела; А1, А2 - работа объёмных и поверхностных сил на вызванных ими перемещениях.

Составляющие функционала энергии Э в декартовой системе координат х, у, 2 записываются следующим образом:

А =| (ХУи + У¥у + ZVw ) ёУ;

V

А2 = | (ХБи + У3у + ZSw)ёБ,

где сх, с , ..., угу, ухг — компоненты тензоров напряжений и деформаций соответственно; ХУ, УУ, ZV — компоненты объёмных нагрузок; ХБ, , ZS — компо-

ненты поверхностных нагрузок; и, V, w — компоненты вектора перемещения в направлении осей координат х, у, г соответственно; У — объём конструкции; Бс — часть поверхности конструкции Б, на которой заданы внешние нагрузки.

Принцип Лагранжа является экстремальным принципом, так как для действительного НДС функционал полной потенциальной энергии системы в состоянии устойчивого равновесия, как это следует из теоремы Лежена Дирихле, достигает минимума по сравнению со всеми возможными состояниями упругого тела. Поэтому задача определения НДС тела при заданных нагрузках состоит в нахождении функций перемещений и(х, у, г), v(х, у, г), w(х, у, г), доставляющих минимум функционалу потенциальной энергии Э.

Механически уравнение (1) в интегральной форме выражает условия равновесия деформированного тела. При реализации минимума потенциальной энергии дифференциальные уравнения равновесия внутри тела и уравнения равновесия на его поверхности (т. е. граничные условия в напряжениях) при вариационной постановке задачи теории упругости являются уравнениями Эйлера функционала энергии Э и выполняются автоматически. При решении задачи теории упругости в дифференциальной постановке удовлетворение уравнений Навье и условий на поверхности представляет значительную трудность.

В основе численных методов решения задач механики деформируемого твёрдого тела лежит замена континуальной расчётной модели с непрерывным распределением параметров и бесконечным числом степеней свободы дискретной моделью, имеющей конечное число неизвестных. Нанесём на исследуемое тело пространственную конечно-разностную сетку размером ( х у х к). Здесь /, у, к - количество узлов сетки по координате х, у, г соответственно. Дискретизация расчётной области балки-стенки проводится с помощью ограниченных четырёхугольными плоскостями многогранников (шестигранников), которые традиционно называются ячейками и из которых составляется так называемый шаблон (рис. 2).

Сеточный аналог функционала энергии реализуется заменой составляющих его интегралов конечными суммами, а производных - их конечноразностными представлениями. Приближённое конечно-разностное выражение сеточного аналога Э(и, vp, wq) получается в предположении, что все

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

k

Рис. 2. Шаблон для случая дискретизации расчётной области с помощью многогранников

Для аппроксимации частных производных от перемещений в трёхмерном пространстве в декартовой системе координат используются выражения производной по объёму от скалярной функции F, которые преобразуются с помощью формулы Остроградского - Гаусса в поверхностные интегралы по замкнутой поверхности:

(££> ¥п^ Е п.dS (££> Е п^

дЕ Ч -я 7

— = Ит—---, — = Ит—---------------------, ------= Ит—----------. (2)

дх V у ду ¥ V дz г V

Здесь пх, пу, пг - косинусы углов между вектором нормали п к замкнутой поверхности S, заключающей в себя объем V, и соответствующей осью координат. При этом функция Е должна быть непрерывной в объёме V и иметь там непрерывные ограниченные первые частные производные. В выражениях (2) под функцией Е следует понимать функции перемещений и, V, ^ .

Использование необходимого условия экстремума функции энергии Э(и, , Vp , )

дЭ(и ,V , ^ ) дЭ(и_,Vр, wя) дЭ(и,Vр, wя)

у р’ я' = о у р’ д' = о у р’ д' = о (3)

ди ’ дv ’ дw ’

* р я

(^ = *1, *1 +1,..., S; р = р1, р1 +1,..., Р; я = д1, д1 +1, ..., О) сводит задачу минимизации функции многих переменных Э(и*, Vр, wq) к решению системы N линейных алгебраических уравнений (СЛАУ) относительно N искомых компонент вектора перемещений и*, vp, wq узлов конечноразностной сетки:

N

Еатп/п = Ьт, т = 1,2 ..., N. (4)

п=1

Здесь коэффициенты атп - элементы матрицы {А} СЛАУ; /п - проекции и*, vp, wq вектора перемещений на оси координат; Ьт - свободные члены, включающие в себя статические и геометрические граничные условия; N = ^ - *1) + (Р - р1) + (О - я1) + 3.

Получающаяся матрица системы (4) имеет ленточную структуру, симметрична относительно главной диагонали и положительно определена, что очень важно при численной реализации. Длина строки матрицы или, что то же, ширина ленты зависит от размеров конечно-разностной сетки и порядка обхода узлов сетки при составлении уравнений (4) и может достигать большой величины. Таким образом, ширина ленты матрицы влияет на размер самой матрицы {А} в (4). Следует отметить, что в строке матрицы подавляющее большинство коэффициентов атп - нули.

Реализация СЛАУ (4) проводилась итерационным методом верхней релаксации с выбором оптимального коэффициента релаксации [3], и поэтому не было необходимости получать и использовать при решении громоздкую матрицу коэффициентов {А}.

К целому ряду представленных в [2] достоинств итерационных методов решения СЛАУ следует отнести тот факт, что итерационные методы «... в отличие от прямых методов имеют тенденцию быть самокорректирующимися и, следовательно, минимизируют ошибки округления» [4]. Таким образом, отдельная ошибка, допущенная в вычислениях, не отражается на окончательном результате, так как ошибочное приближение можно рассматривать как новый начальный вектор.

Метод верхней релаксации получается модификацией метода Зейделя (Гаусса - Зейделя), сходимость которого доказана на основе общей теории разностных схем [5], введением в его итерационную схему коэффициента релаксации ю , который служит для ускорения сходимости итерационного процесса решений СЛАУ. Алгоритм метода верхней релаксации записывается следующим образом:

т -1 N

а г(к+1) = (1 а ) - юЁ а /~('к+1) -шЁ а к) + юЬ .

тт^ т V ' тп^ т тпл п тпл п т

п =1 п=т+1

В работе [6] доказано, что для 0 < ш< 2 метод всегда сходится.

При получении аналитического выражения условия стационарности (3) функции энергии в узле (/, у, к) конечно-разностной сетки (т. е. трёх уравнений СЛАУ) с целью определения соотношений для вычисления коэффициентов атп при искомых перемещениях в уравнениях (4) необходимо пользоваться шаблоном, представленным на рис. 2. В общем случае, в каждом из трёх уравнений содержатся 81 слагаемое вида атп/п с перемещениями 27 узлов конечно-разностной сетки и свободный член. Таким образом, при варьирова-

нии энергии в узле (, у, к) необходимо вычислять 243 коэффициента и три свободных члена.

Структура выражения функции энергии такова, что от смещения в узле (, у, к) сетки с номером 0 (рис. 2) зависит лишь часть энергии Э(и, ур, ),

а именно — энергия ячеек, содержащих этот узел [7]. В зависимости от места положения узла, количество таких ячеек меняется от одной в угловом узле до восьми ячеек во внутреннем узле конечно-разностной сетки. Поэтому в общем случае выражения для производных от функции энергии (3) в узле (, у, к) также состоят из восьми слагаемых по количеству ячеек, которым

принадлежит рассматриваемый узел. По этой же причине не все уравнения СЛАУ (4) содержат 81 слагаемое.

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

и = ( А - А )/2, (5)

или из соотношения, получаемого подстановкой (5) в выражение для функционала энергии

Э = и - А1 - А2 = и - 2и = -и. (6)

Следует отметить, что выражение (6) позволяет в процессе отладки программы и счёта задачи контролировать правильность процесса вычислений.

После определения перемещений из решения СЛАУ (4) с помощью формул (2), представленных через перемещения в узлах ячейки конечноразностной сетки и уравнений Коши

ди ду дм ди ду ду дм дм ди

8х = я, , 8 У =~к , 82 = "я , ^ХУ = ^ + я , ^ У2 = я + "я , ^= я, + я , дх ду д2 ду дх д2 ду дх д2

определяются деформации в ячейках.

Затем по соотношениям обобщённого закона Гука вычисляются напряжения в ячейках конечно-разностной сетки

Сх = °118х + а128у + °1382 , Тху = а41 Уху ,

Су = а218х + а228у + а2382 , Ту2 = %Уу2, >

Сг = аз18х + аз28у + а3382 , Т2х = «61У гх , ,

где а11 = а22 = а33 = К + 4О/ 3; а41 = а51 = а61 = О;

а12 = а13 = а21 = а23 = а31 = а32 = К - 2О/3; К = Е/[3(1 - 2ц,)] — модуль объёмной деформации; Е — модуль упругости при растяжении-сжатии (модуль Юнга); ц — коэффициент поперечной деформации (коэффициент Пуассона);

О = Е [2 (1 + ц)] — модуль упругости при сдвиге.

При анализе НДС конструкций также необходимо знать величины интенсивностей напряжений и деформаций, которые определяются по формулам

°і = ^2 7(с --с г)2 +(с г-с *)2 +К- с -)2+6(т1у + т1+ті),

8. =-—л1(е -в )2 + (в -8 )2 + (в -8 )2 + —(у2 +у2 +у2)

і 3\ - у 2 ху 2Х

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

К достоинствам ВРМ также следует отнести следующие положения:

- простота математической формулировки задачи;

- ясный физический смысл используемого функционала;

- возможность использования метода для расчёта тел сложной формы, в том числе тел, неоднородных по деформационно-прочностным характеристикам материалов;

- сводимость проблемы к решению СЛАУ, для реализации которой существует достаточно надёжный математический аппарат линейной алгебры и вычислительные алгоритмы.

Численные результаты

Расчёт балки-стенки проводился для следующих значений геометрических размеров, физико-механических характеристик материала и внешних нагрузок: длина а = 8 м, высота Ь = 8 м, толщина Н = 0,5 м; модуль Юнга Е = 2-104 МПа, коэффициент Пуассона ц = 0,3; д = 2 МПа, д = 1,6 МПа; х1 = 1,75 м; х2 = 2,25 м; х3 = 3,00 м.

Вследствие наличия двух плоскостей симметрии для балки-стенки, внешней нагрузки и условий опирания рассматривалась только четверть конструкции АВСБ0Р0Н (рис. 3). На поверхностях 0АВР и 0Н0Р реализовывались условия симметрии перемещений и и ^ соответственно.

в н

0

Рис. 3. Расчётная область балки-стенки

Задача решалась при следующих статических и геометрических граничных условиях для четвёртой части объёма балки-стенки:

- нижняя поверхность у = 0 (0АБИ):

свободная поверхность 0 < х < х3 - с = т = тгу = 0, (7)

опорная поверхность х3 < х < а/2 - V = 0, ту=ту = 0; (8)

- верхняя поверхность у = Ь (ЕБСО):

поверхность х1 < х < х2 - су = я + , тху = тгу = 0, (9)

поверхности 0 < х < х1 и х2 < х < а/2 - су = я, тху = тгу = 0; (10)

- левая поверхность х = 0 (0АБЕ): и = 0, тух = тгх = 0;

- свободная правая поверхность х = а/2 (ИВСО): сх = т = тгх = 0;

- задняя поверхность г = 0 (0И0Е): м = 0, тхг = т = 0;

- свободная передняя поверхность

г = И/2 (АБСБ): сг = тхг = т„ =0. (11)

Тестирование программы расчёта НДС проводилось на примере балки-стенки при действии внешней поверхностной нагрузки интенсивностью я .

Опирание конструкции осуществлялось для значения х3 = 0 всей нижней поверхностью у = 0, на которой реализовывались граничные условия равенства нулю касательных напряжений и перемещений по оси :

V = 0, т = т = 0.

’ ху гу

В балке-стенке для этого варианта реализуется одномерное упругое НДС. Во всём объёме конструкции напряжения су = -2 МПа, а величины остальных напряжений по модулю примерно на четыре порядка меньше напряжения су. Линейные деформации, которые по абсолютной величине примерно на три порядка больше угловых, равны ву =-0,1 • 10-3, вх =вг = 0,3 -10-4. Перемещения вдоль оси у на верхней поверхности у = Ь балки-стенки V = —0,8 -10-3 м. Интенсивности напряжений и деформаций равны сi = 2 МПа, в; = 0,87 -10-4. Численно полученные результаты практически точно соответствуют результатам, получаемым при расчёте колонны по формулам сопротивления материалов для случая центрального сжатия. В расчётах выход материала в область пластического деформирования контролировался по величине интенсивности деформаций в;, которая сравнивалась с величиной деформации начала текучести вТ = 0,002.

Учёт действия дополнительной нагрузки qf, представленной на рис. 1,

приводит к адекватным изменениям в распределении параметров НДС в балке-стенке. Наибольшие напряжения и деформации имеют место в зоне приложения нагрузки : сх =-1,22 МПа, су=-3,60 МПа, сг = 0,17 МПа,

вх = 0,47-10-4, ву =-0,17• 10-3, вг = 0,70-10-4, с, = 3,35 МПа, в, = 0,15-10-3. Напряжения с и сх, действующие в координатной плоскости X07 , значительно превышают по величине напряжение сг. Поэтому можно говорить о двумерном, плоском напряжённом состоянии, реализуемом на этом участке конструкции. Что касается остальной, большей части объёма балки-стенки, то её НДС можно считать одномерным, т. к. сжимающие напряжения с значительно превосходят по величине остальные напряжения. Представленные результаты для двух вариантов задачи получены на конечно-разностной сетке х у х к) = (51^93x5) при реализации системы примерно 71000 линейных алгебраических уравнений. Вычисления перемещений проводились с погрешностью 0,1 %.

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

Анализ напряжений и деформаций выявил наиболее нагруженный по всей толщине конструкции участок - зону на нижней поверхности балки-стенки у края опоры х = х3 . Наибольшие значения нормальных напряжений, которые имеют место на оси х (т. е. для г = 0), равны: с = - 38,38 МПа, сх = -16,66 МПа, с г = -11,86 МПа. Касательное напряжение в координатной плоскости X07 тху = 8,80 МПа, и является максимальным по сравнению

с другими касательными напряжениями: оно в 3-6 раз превосходит их. Линейные деформации увеличились примерно на порядок, и у внешней поверхности г = И/ 2 достигают следующих значений: в х = 0,45 -10-3,

в =-0,15 • 10-2, вг = 0,52 -10-3. Величины угловых деформаций возросли на шесть-семь порядков и стали сопоставимы с линейными деформациями: у ху = 0,11 • 10-2, у уг = -0,18 • 10-3, у а = 0,35 • 10-3. Интенсивности напряжений и деформаций по сравнению с тестовым вариантом задачи увеличились примерно на порядок - с, = 26,90 МПа, в, = 0,12 -10-2.

На рис. 4 представлено распределение напряжений с на нижней поверхности у = 0, а на рис. 5 - на оси х (т. е. для г = 0) и пяти значений координаты у : кривая 1 - непосредственно у нижней поверхности для у = 0,02 м; кривая 2 - для у = 0,22 м; кривая 3 - для у = 0,42 м; кривая 4 - для у = 0,82 м; кривая 5 - у верхней поверхности балки-стенки для у = 7,98 м.

Как следует из результатов, представленных на рис. 5, при удалении от нижней поверхности балки-стенки, на которой реализуются граничные условия (7), (8), абсолютная величина напряжений с быстро уменьшается, и уже

для у = 0, 82 м (кривая 4) наибольшее значение напряжения с = - 6,3 МПа,

а интенсивность напряжений с, = 7,0 МПа. На верхней поверхности у = Ь

балки-стенки для напряжений с (кривая 5) выполняются граничные условия (9), (10), а максимальная величина интенсивности напряжений с{ = 3,2 МПа.

Рис. 4. Распределение напряжений ау в плоскости X02 для координаты у = 0

-40

-36

-32

-28

-24

-20

-16

-12

-8

-4

0

4

2 \ -ту

--

А"; ^ ' *" а/2 X

*2

Х3

(Г, , МПа

Рис. 5. Распределение напряжений с по координатам х и у

Распределение нормальных к координатной плоскости X07 напряжений сг по толщине балки-стенки для х = х3 и трёх значений координаты у представлено на рис. 6, где 1 = 2г/И - безразмерная величина.

В связи с тем, что эпюры напряжений строятся для середины ячеек конечно-разностной сетки, представленные на рис. 6 кривые не выходят на внешнюю поверхность г = И/2 конструкции. Если кривые зависимости сг (г) экстраполировать до значения 1 = 1,0, т. е. выйти на поверхность бал-

ки-стенки, то величины напряжений с2 будут удовлетворять граничным условиям (11). Следует отметить, что изменения напряжений сг вдоль оси X у края опоры х = х3 имеют такой же характер, как и напряжений с .

1 О-, МПа

0 0,25 0,50 0,75 1,00

Рис. 6. Распределение напряжений с2 по толщине конструкции для у = 0,02 м (кривая 1); у = 0,06 м (кривая 2); у = 0,10 м (кривая 3)

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

На рис. 7 представлено распределение напряжений сх и су по высоте балки-стенки для значений координат х = 2 = 0. Эпюры напряжений качественно согласуются с результатами других авторов, например с результатами работ [8-13], полученными при реализации задачи в плоской двумерной постановке.

У, м

7

6

1 5

\ 4

\ 3

ч 2

стД

0

-2-1012 з<ХМПа

Рис. 7. Распределение напряжений сх и с по высоте балки-стенки вдоль оси у для значений х = 2 = 0

Из анализа распределения остальных напряжений следует, что во всём объёме балки-стенки, кроме небольшой области на нижней поверхности конструкции у края опоры х = х3, реализуется плоское напряжённое состояние, при котором можно пренебречь напряжениями с2, т, ту2.

Расчёт НДС балки-стенки на конечно-разностной сетке размером (101x201^9) при реализации системы примерно 548000 линейных алгебраических уравнений при вычислении перемещений с погрешностью 20 % был проведён за 7 часов на ПК с процессором ЛМБ ЛШоп с тактовой частотой 2,6 ГГц. Для решения задачи с погрешностью 10 % потребовалось уже 26 часов. При этом максимальные значения перемещений, деформаций и напряжений для этих вариантов имеют место в одних и тех же ячейках сетки, а их численные значения отличаются друг от друга в третьей или четвёртой значащей цифре.

Если для реализации СЛАУ воспользоваться методом исключения Гаусса, то вследствие симметрии матрицы коэффициентов {А} при решении

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

Выводы

Результаты численного моделирования упругого напряжённо-деформированного состояния балки-стенки в трёхмерной постановке подтвердили очевидные предположения о наличии опасных с точки зрения прочности зон при деформировании конструкции: под нагрузкой qf и особенно

у края опоры х = х3. В первой из перечисленных зон реализуется плоское напряжённое состояние, во второй - существенно трёхмерное НДС.

В рассмотренном варианте задачи максимальная величина напряжения сг во второй опасной зоне балки-стенки составляет 31 % от величины максимального напряжения с . Поэтому плоская постановка задачи о НДС балки-стенки, связанная с пренебрежением напряжениями с2, тх2, т , приведёт

к существенной погрешности при определении параметров НДС, в частности интенсивности напряжений сi во второй зоне. Именно эта величина часто используется для оценки уровня напряжённости материала конструкции и её прочности. Расчёты, проведённые для других значений толщины Н конструкции, показали, что при её уменьшении отношение напряжений с^су также

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

Библиографический список

1. Барашков, В.Н. Алгоритм реализации задачи теории упругости и пластичности вариационно-разностным методом. Ч. I / В.Н. Барашков // Известия Томского политехнического университета. - 2003. - Т. 306. - № 3. - С. 23-28.

2. Барашков, В.Н. Алгоритм реализации задачи теории упругости и пластичности вариационно-разностным методом. Ч. II / В.Н. Барашков // Известия Томского политехнического университета. - 2003. - Т. 306. - № 4. - С. 23-27.

3. Деруга, А.П. Об одном способе определения оптимального коэффициента сверхрелаксации / А.П. Деруга, А.А. Ларионов // Пространственные конструкции в Красноярском крае: сб. статей. - Красноярск : Изд-во кПи, 1978. - Вып. 11. - С. 175-178.

4. Вазов, В. Разностные методы решения дифференциальных уравнений в частных производных / В. Вазов, Дж. Форсайт. - М. : Иностр. лит-ра, 1963. - 488 с.

5. Самарский, А.А. Теория разностных схем / А.А. Самарский. - М. : Наука, 1977. - 656 с.

6. Varga, R.S. Matrix iterative analysis / R.S. Varga. - Prentice Hall, Englewood Cliffs, New Jer-sy, S. n., 1962.

7. Теоретические и экспериментальные исследования высокоскоростного взаимодействия тел / под ред. А.В. Герасимова. - Томск : Изд-во Том. ун-та, 2007. - 572 с.

8. Синицын, А.П. Решение задачи о плоском напряженном состоянии в конечных разностях / А.П. Синицын // Вестник Военно-инженерной академии РККА, 3. - М. : Изд-во ВИА РККА, 1934. - С. 15-30.

9. Гольденблат, И.И. Расчет и конструирование железобетонных балок / И.И. Гольденб-лат. - М. : Стройиздат, 1940. - 84 с.

10. Варвак, П.М. Развитие и приложение метода сеток к расчету пластинок. Ч. 1 / П.М. Вар-вак. - Киев : Изд-во АН УССР. - 1949. - 136 с.

11. Варвак, П.М. Развитие и приложение метода сеток к расчету пластинок. Ч. 2 / П.М. Варвак. - Киев : Изд-во АН УССР. - 1952. - 116 с.

12. Жемочкин, Б.Н. Теория упругости / Б.Н. Жемочкин. - М. : Госстройиздат, 1957. - 256 с.

13. Вайнберг, Д.В. Пластины, диски, балки-стенки / Д.В. Вайнберг, Е.Д. Вайнберг. - Киев : Изд-во АН УССР, 1959. - 1049 с.