Колпаков_P. # 02, февраль 2012
 Скачать 1.44 Mb.
|
|
http://technomag.edu.ru/doc/334177.html 1 Математическое моделирование функционирования взрывных устройств # 02, февраль 2012 Колпаков В. И. УДК 519.63:532:539.5 МГТУ им. Н.Э. Баумана Взрывные или импульсные устройства в настоящее время широко используются в боеприпасах и средствах поражения, в ракетной технике, в работах по интенсификации скважин при газо- и нефтедобыче. К ним, например, можно отнести фугасные, осколочные и кумулятивные заряды разного назначения. Однако, несмотря на то, что экспериментальные исследования играют ключевую роль в изучении импульсных устройств и технологий, использующих в своей основе взрывчатые вещества (ВВ), без глубокого теоретического анализа, как правило, не удается достичь требуемого результата. Кроме того, современные условия диктуют необходимость существенного сокращения количества испытаний при их отработке. В этой связи существенно возросли значение и практическая ценность исследований, проводимых на основе численных методов механики сплошной среды, которые, в свою очередь, предъявляют повышенные требования не только к качеству физико-математических моделей, но и уровню разработанных на их основе алгоритмов расчета. В рассмотренном ключе в настоящей работе излагается методика численного расчета, хорошо зарекомендовавшая себя при моделировании функционирования взрывных устройств различного назначения в течение последних пятнадцати лет. 1. ФИЗИКО-МАТЕМАТИЧЕСКАЯ ПОСТАНОВКА ЗАДАЧИ Изучение процессов метания оболочки осколочных макетов или обжатия кумулятивной облицовки (КО) продуктами детонации (ПД) с последующим формированием кумулятивной струи (КС) или поражающего элемента (ПЭ) как для осесимметричных, таки для удлиненных кумулятивных зарядов (КЗ) целесообразно проводить в двумерной постановке (в цилиндрической или декартовой системе координат соответственно. Есте- 77-30569/334177 , №02 февраль 2012 г. http://technomag.edu.ru 2 ственно, в первом случае не учитывается асимметрия в изготовлении и инициировании заряда ВВ, во втором — протяженность и пространственная кривизна реальных КЗ. Кумулятивные заряды, формирующие ПЭ, принято также называть снарядоформирующими зарядами (СФЗ). Типовые расчетные схемы рассматриваемых взрывных устройств показаны на рис. 1. Здесь 1 – ВВ, 2 – корпус КЗ или осколочного макета, 3 – точка инициирования (ТИ), 4 – КО, 5 – фронт детонационной волны (ДВ), 6 – ПД, 7 – линза (или экран, 8 – реперные точки (маркеры, используемые для идентификации течения в КО или осколочной оболочке. а) б) Рис. 1. Расчетные схемы взрывных устройства) КЗ; (б) осколочный макет, где 1 – ВВ, 2 – корпус КЗ или осколочного макета, 3 – точка инициирования, 4 – КО, 5 – фронт детонационной волны, 6 – ПД, 7 – линза (экран, 8 – реперные точки (маркеры) Так как для обеих схем постановка задачи примерно идентична, рассмотрим ее на примере КЗ. Для этого будем полагать, что в начальный момент времени (t = 0) в точке 3 точка инициирования) осуществляется подрыв заряда ВВ с начальной плотностью ρ вв и теплотой взрывчатого превращения Q. От точки инициирования начинает распространяться фронт ДВ (кривая 5) со скоростью D вв с образованием ПД (зона 6). Стечением времени ДВ начинает отражаться от поверхностей корпуса (позиция 2) и КО (позиция 4), на которые действуют давления порядка 20 … 60 ГПа. Его величина зависит от свойств ВВ, угла подхода фронта ДВ к поверхности облицовки, материала и толщины облицовки. Под действием ПД кумулятивная облицовка начинает обжиматься с образованием КС или ПЭ. При этом для получения общих закономерностей или особенностей их формирования для конкретного КЗ, обусловленных формой облицовки, геометрией заряда, формой и месторасположением линзы и ТИ, физико-механическими свойствами используемого http://technomag.edu.ru/doc/334177.html 3 состава ВВ или материалов облицовки и корпуса, целесообразно использовать модель сжимаемой идеальной упругопластической среды с баротропным уравнением состояния [1 – 4]. Вязкими свойствами материала облицовки в процессе формирования КС или ПЭ можно пренебречь. Использование баротропного уравнения состояния позволяет избежать интегрирования уравнения энергии в системе соотношений, описывающих поведение взаимодействующих сред. Причем если при формировании КС можно ограничиться газодинамической моделью, то для определения формы и кинематических характеристик ПЭ, особенно на поздних стадиях движения, является важным учет проявления упругопласти- ческих свойств материала. Для зоны ПД течение считается изэнтропическим, так как даже отражение фронта ДВ от абсолютно жесткой поверхности не приводит к заметному росту энтропии. Распространение детонации рассматривается вне общей системы уравнений и задает границу области, охваченную течением. С учетом сделанных допущений система уравнений, описывающая двумерное течение в переменных Эйлера, имеет вид 0 = ρ α + ∂ ρ ∂ + ∂ ρ ∂ + ∂ ρ ∂ r v z v r v t r z r ) ( ) ( (1) ( ) zz rr rz rr r z r r r r D D r z D r z v v r v v t v t d v d σ σ σ + α + ∂ ∂ + ∂ σ ∂ = ∂ ∂ + ∂ ∂ + ∂ ∂ ρ = ρ 2 (2) ( ) rz rz zz z z z r z z D r r D z z v v r v v t v t d v d σ σ α + ∂ ∂ + ∂ σ ∂ = ∂ ∂ + ∂ ∂ + ∂ ∂ ρ = ρ (3) ( ) [ ] ( ) [ ] − ρ ρ ρ + ρ < ρ ρ ⇒ σ − ≥ ρ ρ ⇒ − ρ ρ ⋅ = ρ = γ 1 C 1 1 1 ) ( 3 ВВ ВВ 0 * P 0 для КО, корпуса, линзы, (4) для КО, корпуса, линзы, для ПД, для ВВ; p D rr rr − = σ σ , p D zz zz − = σ σ , rz rz D σ = σ ; (5) ρ ρ + ∂ ∂ = σ dt d r v G Dt DD r rr 3 1 2 ; ρ ρ + ∂ ∂ = σ dt d z v G Dt DD z zz 3 1 2 ; ∂ ∂ + ∂ ∂ = σ r v z v G Dt DD z r rz ; (6) ( ) 2 2 2 2 3 2 2 Y D D D D D f zz rr zz rz rr ≤ ⋅ + + + ⋅ = σ σ σ σ σ (7) Здесь ρ – плотность p – давление t – текущее время v r , v z – компоненты вектора скорости по осями выбранной системы координат (рис. 1); σ rr , σ zz – нормальные и σ rz – касательные напряжения D σrr , D σrz , D σzz – компоненты девиатора напряжений 77-30569/334177 , №02 февраль 2012 г. http://technomag.edu.ru 4 D(…) / Dt – производная в смысле Яуманна, учитывающая поправки к составляющим де- виатора напряжений, обусловленные поворотом фиксированного элемента среды как целого текущие значения модуля сдвига и динамического предела текучести материала среды, принимаемых в областях разгрузки в виде функциональных зависимостей от их начальных значений и безразмерной плотности γ = ρ /ρ 0 , причем при G 0 = 0 и Y 0 = 0 исходная система уравнений естественным образом трансформируется в уравнения идеального газа, описывающие поведение ПД; α – коэффициент симметрии для плоского случая α = 0, для осесимметричного случая α = 1). В приведенной системе уравнений соотношения (1) – (3) представляют собой, соответственно, законы сохранения массы и импульса (4) – баротропные уравнения состояние взаимодействующих сред (КО, корпуса, линзы, ВВ и ПД), конкретизируемые ниже. Соотношения (5) связывают компоненты тензора полных напряжений с шаровой и девиа- торной составляющими (6) – закон Гука в дифференциальной форме (7) – условие пластического течения Мизеса σ i ≤ Y (σ i – интенсивность напряжений. Если условие (7) нарушается, те. (f > 2Y 2 /3) и материал деформируемой среды находится в состоянии пластического течения, то компоненты D σrr , D σrz , умножаются на множитель процедура приведения напряжений к кругу текучести [6]). Для описания разрушения металлических элементов конструкции зарядов под действием ПД использовалась комбинация критериев откольной прочности p = – σ * p ( σ * p – от- кольная прочность) и предельных пластических деформаций р р (р предельная пластическая деформация. Первый из них использовался на этапе нагружения элементов конструкции ДВ, второй – при последующем пластическом деформировании. Причем при ρ/ρ 0 < 1 и p = – σ * p – величина давления ограничивалась в соответствии с условиями [7, 8] ( ) ( ) γ Φ ρ = p p ; ( ) γ Φ = 0 Y Y ; ( ) γ Φ = 0 G G , (8) где ( ) ( ) ( ) γ ≤ γ γ ≤ γ < γ γ γ γ γ γ > γ = γ Φ ; для , для - - для 2 2 1 2 1 2 1 0 , (9) γ 1 , γ 2 – постоянные материала, числовые значения которых, например, для стали поданным работ [9] составляют – γ 1 = 0.95, γ 2 = 0.90 для условий нагружения пластин и γ 1 = 0.975, γ 2 = 0.95 – для цилиндрических оболочек. В качестве уравнений состояния материалов облицовки, линзы и корпуса заряда обычно используются ударные адиабаты в форме http://technomag.edu.ru/doc/334177.html 5 ( ) ( ) [ ] 1 0 − ρ ρ = ρ = b n b A p p , где ρ 0 – начальная плотность А, n b – эмпирические константы [1], причем при n b ≠ 1 это уравнение называется ударной адиабатой Тэта, а при n b = 1 оно вырождается в линейную баротропную зависимость, в которой A b = K 0 , где K 0 – модуль объемного сжатия. Для продуктов детонации в качестве уравнения состояния рекомендуется использовать изоэнтропу в форме степенного двучлена Коэффициенты изоэнтропы В, Си) определяются по параметрам в точке Чепмена – Жуге [1 – 3] ( ) , 1 ) ( 1 CJ CJ CJ CJ e p p k n ρ − γ − γ − + = ( ) 1 1 CJ CJ CJ CJ γ − − ⋅ ρ ρ − γ − = n n e p B n , γ ρ ρ − = CJ CJ CJ n B p C , где , CJ p , CJ ρ CJ e – давление, плотность и энергия на фронте ДВ, k – показатель адиабаты для конденсированных ВВ k ≈ 3): , 1 2 BB BB CJ + ρ = k D p , 1 BB CJ ρ + = ρ k k 1 Для расчета детонации используется искусственный прием, согласно которому наиболее полное выполнение условий на фронте ДВ достигается подбором изэнтропы ПД и уравнения состояния непрореагировавшего ВВ [ ] 1 ) ( 3 BB BB − ρ ρ = C p , в котором коэффициент С вв определяется из условия 75 0 4 Область разностной сетки, заполненная воздухом, в представленной постановке рассчитывается приближенно. Для этого в качестве изэнтропы условного воздуха используется изэнтропа ПД. Граничные условия для рассматриваемых задач в рамках идеализированных расчетных схем (рис. 1), задаются на участках поверхностей взаимодействующих сред. Так, например, поверхности КО и корпуса, контактирующих с воздухом, и поверхность продуктов детонации, истекающих из корпуса, допустимо рассматривать как свободные от действия внешних поверхностных сил. То есть, пренебрегая силами атмосферного давления для ПД), где n j – вектор единичной нормали. На поверхностях КО, линзы и корпуса, контактирующих с ПД, накладываются ограничения на скорости движения индивидуальных точек в соответствии с условием непроницаемости 77-30569/334177 , №02 февраль 2012 г. http://technomag.edu.ru 6 ; ) ПД ( (KO) i i i i n v n v = ; ) ПД ( ) линза ( i i i i n v n v = , ) ПД ( ) корп ( i i i i n v n v = а также на напряженное состояние, реализующееся в этих точках в соответствии с третьим законом Ньютона ; ) ПД ( ) KO ( ni j ij p n = σ ; ) ПД ( ) линза ( ni j ij p n = σ ) ПД ( ) корп ( ni j ij p n = σ При формулировке граничных условий на оси симметрии (ось 0z) необходимо учитывать, что при r = 0 частицы среды движутся только в осевом направлении (v r = 0), а осевые ускорения этих частиц должны быть ограничены. Из уравнений движения (2) – (3) следует, что это может быть реализовано только при отсутствии касательных напряжений на оси симметрии (σ rz = 0). Давление, плотность и энергия во фронте ДВ задаются равными параметрам в точке Чепмена – Жуге. Выделение контактных разрывов типа ПД – корпус или ПД – КО осуществлялось методом концентраций [2, 3], разработанного В.В. Кореньковым. Для этого в дополнение к основным характеристикам течения (ρ, p, v r , v z , D σrr , D σrz , D σzz ) для неоднородной системы ПД (ВВ) – КО (корпус) определяют массовую концентрацию веществ следующим образом = = корпуса. для воздуха, или ПД ВВ, для 0 Здесь M 1 – масса ВВ, ПД или воздуха в ячейке, а M – масса ячейки. Таким образом, для однородных ячеек, содержащих только первую компоненту w = 1, для однородных ячеек со второй компонентой w = 0 и, наконец, для смешанных ячеек 0 < w < 1. Локализация контактных границ осуществляется из анализа текущего распределения концентрации взаимодействующих веществ, которое определяется законом сохранения концентрации 0 = ρ α + ∂ ρ ∂ + ∂ ρ ∂ + ∂ ρ ∂ r v w z v w r v w t w r z r ) ( ) ( ) ( (10) Давление в смешанной ячейке вычисляется из условия аддитивности удельных объемов и равенства давлений содержащихся в ней компонент ( ) ρ = ρ ρ + ρ − = ρ = = = = ), ( ) ( ; 1 1 1 1 0 0 0 1 w w w w p p w w (11) где ρ, ρ w = 0 , ρ w = 1 – текущая плотность смешанной ячейки и плотности содержащихся в ней компонент p 0 (ρ w = 0 ) – ударная адиабата материала корпуса p 1 (ρ w = 1 ) – изэнтропа ПД или уравнение состояния ВВ. По значениям поля концентраций определяется вид контактной границы и на его основе рассчитывается поток из смешанной ячейки. На разрывах типа http://technomag.edu.ru/doc/334177.html 7 ПД – КО или ПД – корпус в процессе решения обеспечиваются условия равенства нормальных составляющих скоростей и напряжений. Выделения границ ПД – воздух не производится. Начальные условия конкретной задачи задаются распределением параметров ρ, p, w, ив поле течения. Компоненты напряжений принимаются равными нулю. Для получения более полной информации о процессах схлопывания КО или метания корпуса последние обычно маркируют реперными точками — маркерами или трассе- рами (позиция 8 на рис, в которых дополнительно вычисляют параметры текущего состояния среды P m = {ρ, p, v r , v z …} по формуле ( ) [ ] ( ) z 4 1 ∆ ∆ ⋅ = ∑ = r A k P P k k m m (12) Здесь k = 1 2, 3, 4; A k – площади, определяемые текущим положением маркера m [2 – 5], а P m (k) – параметры состояния среды в ячейках (k), окружающих маркер. Система уравнений (1) – (11) интегрировалась модифицированным методом крупных частиц, подробно изложенным в следующем разделе (разд. 2), на программном комплексе, позволяющем проводить вычисления во всей области интегрирования на неподвижной (эйлеровой) сетке без предварительного анализа течения. |