Numerical simulation of the interaction of a shock wave with a permeable deformable granulated layer

Full Text

Введение. Перспективным элементом, защищающим конструкции от импульсных воздействий, являются насыпные слои [1-3]. В научной литературе существует ряд работ, посвященных исследованию деформируемых пористых преград [4-22], где показано, что размещение слоя пористых насыпных слоев перед стенкой мишени уменьшает величины остаточного напряжения и энергии мишени и существенно снижает скорость приложенной нагрузки. Полученные результаты ограничиваются малыми нагрузками, которые не вызывают развитое пластическое течение или разрушение материала шариков [14, 15]. При этом предполагалось, что проницаемые элементы в процессе взаимодействия с ударными волнами испытывают малые деформации и их проницаемость не изменяется. В данной статье рассматривается процесс взаимодействия плоских ударных волн с деформируемыми упругопластическими гранулированными слоями, которые могут испытывать большие деформации, влияющие на их проницаемость.

1. Уравнения в одномерном приближении. Динамическое поведение гранулированного слоя с содержащимся поровым газом описывается на основе уравнений динамики двух взаимопроникающих континуумов, каждый из которых имеет свои скорости, напряжения (давления) и температуры. Одномерные уравнения динамики порового газа в форме законов сохранения массы, импульса и энергии имеют вид:

$\frac{\partial \left({\alpha }_1{\rho }_1^0\right)}{\partial t}+\frac{\partial }{\partial х}\left(\beta {\rho }_1^0{u}_1\right)=\mathrm{0,}\frac{\partial \left({\alpha }_1{\rho }_1^0{u}_1\right)}{\partial t}+\frac{\partial }{\partial х}\left(\beta \right(p_1+{\rho }_1^0{u}_1^2\left)\right)=\mathrm{0,}$

$\frac{\partial }{\partial t}\left({\alpha }_1{\rho }_1^0\right(e_1+\frac{u_1^2}{2}\left)\right)+\frac{\partial }{\partial x}\left\{\beta \right({\rho }_1^0{u}_1(e_1+\frac{u_1^2}{2})+p_1{u}_1\left)\right\}=\mathrm{0,}$

$p_1=(\gamma -1){\rho }_1^0{e}_1,T_1=e_1/c_{v_1}\left(T_1\right).$

Здесь и далее нижний индекс 1 относится к газу; 2 — к твердой компоненте; $t$ — время; $x$ — координата; ${\rho }_1^0$ — истинная плотность газа; $u$ — скорость; $e$ — внутренняя энергия; $c_{v_1}$ — удельная теплоемкость; $\gamma$ — показатель адиабаты газа; $f$, $q$ — межфазные силы и тепловой поток, действующие на "элементарный узел" твердой компоненты со стороны порового газа; $n$ — количество таких узлов в единице объема смеси; ${\alpha }_i$ — объемные концентрации компонент; $\beta$ — коэффициент проницаемости сечений. Коэффициент проницаемости принимается в виде отношения площади пор к общей площади элемента среды. Из системы (1) при ${\alpha }_1=\beta =1$ и $n=0$ следуют уравнения газовой динамики для однородной среды.

Уравнения динамического деформирования гранулированного слоя в одномерном приближении как скелета двухфазной среды имеют вид [14, 23]:

$\frac{\partial \left({\rho }_2\right)}{\partial t}+\frac{\partial }{\partial x}\left({\rho }_2{u}_2\right)=\mathrm{0,}$

$\frac{\partial \left({\rho }_2{u}_2\right)}{\partial t}+\frac{\partial }{\partial x}\left({\rho }_2{u}_2^2\right)=nf,$

$\frac{\partial \sigma }{\partial t}+u_2\frac{\partial \sigma }{\partial x}=\frac{\partial f_n}{\partial \epsilon }\frac{\partial u_2}{\partial x}.$ (2)

Здесь $\sigma$, $\epsilon$ — напряжение и деформация, ${\rho }_2$ — приведенная плотность гранулированного слоя ( ${\rho }_2={\alpha }_2{\rho }_2^0$, ${\alpha }_1+{\alpha }_2=1$). Закон связи между напряжением и деформацией имеет вид $\sigma =f_n(\epsilon ,{\sigma }^{*})$. Параметры, отмеченные верхним индексом , представляют собой максимальные значения, достигнутые частицей при нагружении в соответствующем направлении, они необходимы для описания разгрузки частиц среды, испытывающей необратимые деформации. Конкретный вид уравнения состояния определяется по результатам экспериментальных исследований сжатия гранулированных слоев. При взаимодействии твердого и газового компонентов в качестве межфазных сил учитываются: силы сопротивления частичек твердой фазы при их обтекании поровым газом; силы Стокса вязкого трения, а также конвективный теплообмен через межфазную поверхность [23]. Твердая фаза как скелет высокопористой среды может сильно сжиматься в процессе деформации, поэтому параметры ${\alpha }_1$, $\beta$ будут зависеть от степени ее сжатия. Конкретный вид этих зависимостей определяется численными исследованиями деформирования фрагментов гранулированного слоя в трехмерной постановке [24].

Поверхности контакта чистого газа с пористым гранулированным слоем являются поверхностью разрыва пористости и проницаемости. Как показали исследования [25-27], на контактных границах "чистый газ" – "поровый газ" должны выполняться специальные условия, как на скачке пористости:

${\rho }_5^0{u}_5=\beta {\rho }_4^0{u}_4,$

${\rho }_5^0{u}_5^2+\beta p_5=\beta ({\rho }_4^0{u}_4^2+p_4)$

${\rho }_5^0{u}_5(\frac{\gamma }{\gamma -1}\frac{p_5}{{\rho }_5^0}+\frac{u_5^2}{2})=\beta {\rho }_4^0{u}_4(\frac{\gamma }{\gamma -1}\frac{p_4}{{\rho }_4^0}+\frac{u_4^2}{2}).$ (3)

Уравнения (3) отражают законы сохранения массы, импульса и энергии на скачке. Индекс 5 относится к параметрам со стороны "чистого газа", индекс 4 — к параметрам со стороны порового газа.

Построенная нелинейная математическая модель описывает процессы взаимодействия проницаемых деформируемых гранулированных слоев с ударными волнами. Модель учитывает изменение пористости среды и ее проницаемости от степени деформации гранулированного слоя. Численное решение уравнений (1), (2) производится с помощью схемы С. К. Годунова [28], адаптированной к задачам динамики взаимопроникающих сред. Численная реализация контактных условий (3) производится на основе решения задачи распада разрыва на скачке пористости [25, 27].

2. Постановка задачи. Постановка задачи показана на рис. 1. Расчетная область включает в себя три подобласти: воздух, гранулированный слой, воздух. Подвижный деформируемый гранулированный слой размещается в подобласти 2 от $x=x_2$ до $x=x_3$, таким образом, начальная толщина гранулированного слоя составляет $H=x_3-x_2=0.017$ м. Коэффициент проницаемости слоя принимается равным $\alpha =0.215$, начальные параметры порового газа $p_0=0.1$ МПа, ${\rho }_0=1.23$ кг/м , $u_0=0$.

 

Рис. 1. Постановка задачи численного моделирования

 

Координаты границ подобластей имеют следующие значения: $x_1=-0.2$ м, $x_2=0$ м, $x_3=0.017$ м, $x_4=0.217$ м. В первой подобласти задаются параметры, соответствующие параметрам газа (воздуха) за фронтом набегающей плоской ударной волны. На искусственных границах подобластей газа $x=x_1$ и $x=x_4$ ставятся условия по давлению, соответствующие начальным условиям по этим подобластям. Размеры расчетных областей (подобласти 1 и 2) выбраны из условий, чтобы волновые возмущения от гранулированного слоя не отразились от искусственных границ и не исказили численное решение в течение интервала времени 0.4 мс, когда уже сформированы отраженные и прошедшие ударные волны в газе.

Коэффициент проницаемости принимается в виде отношения площади пор к общей площади элемента среды. В третьей подобласти — покоящийся газ с начальными параметрами, как и во второй подобласти. Показатель адиабаты газа во всех подобластях $\gamma =1.4$. Решение получено на разностной сетке с размером ячеек 0.0005 м. Гранулированный слой предполагается деформируемым. Кривые одноосного сжатия имеют вид, характерный для пористых материалов [29], и в переменных "давление – плотность" показаны на рис. 2. Пунктиром показана разгрузочная ветвь, принимаемая в виде прямой линии, с тангенсом угла наклона (квадрат скорости звука) равным $2.9\cdot {10}^5$ м /с . В начальный момент времени слой покоится, деформации слоя отсутствуют, его начальная плотность равна 680 кг/м , давление — 0.1 МПа. Уменьшение проницаемости слоя при его сжатии с начального значения $\alpha =0.215$ происходит по линейному закону (рис. 3). Вид зависимости изменения проницаемости от плотности установлен по результатам численного моделирования трехмерных задач упругопластического сжатия симметричных фрагментов гранулированных слоев [30].

 

Рис. 2. Изменение плотности относительно давления

 

Рис. 3. Изменение проницаемости

 

В области $x<0$ в качестве начальных условий задаются постоянные параметры за фронтом набегающей плоской ударной волны (УВ) в направлении оси $Ox$. Эти параметры определяются числом Маха ударного фронта $M_0$ и вычисляются согласно формулам [31]:

$\frac{p_{sw}}{p_0}=\frac{2\gamma M_0^2-(\gamma -1)}{\gamma +1},\frac{{\rho }_{sw}}{{\rho }_0}=\frac{(\gamma +1)M_0^2}{(\gamma -1)M_0^2+2},$

$\frac{T_{sw}}{T_0}=\frac{(2\gamma M_0^2-(\gamma -1\left)\right)\left(\right(\gamma -1)M_0^2+2)}{{(\gamma +1)}^2{M}_0^2},$

$\frac{u_{sw}}{a_0}=\frac{2}{\gamma +1}(M_0-\frac{1}{M_0}),$

где $\gamma$ — показатель адиабаты, $M_0=u_0/a_0$, $a_0$ — скорость звука перед УВ.

Рассматривались три варианта задания интенсивности набегающей УВ. В первом варианте полагалось $p_{sw}=0.3$ МПа, ${\rho }_{sw}=2.596$ кг/см , $u_{sw}=292.5$ м/с. На рис. 4–6 представлены распределения давлений, плотностей и скоростей в газе по расчетной области задачи ($-0.1м<x<0.1м$) в момент времени $t=0.2$ мс, когда формируются отраженные и проходящие через гранулированный слой ударные волны. Видно, что амплитуда давления отраженной волны в газе (кривые 1) более чем в 2 раза превышает амплитуду падающей ударной волны. Проходящая волна также является нелинейной с амплитудой 0.136 МПа (кривая 3). Наибольшие давления образуются в поровом газе гранулированной среды с амплитудами, на порядок превышающими амплитуду набегающей, отраженной и проходящей волны (кривая 2). Распределения скорости и плотности подтверждают нелинейный характер протекающих процессов (рис. 5, 6). Хорошо видны характерные скачки всех параметров на границах гранулированного слоя, что соответствует физике протекающих явлений [14].

 

Рис. 4. Распределение давления по расчетной области при $t=0.2$мс

 

Рис. 5. Распределение скорости по расчетной области при $t=0.2$мс

 

Рис. 6. Распределение плотности по расчетной области при $t=0.2$мс

 

На рис. 7, 8 представлены распределения давлений и плотностей по гранулированному слою (твердой фазе) в моменты времени $t=\{0.1;0.2;0.3\}$ мс. Здесь цифрой 1 отмечено распределение полей в момент времени $t=0.1$ мс; цифрой 2 — $t=0.2$ мс; цифрой 3 — $t=0.3$ мс.

 

Рис. 7. Распределение давления по расчетной области в твердой фазе слоя

 

Рис. 8. Распределение плотности по расчетной области в твердой фазе слоя

 

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

Для этого проведены расчеты взаимодействия гранулированного слоя с набегающими УВ различной амплитуды: $p_{sw}=0.3$ МПа, $p_{sw}=0.2$ МПа, $p_{sw}=0.15$ МПа.

На рис. 9, 10 представлены временные зависимости давлений проходящих (рис. 9) и отраженных (рис. 10) волн в точках на удалении  и  от границ слоя.

 

Рис. 9. Изменение давления проходящих волн

 

Рис. 10. Изменение давления отраженных волн

 

Цифрами 1 и 1' отмечено решение, полученное при параметрах газа за фронтом плоской ударной волны: $p_{sw}=0.3$ МПа. Цифрами 2 и 2' — решение, полученное при параметрах газа за фронтом плоской ударной волны: $p_{sw}=0.2$ МПа. Цифрами 3 и 3' — решение, полученное при параметрах газа за фронтом плоской ударной волны: $p_{sw}=0.2$ МПа. Сплошные линии (1, 2, 3) соответствуют решению с учетом изменения коэффициента проницаемости, пунктирные линии (1', 2', 3') "— решению с постоянным коэффициентом проницаемости 0.215. Влияние учета изменения проницаемости особо заметно на параметрах отраженных волн, и это влияние возрастает с ростом амплитуды набегающей волны. Эта же закономерность проявляется и для проходящих волн, но в меньшей мере.

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

Конкурирующие интересы. Заявляем, что в отношении авторства и публикации этой статьи конфликта интересов не имеем.

Авторский вклад и ответственность. И. А. Модин — проведение численных расчетов и экспериментальных исследований, обработка и анализ результатов, работа с черновиком и переработанным вариантом рукописи. А. В. Кочетков — идея исследования, формулировка целей и задач исследования, визуализация и верификация результатов Е. Г. Глазова — проведение численных расчетов, визуализация и верификация результатов, работа с черновиком и переработанным вариантом рукописи. Авторы несут полную ответственность за предоставление окончательной рукописи в печать. Окончательная версия рукописи была одобрена всеми авторами.

Финансирование. Работа выполнена при финансовой поддержке гранта Российского научного фонда (РНФ № 20–79–00108).

About the authors

Ivan A. Modin

Research Institute of Mechanics, National Research Lobachevsky State University of Nizhny Novgorod

Email: mianet@mail.ru
ORCID iD: 0000-0002-3561-4606
SPIN-code: 4839-8129
Scopus Author ID: 57192279101
ResearcherId: E-9088-2019
http://www.mathnet.ru/rus/person138504

Cand. Techn. Sci.; Researcher; Researcher; Lab. of Simulation of Physical and Mechanical Processes

Russian Federation, 23, korp. 6, pr. Gagarina, Nizhny Novgorod, 603022

Anatoliy V. Kochetkov

Research Institute of Mechanics, National Research Lobachevsky State University of Nizhny Novgorod

Author for correspondence.
Email: kochetkov@mech.unn.ru
ORCID iD: 0000-0001-7939-8207
SPIN-code: 9964-3450
Scopus Author ID: 23004869700
http://www.mathnet.ru/person32889

Dr. Phys. & Math. Sci.; Professor; Dept. of Theoretical, Computer and Experimental Mechanics; (Institute of Information Technology, Mathematics and Mechanics UNN)

Russian Federation, 23, korp. 6, pr. Gagarina, Nizhny Novgorod, 603022

Elena G. Glazova

Research Institute of Mechanics, National Research Lobachevsky State University of Nizhny Novgorod

Email: glazova@mech.unn.ru
ORCID iD: 0000-0003-4351-889X
SPIN-code: 5314-3012
Scopus Author ID: 55248346700
http://www.mathnet.ru/person163935

Cand. Phys. & Math. Sci.; Senior Researcher; Lab. of Dynamics of Multicomponent Media (Research Institute of Mechanics UNN)

Russian Federation, 23, korp. 6, pr. Gagarina, Nizhny Novgorod, 603022

References

Supplementary files