Рассматриваются вопросы численного моделирования взрывного воздействия на упругую полуплоскость с полостью (соотношение ширины к высоте один к двенадцати). Поставленная задача решается с помощью численного моделирования уравнений динамики сплошных сред.
Волны напряжений различной природы, распространяясь, в деформируемом теле взаимодействуют, друг с другом, что приводит к образованию новых областей возмущений, перераспределению напряжений и деформаций.
При интерференции волн напряжений их интенсивности складываются. Они могут достигать значений, превосходящих предел прочности материала. В этом случае наступает разрушение материала.
После трехкратного или четырехкратного прохождения и отражения волн напряжений в теле процесс распространения возмущений становится установившимся, напряжения и деформации усредняются, тело находится в колебательном движении.
Расчеты проводились при следующих единицах измерения: килограмм-сила (кгс); сантиметр (см); секунда (с). Для перехода в другие единицы измерения были приняты следующие допущения: 1 кгс/см2 ≈ 0,1 МПа; 1 кгс•с2/см4 ≈ 109 кг/м3.
В работах [1–10] приведена информация о моделировании нестационарных волн напряжений в деформируемых телах сложной формы с помощью применяемого численного метода, алгоритма и комплекса программ.
Некоторая информация о физической достоверности и математической точности рассматриваемого численного метода, алгоритма и комплекса программ приведена в следующих работах [1, 5, 8–10].
Для решения задачи о моделировании упругих нестационарных волн напряжений в деформируемых областях сложной формы рассмотрим некоторое тело Γ в прямоугольной декартовой системе координат XOY, которому в начальный момент времени t = 0 сообщается механическое нестационарное импульсное воздействие.
Предположим, что тело Γ изготовлено из однородного изотропного материала, подчиняющегося упругому закону Гука при малых упругих деформациях.
Точные уравнения двумерной (плоское напряженное состояние) динамической теории упругости имеют вид
(1)
где σx, σy и τxy – компоненты тензора упругих напряжений; εx, εy и γxy – компоненты тензора упругих деформаций; u и v – составляющие вектора упругих перемещений вдоль осей OX и OY соответственно; ρ – плотность материала; 
– скорость продольной упругой волны; 
– скорость поперечной упругой волны; ν – коэффициент Пуассона; E – модуль упругости; 
– граничный контур тела Γ.
Систему (1) в области, занимаемой телом Γ, следует интегрировать при начальных и граничных условиях.
Для решения двумерной нестационарной динамической задачи математической теории упругости с начальными и граничными условиями (1) используем метод конечных элементов в перемещениях.
Принимая во внимание определение матрицы жесткости, вектора инерции и вектора внешних сил для тела Γ, записываем приближенное значение уравнения движения в теории упругости
(2)
где 
– диагональная матрица инерции; 
– матрица жесткости; 
– вектор узловых упругих перемещений; 
– вектор узловых упругих скоростей перемещений; 
– вектор узловых упругих ускорений; 
– вектор внешних узловых упругих сил.
Интегрируя уравнения (2) конечноэлементным вариантом метода Галеркина, получим явную двухслойную конечноэлементную линейную схему в перемещениях для внутренних и граничных узловых точек
(3)
Шаг по временной переменной координате Δt выбирается из следующего соотношения
(i = 1, 2, 3, ...), (4)
где Δl – длина стороны конечного элемента.
На основе метода конечных элементов в перемещениях разработаны алгоритм и комплекс программ для решения линейных плоских двумерных задач, которые позволяют решать задачи при нестационарных волновых воздействиях на сложные системы. При разработке комплекса программ использовался алгоритмический язык Фортран-90.
Рис. 1. Постановка задачи о воздействии сосредоточенной взрывной волны на свободной поверхности упругой полуплоскости с полостью (соотношение ширины к высоте один к двенадцати)
Исследуемая область разбивается по пространственным переменным на треугольные конечные элементы с тремя узловыми точками с линейной аппроксимацией упругих перемещений и на прямоугольные конечные элементы с четырьмя узловыми точками с билинейной аппроксимацией упругих перемещений. По временной переменной исследуемая область разбивается на линейные конечные элементы с двумя узловыми точками с линейной аппроксимацией упругих перемещений.
Рис. 2. Изменение упругого контурного напряжения 
во времени t/Δt в точке A1: 1 – в задаче без полости; 2 – в задаче с полостью (соотношение ширины к высоте один к двенадцати)
Рис. 3. Изменение упругого контурного напряжения 
во времени t/Δt в точке A2: 1 – в задаче без полости; 2 – в задаче с полостью (соотношение ширины к высоте один к двенадцати)
Рис. 4. Изменение упругого контурного напряжения 
во времени t/Δt в точке A3: 1 – в задаче без полости; 2 – в задаче с полостью (соотношение ширины к высоте один к двенадцати)
Рис. 5. Изменение упругого контурного напряжения 
во времени t/Δt в точке A4: 1 – в задаче без полости; 2 – в задаче с полостью (соотношение ширины к высоте один к двенадцати)
Рассмотрим задачу о воздействии сосредоточенной взрывной волны на свободной поверхности упругой полуплоскости с полостью (соотношение ширины к высоте один к двенадцати) (рис. 1).
В точке F перпендикулярно свободной поверхности ABEFG приложено сосредоточенное нормальное напряжение σy (рис. 1), которое при 0 ≤ n ≤ 10 (n = t/Δt) изменяется линейно от 0 до P, а при 10 ≤ n ≤ 20 от P до 0 (P = σ0, σ0 = –1 МПа (–1 кгс/см2)).
Граничные условия для контура GHIA при t > 0 
. Отраженные волны от контура GHIA не доходят до исследуемых точек при 0 ≤ n ≤ 200. Контур ABCDEFG свободен от нагрузок, кроме точки F, где приложено сосредоточенное упругое нормальное напряжение σy.
Расчеты проведены при следующих исходных данных: H = Δx = Δy; Δt = 1,393•10–6 с; E = 3,15•104 МПа (3,15•105 кгс/см2); ν = 0,2; ρ = 0,255•104 кг/м3 (0,255•10–5 кгс•с2/см4); Cp = 3587 м/с; CS = 2269 м/с. Решается система уравнений из 59048 неизвестных.
Результаты расчетов для контурного напряжения 
во времени n получены в точках A1–A4 (рис. 1), находящихся на свободной поверхности упругой полуплоскости.
На рис. 2–5 приведены контурные напряжения 
во времени n в точках A1–A4: 1 – в задаче без полости; 2 – в задаче с полостью (соотношение ширины к высоте один к двенадцати).
Выводы
1. Для прогноза безопасности объекта, находящегося в твердой деформируемой среде, при волновых воздействиях применяется численное моделирование.
2. Разработаны методика, алгоритм и комплекс программ для решения линейных двумерных плоских задач, которые позволяют решать сложные задачи при волновых воздействиях.
3. Линейная динамическая задача с начальными и граничными условиями в виде дифференциальных уравнений в частных производных приведена к системе линейных обыкновенных дифференциальных уравнений с начальными условиями, которая решается по явной двухслойной схеме.
4. Рассмотрена задача о воздействии сосредоточенной взрывной волны на свободной поверхности упругой полуплоскости с полостью (соотношение ширины к высоте один к двенадцати).
5. Полученные результаты можно оценить как первое приближение к решению сложной комплексной задачи о воздействии сосредоточенной взрывной волны на свободной поверхности полуплоскости с полостью (соотношение ширины к высоте один к двенадцати), с помощью численного моделирования волновых уравнений теории упругости.