ЖЭТФ, 2020, том 157, вып. 3, стр. 406-427
© 2020
ДОПЛЕРОВСКИЕ ГЕТЕРОДИННЫЕ ИЗМЕРЕНИЯ И
МОДЕЛИРОВАНИЕ ВЫБРОСА ЧАСТИЦ С ПОВЕРХНОСТИ
УДАРНО-НАГРУЖЕННЫХ ОБРАЗЦОВ
А. В. Андрияшa, С. А. Дьячковa,b, В. В. Жаховскийa,b, Д. А. Калашниковc,
А. Н. Кондратьевa*, С. Е. Куратовa, А. Л. Михайловc, Д. Б. Рогозкинa,d,
А. В. Федоровc, С. А. Финюшинc, Е. А. Чудаковc
a Всероссийский научно-исследовательский институт автоматики им. Н. Л. Духова
127055, Москва, Россия
b Объединенный институт высоких температур Российской академии наук
125412, Москва, Россия
c Российский федеральный ядерный центр —
Всероссийский научно-исследовательский институт экспериментальной физики
607188, Саров, Нижегородская обл., Россия
d Национальный исследовательский ядерный университет «МИФИ»
115409, Москва, Россия
Поступила в редакцию 4 июня 2019 г.,
после переработки 13 ноября 2019 г.
Принята к публикации 14 ноября 2019 г.
Представлены результаты доплеровских гетеродинных измерений выброса пылевых частиц при ударно-
волновом нагружении металлических образцов. Эксперименты проводились на образцах из олова и свин-
ца заданной толщины и качества обработки поверхности. Для условий, близких к экспериментальным,
методом SPH (smooth particle hydrodynamics) выполнено прямое численное моделирование процесса вы-
броса массы с поверхности ударно-нагруженных образцов. Определены поверхностная плотность и на-
чальное распределение объемной плотности выброшенной массы по скорости. На основе этих результатов
рассчитана временная зависимость профиля объемной плотности при расширении образовавшегося пы-
левого облака в воздушную среду. По спектральным данным доплеровских гетеродинных измерений с
помощью подхода, основанного на транспортном уравнении для корреляционной функции рассеянного
поля, восстановлены основные параметры распределения пылевых частиц по скорости, поверхностная
плотность пыления и др. При соответствующем выборе разброса пылевых частиц по размерам удается
описать наблюдаемую в экспериментах временную динамику спектров, обусловленную торможением пы-
ли в воздухе. Восстановленные из эксперимента значения массы выброшенного вещества согласуются с
результатами SPH-моделирования.
DOI: 10.31857/S0044451020030037
ми, фольгами Эсея, теневой протоно- и рентгено-
графией и др.) [1-16], стали широко использовать-
1. ВВЕДЕНИЕ
ся доплеровские гетеродинные измерения скорости
(photon doppler velocimetry, PDV) [17] (см. также
В последнее десятилетие для изучения динами-
[9, 10, 12-16, 18-31]). Этот метод позволяет одновре-
ческих процессов, сопровождающих выход ударной
менно регистрировать скорости движения несколь-
волны на поверхность материалов, наряду с тра-
ких или даже большого числа объектов, попадаю-
диционными методами измерений физических ха-
щих в область зондирования. Именно такая ситуа-
рактеристик продуктов разрушения (пьезодатчика-
ция реализуется при выходе ударной волны на по-
верхность образца, когда происходит кумулятивный
* E-mail: an.kondratev@physics.msu.ru
406
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
выброс твердотельных или жидко-капельных мик-
стик пыления с результатами численного моделиро-
рочастиц из неоднородностей поверхности.
вания выброса массы методом сглаженных частиц
При измерении скорости гетеродинным методом
(smoothed particle hydrodynamics, SPH [32]). Допле-
используется смешение двух волн, различающихся
ровские гетеродинные измерения выброса пыли про-
по частоте [17, 18]. Одна, опорная, волна представ-
водились при ударно-волновом нагружении образ-
ляет собой сигнал на исходной частоте, другая, от-
цов Pb и Sn заданной толщины и качества обработ-
раженная от движущегося объекта, имеет частоту,
ки поверхности (шероховатость Rz = 7 мкм и Rz =
смещенную из-за эффекта Доплера. Интерферен-
= 20 мкм для образцов соответственно Pb и Sn). В
ция волн в детекторе приводит к временным би-
зависимости от схемы ударного воздействия давле-
ениям интенсивности. По частоте биений, которая
ние в ударной волне варьировалось от 13.5 до 23 ГПа
совпадает с доплеровским сдвигом ω = 2(v/c)ω0
для Pb и составляло 16 ГПа для Sn. В эксперимен-
(v — скорость движения объекта, c — скорость света,
тах продемонстрирована качественная зависимость
ω0 — частота зондирующего излучения), определя-
процесса пыления от давления в ударной волне.
ется значение скорости в заданный момент времени.
Для условий, близких к экспериментальным,
При отражении от многочастичной системы, та-
проведено прямое численное SPH-моделирование
кой как пылевой выброс с поверхности ударно-на-
процесса выброса массы с поверхности ударно-на-
груженного образца, гетеродинный сигнал состоит
груженных образцов. Выполнены расчеты поверх-
из большого числа гармоник с различными амп-
ностной плотности и начального распределения объ-
литудами и фазами. Результаты обработки вре-
емной плотности выброшенной массы по скорости в
менной развертки биений с помощью дискретного
зависимости от параметров профиля неровности по-
фурье-преобразования собираются в так называе-
верхности образца. На основе этих результатов рас-
мые спектрограммы [17-19] — двумерные распреде-
считана временная эволюция распределения объем-
ления на плоскости «частота (или скорость) — вре-
ной плотности пылевых частиц при их торможении
мя», где каждому значению фурье-спектра ставится
в воздухе. Для теоретического расчета доплеровско-
в соответствие точка с яркостью, пропорциональной
го спектра обратнорассеянного сигнала разработан
спектральной амплитуде. Такое представление ре-
численный метод решения транспортного уравне-
зультатов PDV-измерений является наглядным, но в
ния для корреляционной функции поля, обобщаю-
случае отражения от многочастичной системы недо-
щий методы работ [28,30] с учетом перераспределе-
статочно информативным для количественного ана-
ния частиц по скоростям и координатам в процессе
лиза данных эксперимента. В ситуации, когда от-
торможения. Входящие в решение оптическая тол-
раженный сигнал формируется многократно рассе-
щина и параметры начального профиля объемной
янными волнами, спектральное распределение уже
плотности пылевого выброса, распределение частиц
нельзя однозначно связать со скоростями отдельных
по размерам определены из условия наилучшего со-
частиц и требуется более сложная статистическая
гласия теоретически рассчитанного спектра с зави-
интерпретация экспериментальных данных.
сящими от времени спектральными данными допле-
Для анализа данных доплеровских гетеродин-
ровских гетеродинных измерений. Восстановленные
ных измерений в работе [28] был предложен под-
таким способом из эксперимента значения поверх-
ход, основанный на использовании транспортного
ностной плотности пыления находятся в согласии с
уравнения для корреляционной функции рассеянно-
результатами SPH-моделирования.
го поля. Такой подход позволяет связать доплеровс-
кий спектр обратнорассеянного сигнала с исходны-
ми параметрами пылевых выбросов и открывает но-
2. ЭКСПЕРИМЕНТ
вые возможности для восстановления по спектраль-
ным PDV-данным различных физических характе-
При ударном воздействии на образец происхо-
ристик пыления ударно-нагруженных образцов (по-
дит разрушение его поверхности, которое сопро-
верхностной плотности пыления, распределения пы-
вождается кумулятивным выбросом микрочастиц.
левых частиц по скоростям и размерам и др.).
В результате перед образцом возникает баллисти-
Ниже представлены результаты применения под-
чески расширяющееся облако металлической пы-
хода работы [28] к анализу новых эксперименталь-
ли. В настоящей работе для исследования пыления
ных данных по выбросу частиц при ударном воздей-
ударно-нагруженных образцов использовались до-
ствии на образцы свинца (Pb) и олова (Sn) и срав-
плеровские гетеродинные измерения скорости дви-
нения восстановленных из эксперимента характери-
жения фрагментов разрушения поверхности. Схе-
407
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
z, мкм
a
7
4
2
6
0
-2
5
-4
50
100
150
200
250
300
350
x, мкм
4
z, мкм
3
15
б
2
10
5
1
0
Рис. 1. Схема эксперимента: 1 — капсюль электродетона-
тора; 2 — пластическое ВВ; 3 — плоский ударник (в экспе-
-5
риментах с образцами Pb — стальная пластина диаметром
30 мм, толщиной 2 мм (13.5 ГПа) и 1.7 мм (23 ГПа)); 4
-10
исследуемый образец (пластина Pb диаметром 30 мм и
0
50
100
150
200
250
300
350
толщиной 1.7 мм, пластина Sn диаметром 40 мм и тол-
x, мкм
щиной 2.4 мм); 5 — движущиеся фрагменты разрушения
Рис. 2. Профили поверхности образцов Pb (а) и Sn (б)
поверхности; 6 — зондирующий луч; 7 — коллиматор
ма микросборки, с помощью которой выполнялись
PDV-сигнала на этом интервале описываются мас-
экспериментальные исследования, проиллюстриро-
сивом из 5 · 106 точек. Для получения частотного
вана на рис. 1. Ударная волна в исследуемом образ-
спектра PDV-сигнала временной интервал был раз-
це формируется при воздействии на него ударником
делен на равные части — «окна». В каждом окне
различной толщины из стали или непосредственно
проводилось дискретное преобразование Фурье. В
взрывчатым веществом (ВВ).
результате был получен набор PDV-спектров, отве-
чающих различным моментам времени разлета пы-
В экспериментах исследовалось пыление образ-
левого облака. Результаты измерений скорости дви-
цов из свинца и олова в воздушной среде при нор-
жения продуктов разрушения ударно-нагруженных
мальных условиях. Пластины свинца имели толщи-
образцов показаны в виде спектрограмм на рис. 3,
ну 1.7 мм, шероховатость поверхности Rz = 7 мкм,
на которых изображен квадрат амплитуды спектра
пластины олова — толщину 2.4 мм, шероховатость
гетеродинных биений |J(ω)|2.
поверхности Rz = 20 мкм. Примеры профилограмм
поверхности исследуемых образцов показаны на
Значение давления в ударной волне при выходе
рис. 2. Исследование пыления ударно-нагруженных
на свободную поверхность напрямую связано с ма-
образцов проводилось PDV-методом в составе ком-
териалом образца, его толщиной, а также со схемой
плекса гетеродин-интерферометр
[25-27]. Грани-
возбуждения ударной волны (ударниками различ-
цы полосы пропускания фотодиода и осциллографа
ной толщины или непосредственно ВВ). Использо-
позволяли регистрировать скорости до 10 км/с. В
валась схема с инициированием ВВ в одной точке от
экспериментах измерялась временная зависимость
электродетонатора. Диаметр ВВ в экспериментах с
фототока, возникающего при интерференции опор-
образцами Pb составлял 15 мм, толщина варьиро-
ного и обратнорассеянного от пылевого облака вол-
валась от 10 до 20 мм. В эксперименте с образцом
новых полей. Полный временной интервал измере-
Sn диаметр и толщина ВВ равнялись 10 мм. Рассто-
ний составлял 100 мкс. Временные биения фототока
яние от коллиматоров до исследуемой поверхности
408
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
v, км/с
v, км/с
|J(
)|2
|J(
)|2
1.6
1.6
2.5
1.00
2.5
1.50
а
1.4
б
1.4
1.2
1.2
1.0
0.40
1.0
0.84
2.0
2.0
0.8
0.8
0.6
0.6
0.16
0.47
4.5
5.0
5.5
6.0
6.5
4.5
5.0
5.5
6.0
6.5
1.5
1.5
t, мкс
t, мкс
0.063
0.27
1.0
1.0
0.025
0.15
0.5
0.5
0
10
20
30
40
50
0
10
20
30
40
50
Время, мкс
Время, мкс
v, км/с
v, км/с
2
|J(
)|2
|J( )|
2.5
2.5
2.5
1.00
2.5
в
2.0
г
2.0
0.050
1.5
1.5
0.40
2.0
1.0
2.0
1.0
0.024
0.16
5.5
6.0
6.5
7.0
7.5
9.0
9.5
10.0
10.5
0.011
t, мкс
t, мкс
1.5
0.063
1.5
0.0063
0.025
0.0025
1.0
1.0
0
10
20
30
40
0
10
20
30
40
50
Время, мкс
Время, мкс
Рис. 3. (В цвете онлайн) Спектрограммы движения продуктов разрушения ударно-нагруженных образцов Pb (канал 1
(а), канал 2, давление в ударной волне 13.5 ГПа (б), давление в ударной волне 23 ГПа (в)) и Sn (давление в ударной волне
16 ГПа (г)). На вставках показаны участки спектрограмм, отвечающие первым 1.5 мкс после выхода ударной волны на
поверхность
Таблица
1. Параметры образцов и результаты
составляло 40 мм. Давление в ударной волне опре-
PDV-измерений
делялось по скорости свободной поверхности [33],
значение которой оценивалось по положению пика
№ Веще- Rz, λ, P, vfs, vmax,
в фурье-спектре в момент времени, соответствую-
щий выходу ударной волны на поверхность образца.
опыта ство мкм мкм ГПа км/с км/с
Для свинцовых пластин давление в двух опытах со
1
Pb
7
60
13.5
0.88
1.50
стальными ударниками толщиной 2 мм и 1.7 мм со-
2
Pb
7
60
23.0
1.32
2.30
ставляло соответственно 13.5 ГПа и 23 ГПа. В опы-
3
Sn
20
80
16.0
1.23
2.10
те с образцом Sn ударная волна возбуждалась непо-
средственно ВВ, и давление в ударной волне состав-
ляло 16 ГПа. При таких давлениях можно считать,
что образцы после нагружения оставались в твер-
характеристики поверхности исследуемых образ-
дом состоянии.
цов: средняя глубина неровностей Rz и среднее
расстояние между их вершинами — «длина волны»
Согласно результатам, представленным на
λ. Погрешность измерения скорости vfs методом
рис. 3, в экспериментах зарегистрированы сигналы
PDV в экспериментах не превышала 20 м/с.
как от выбрасываемых с поверхности пылевых
частиц, так и от свободной поверхности. Значения
На рис. 3а,б представлены результаты регист-
скорости свободной поверхности (vfs) и оценки
рации скорости свободной поверхности и облака час-
максимальной скорости движения частиц (vmax)
тиц в опыте 1 (P = 13.5 ГПа) для двух каналов реги-
для трех описанных выше опытов приведены в
страции, соответствующих зондированию точек на
табл. 1, в которой также представлены соответ-
поверхности, разнесенных на 5 мм. Плотность по-
ствующие значения давления в ударной волне и
тока частиц от участка поверхности, соответствую-
409
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
v, км/с
щего каналу 2, оказывается ниже плотности пото-
2.0
ка в канале 1. Поэтому в канале 2 отчетливо виден
сигнал от свободной поверхности на сравнительно
1.8
большом временном участке в начале записи реги-
1.6
страции. На спектрограмме опыта 2 (P = 23 ГПа),
представленной на рис. 3в, сигнал от свободной по-
1.4
2
верхности виден в течение первой микросекунды с
1.2
момента выхода ударной волны на поверхность об-
1.0
разца. Затем плотное облако пылевых частиц экра-
нирует свободную поверхность, и сигнал от нее ис-
0.8
1
чезает.
0.6
Максимальные значения скорости пыли оцени-
6
8
10
12
14
16
18
20
вались по фурье-спектру в первые моменты после
Время, мкс
выхода ударной волны на поверхность и составля-
ли в первых двух опытах соответственно около 1.5
Рис. 4. Наложение временных зависимостей скорости час-
и 2.3 км/с. Отношения vmax/vfs в этих опытах ока-
тиц v(tm) для различных v0 на спектрограмму опыта 1 (об-
зываются близкими, что можно объяснить одинако-
разец Pb, 13.5 ГПа, канал 1). Размер частиц d = 2.3 мкм.
Скорости свободной поверхности vfs = 0.88 км/с (линия
вым качеством обработки поверхности.
1) и ударной волны в воздухе vsh = 1.15 км/с (линия 2)
В опыте 3 с образцом олова, нагруженным до
давления P = 16 ГПа, слабый сигнал от свободной
газа, d — диаметр частицы, ρ0 — плотность материа-
поверхности виден менее 1 мкс (рис. 3г). Максималь-
ла частицы, ρg — плотность газа. Коэффициент со-
ная скорость пылевых частиц составляла прибли-
противления среды CD зависит от числа Рейнольдса
зительно 2.1 км/с. Регистрируемую максимальную
Re = ρg(v-vg)d/μg, где μg — динамическая вязкость
скорость пыли при подлете к коллиматору можно
среды [36, 37]. Наложение кривых v(tm) для частиц
оценить как 1.3 км/с.
размером d = 2.3 мкм с различной начальной скоро-
Тот факт, что на спектрограммах при малых вре-
стью v0 на верхнюю границу спектрограммы опыта
менах после выхода ударной волны на поверхность
1 (канал 1) проиллюстрировано на рис. 4. Такой спо-
образца удается различить сигнал от свободной по-
соб позволяет, подбирая значения d, определить ха-
верхности (см. вставки на рис. 3), можно объяс-
рактерный размер частиц в пылевом выбросе. Для
нить неоднородным в направлении, параллельном
опытов 1, 2 и 3 размеры частиц можно оценить со-
поверхности, распределением пыли. Выбросы пы-
ответственно как 2.3 мкм, 1.5 мкм и около 3 мкм.
левых частиц в виде плоских струй перемежаются
Как следует из показанных на рис. 3 спектро-
относительно прозрачными областями [34, 35]. По-
грамм, после нескольких первых микросекунд сиг-
этому, если ось зондирования образует малый угол
нал от лидирующих микрочастиц облака кумуля-
с направлением выбросов, пылевые частицы в пер-
тивной пыли из-за их сильного торможения исче-
вые моменты не экранируют сигнал от свободной
зает. В результате для зондирующего излучения от-
поверхности (наиболее отчетливо этот эффект на-
крывается часть пылевого облака, состоящая из час-
блюдается в опыте с образцом Pb, 13.5 ГПа, канал
тиц, движущихся медленнее ударной волны в возду-
2). Со временем поток пылевых частиц закрывает
хе и поэтому испытывающих относительно слабое
поверхность, и сигнал от нее пропадает.
торможение. Сравнение спектрограмм, отвечающих
По наклону начальных участков спектрограмм
опытам 1 и 2, свидетельствует о том, что с ростом
можно сделать оценку характерных размеров от-
давления в ударной волне увеличивается доля вы-
дельных частиц, возникающих в процессе разруше-
брошенных с поверхности пылевых частиц, движу-
ния поверхности при ударно-волновом нагружении.
щихся с относительно большими скоростями (со ско-
Такая оценка основана на зависимости торможения
ростями, превышающими скорость ударной волны в
частиц в воздухе от их размера вследствие гидроди-
воздухе, vsh).
намического сопротивления [36],
dv
3 CDρg
=-
Δv2,
(1)
3. МОДЕЛЬ СРЕДЫ
dtm
4
ρ0d
где tm — время движения, Δv = v - vg — скорость
Для анализа данных PDV-измерений и восста-
движения частицы относительно газа, vg — скорость
новления по этим данным параметров пылевых вы-
410
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
бросов необходимо иметь представления об основ-
1.5-2 раза. Экспериментальные данные часто ап-
ных закономерностях пыления ударно-нагружен-
проксимируют экспоненциальной зависимостью [12]
ных образцов. Сведения о них можно получить
F (x) = β exp [(x - 1)] ,
(4)
на основе экспериментальных данных, полученных
различными контактными методами (с помощью
где параметр β принимает в зависимости от условий
пьезодатчиков, фольг Эсея и др.) [2-6,9-16], анали-
эксперимента значения порядка 5-20 [12, 29]. Наи-
тических расчетов [38, 39], а также данных числен-
большие отклонения от зависимости (4) наблюдает-
ного моделирования методами молекулярной дина-
ся при скоростях, близких к максимальной скорости
мики и SPH-методом [34,35,40-42]. В отличие от мо-
движения частиц и к скорости свободной поверхнос-
лекулярной динамики гидродинамическое SPH-мо-
ти [10,35,40].
делирование позволяет определить при заданных
Соотношение
(2) справедливо при движении
условиях эксперимента (давление в ударной волне,
фрагментов разрушения в вакууме. При разлете пы-
профиль шероховатости поверхности образцов) аб-
левого облака в газовой среде это соотношение при-
солютные значения таких характеристик пылеобра-
менимо только в первые моменты времени после вы-
зования, как максимальная скорость пылевых час-
хода ударной волны на поверхность, когда можно
тиц, распределение объемной плотности по скоро-
пренебречь эффектом торможения частиц.
сти, поверхностная плотность пыления.
Ниже приведены результаты SPH-расчетов
В первом приближении можно считать, что пы-
основных характеристик пылевых выбросов для
левое облако представляет собой плоский слой час-
конкретных условий экспериментов, описанных в
тиц, движущихся без столкновений с определенным
предыдущем разделе.
разбросом по скоростям [3-6,10,12]. Временная эво-
люция распределения объемной плотности вещества
3.1. Моделирование начальной стадии
в пылевом облаке описывается автомодельной зави-
симостью [6]
процесса пыления
(
)
Пределы текучести исследуемых материалов
ρa
z
ρ(z, tm) =
F
,
(2)
(свинца и олова) на три порядка меньше, чем
vfstm
vfstm
амплитуда реализуемых в эксперименте ударных
где ρa — поверхностная плотность пыления, vfs
волн. Поэтому будем предполагать, что имеет место
скорость свободной поверхности, координата z
только пластическое течение материала образца.
расстояние от ее исходного положения. Время tm от-
Движение вещества при численном моделировании
считывается от момента выхода ударной волны на
будем описывать уравнениями сохранения массы,
поверхность образца. Функция F (x) предполагается
импульса и энергии:
нормированной на единицу.
Наряду с объемной плотностью (2) во многих ра-
= -ρ∇ · U,
(5)
ботах [2-6,10,12] используется так называемая ди-
dt
dU
1
намическая поверхностная плотность
=-
∇P,
(6)
dt
ρ
de
P
ρa(z, tm) = dzρ(z, tm),
(3)
=-
∇ · U,
(7)
dt
ρ
z
где ρ — плотность вещества, U — его скорость, P
которая является функцией отношения z/vfstm
давление, e — внутренняя энергия единицы объема.
(или v/vfs). При z
= vfstm величина (3) дает
Так как процесс пыления характеризуется высоки-
значение поверхностной плотности ρa выброшенной
ми скоростями движения вещества, силы вязкости и
массы.
поверхностного натяжения не оказывают заметного
Прямые измерения (с помощью пьезодатчиков,
влияния на течение. Ниже эти силы не учитывают-
фольг Эсея и др.) и численное моделирование пыле-
ся.
ния металлических образцов под воздействием удар-
Чтобы замкнуть систему уравнений
(5)-(7),
ной волны позволяют приближенно установить вид
необходимо задать уравнение состояния P = P (ρ, e).
функции F(x) [3-6,10,12,35,40]. Согласно этим дан-
В данной работе воспользуемся уравнением состоя-
ным, плотность вещества в пылевом облаке рез-
ния в форме Ми - Грюнайзена [43]:
ко падает с ростом скорости частиц, которая мо-
жет превышать скорость свободной поверхности в
P - Pr = γρ(e - er),
(8)
411
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
Таблица 2. Характеристики свинца и олова, ис-
пользуемые при SPH-моделировании
Характеристика материала Свинец Олово
Начальная плотность ρ0, кг/м3
11350
7206
Коэффициент Грюнайзена γ
1.7
1.7
Коэффициент a
1.26
1.55
Рис. 5. Постановка задачи для расчетов методом SPH: об-
Скорость звука c, км/с
2.58
2.37
разец толщиной L с гофрированной поверхностью ударя-
ется о жесткую стенку со скоростью up. В результате удара
возникает ударная волна, которая распространяется к гоф-
где Pr(ρ) и er(ρ) — опорные кривые давления и внут-
рированной поверхности
ренней энергии, γ — коэффициент Грюнайзена. При
плотности вещества выше начальной (ρ0/ρ < 1) со-
стояние вещества описывается ударной адиабатой
кая структура поверхности в реальных образцах мо-
вида us = c + aup, где us — скорость ударной волны,
жет варьироваться. Поэтому рассмотрены два вари-
up — скорость вещества, c — скорость звука, a — чис-
анта моделирования: величина Rz фиксирована, а λ
ленный коэффициент. В этом случае опорные кри-
распределена случайно; величина λ фиксирована, а
вые давления и энергии можно записать в виде [44]
Rz распределена случайно.
Скорость up удара образца о стенку подбиралась
1 - ρ0
Pr 10
Pr = ρ0c2
,
er =
(9)
так, чтобы скорость свободной поверхности vfs
[1 - a(1 - ρ0)]2
ρ0
2
2up, и амплитуда ударной волны соответствова-
При плотности вещества ниже начальной (ρ0/ρ > 1),
ли экспериментальным значениям. Поскольку ско-
т. е. при растяжении вещества, предполагается, что
рости струй определяются углами раствора гофри-
давление линейно связано с удельным объемом:
ровки, диапазон последних был подобран так, что-
бы скорости струй соответствовали наблюдаемым в
Pr = ρ0c2(1 - ρ0).
(10)
эксперименте на начальной стадии (длительностью
менее 1 мкс) процесса пыления.
Значения констант, входящих в соотношения
Для решения уравнений (5)-(7) в настоящей
(8)-(10), приведены в табл. 2.
работе используется бессеточный лагранжев метод
У образцов, используемых в экспериментах,
SPH (см., например, работу [32] и Приложение A).
структура неоднородностей поверхности пред-
Результаты расчета процесса формирования
ставляет собой систему параллельных бороздок.
струй через 0.5 мкс после удара образца о жесткую
Параметры шероховатости характеризуются сред-
стенку показаны на рис. 6. Как следует из расче-
ним расстоянием между бороздками λ и их
тов, скорости струй определяются углом раствора
глубиной Rz . При выходе ударной волны на та-
гофрировки, т. е. соотношением между Rz и λ,
кую поверхность из углублений выходит ансамбль
а количество выброшенной массы существенно
струй с разными скоростями [1-16]. Это отчетливо
зависит только от λ.
видно на PDV-спектрограммах. Причиной разброса
Скорости истечения струй vsp — практически ли-
скоростей струй является разброс значений λ и
нейные функции угла раствора ϕ (см. рис. 7),
Rz и, соответственно, угла раствора гофрировки
ϕ = 2arctg(λ/2Rz). Распределение выброшенной
массы по скорости получается усреднением по
vsp = [1 + α(180 - ϕ)]vfs,
(11)
ансамблю струй и отличается от соответствующего
распределения для одиночной струи.
где численный коэффициент α слабо зависит от дав-
Постановка задачи для SPH-моделирования про-
ления в ударной волне и материала образца. Линей-
иллюстрирована на рис. 5. Используются периоди-
ная аппроксимация дает хорошие результаты для
ческие граничные условия по осям x и y, образцу за-
небольшой амплитуды ударной волны. При более
дается скорость up в направлении жесткой стенки.
сильном воздействии становятся заметными нели-
Угол гофрировки варьируется от некоторого зна-
нейные эффекты.
чения ϕmin до ϕmax в соответствии с наблюдае-
Рассчитанные профили объемной плотности вы-
мым диапазоном скоростей струй. Микроскопичес-
брошенной массы показаны на рис. 6б. При v/vfs >
412
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
> 1 их можно аппроксимировать зависимостью
а
vsp = 1494 м/с
t = 0.5 мкс
xsp - 1 + C
ρ(x) = B + A
+ δρ(x),
(12)
120
x-1+C
. Слагаемое
где x = v/vfs, xsp = vsp/vfs, x < xsp
vsp = 1369 м/с
δρ(x) в (12) описывает утолщение в головной части
130
струи,
(
)
vsp = 1241 м/с
xsp - x
(xsp - x)2
δρ(x) = M
exp
-
,
(13)
140
Δx
2S2Δx2
где Δx — ширина утолщения конца струи на оси ско-
vsp = 1106 м/с
рости. Заметим, что во всех расчетах с высокой точ-
150
ностью можно предполагать, что скорость и коорди-
ната связаны линейно: v = z/tm.
104
Для каждого расчета при заданном значении уг-
б
ла ϕ (или Rz и λ) подбирались коэффициенты A, B,
C, M, S. Далее по ним строились соответствующие
полиномиальные интерполяции A(ϕ), B(ϕ), C(ϕ),
103
M (ϕ), S(ϕ) внутри выбранного диапазона углов
150
140
130
120
[ϕmin, ϕmax]. В результате профиль распределения
массы ρ(x, ϕ) может быть получен для любого уг-
ла ϕ внутри этого диапазона. Далее аналитическая
102
функция ρ(x, ϕ) (12) усреднялась в предположении,
что разброс значений угла ϕ подчиняется гауссову
1.0
1.2
1.4
1.6
1.8
распределению с центром при ϕ0 = (ϕmax + ϕmin)/2
v/vfs
с дисперсией σ = 0.17(ϕmax - ϕmin). В результате
такого усреднения было определено распределение
Рис. 6. (В цвете онлайн) а) Результаты расчетов пыле-
массы в ансамбле струй.
ния свинца для опыта 1 (13.5 ГПа и up = 0.4 км/с) для
Результаты расчетов усредненных таким обра-
углов от ϕ от 120 до 150 на момент времени 0.5 мкс.
зом профилей объемной плотности ρ(x) массы для
б) Результаты SPH-расчетов профилей объемной плотно-
опытов 1, 2 и 3 представлены на рис. 8. Видно, что
сти пыления ρ(v/vfs). Аналитическая аппроксимация (12)
утолщения на концах струй исчезают с усреднени-
показана сплошными черными линиями
ем по ансамблю. Также стоит отметить, что суще-
ственной разницы между средними профилями, по-
v /vspfs -1
лученными при вариации Rz или λ, не наблюдается,
1.2
поэтому на практике для моделирования ансамбля
может использоваться вариация какой-либо одной
1.0
из этих величин.
Выбор диапазона углов для построения средних
0.8
профилей ансамбля струй заметно влияет на наблю-
даемую максимальную скорость частиц. На примере
0.6
олова показано, что усреднение в диапазоне углов от
100 до 140 дает профиль с наибольшей скоростью;
0.4
при этом головные части струй с утолщением дают
Pb, 13.5 ГПа
заметный вклад в среднюю плотность. Это приво-
0.2
Pb, 23 ГПа
дит к небольшому изгибу кривой зависимости ρ(x).
Sn, 16 ГПа
Усреднение в диапазонах 110-150 и 120-160 пол-
0
20
40
60
80
100
ностью подавляет этот эффект.
180
-
Интеграл от усредненной объемной плотности
ρ(x) дает значение динамической поверхностной
Рис. 7. Зависимости скоростей струй от угла гофрировки
плотности пыления (см. соотношение (3)). Результа-
ты интегрирования ρ(x) для трех рассматриваемых
413
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
а
а
101
103
100
102
10-1
1
2
3
1
2
3
101
б
б
101
103
100
102
10-1
1
2
3
1
2
3
101
в
101
в
103
100
102
10-1
1
2
3
1
2
3
10-2
101
1.0
1.2
1.4
1.6
1.8
2.0
1.0
1.2
1.4
1.6
1.8
2.0
v/vfs
v/vfs
Рис. 9. (В цвете онлайн) Результаты расчета профи-
Рис. 8. (В цвете онлайн) Результаты расчета профиля объ-
ля динамической поверхностной плотности ρa(x) ансам-
емной плотности ρ(x) ансамбля струй для опытов 1 (а), 2
бля струй для опытов 1 (а), 2 (б) и 3 (в) при фикси-
(б) и 3 (в) при фиксированных значениях λ (сплошные
рованном значении λ. Диапазон углов усреднения и по-
кривые) и Rz (штриховые). Диапазоны углов усреднения
верхностная плотность выброшенной массы для опыта 1
для опыта 1 (Pb, 13.5 ГПа): кривые 1 140-160; 2
(Pb, 13.5 ГПа): кривые 1 140-160; 2 130-160 ;
130-150 ; 3 120 -150; опыта 2 (Pb, 23 ГПа): 1
3 120-150, ρa(v/vfs = 1) = 13.7 мг/см2; опыта
140-160 ; 2 120-160 ; 3 90-150 ; опыта 3 (Sn,
2 (Pb, 23 ГПа): 1 140-160; 2 120-160; 3
16 ГПа): 1 120-160; 2 110 -150; 3 100-140 .
90-150 , ρa(v/vfs = 1) = 12.1 мг/см2; опыта 3 (Sn,
Расчеты относятся к моменту времени 0.5 мкс после вы-
16 ГПа): 1 120-160 ; 2 110 -150; 3 100-140 ,
хода ударной волны на поверхность
ρa(v/vfs = 1) = 10.8 мг/см2. Расчеты относятся к мо-
менту времени 0.5 мкс после выхода ударной волны на
поверхность
опытов представлены на рис. 9. Численные значе-
ния полной поверхностной плотности ρa выброшен-
0.3-0.5 мкс [34,35] распадаются струи выброшен-
ной массы приведены в подписи к рисунку.
ной из образца массы. Определение размеров пы-
Полученные начальные профили распределения
левых частиц выходит за рамки SPH-моделирова-
плотности далее используются для моделирования
ния [45]. В относительных единицах профиль рас-
торможения пылевого облака в газовой среде, а так-
пределения по размерам удается получить с помо-
же для расчета доплеровского спектра сигнала об-
щью MД [34, 35]. Наиболее полную информацию о
ратного рассеяния.
распределении частиц по размерам на сегодняшний
день дают прямые оптические измерения (гологра-
3.2. Распределение частиц по размерам
фия [46, 47], рассеяние Ми [48, 49]).
Чтобы связать плотностные характеристики пы-
Для заданных экспериментальных условий (ма-
левых выбросов с их оптическими характеристи-
териала образца, профиля неровностей поверхнос-
ками, необходимо иметь информацию о размерах
ти, давления в ударной волне и др.) характерные
пылевых частиц, на которые на временах tm
размеры фрагментов разрушения (микрочастиц ку-
414
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
мулятивной природы, осколков поверхности и др.)
1.0
могут зависеть от их относительной скорости v/vfs
[46]. Большинство имеющихся данных о размерах
относится к наиболее быстрым, лидирующим час-
0.8
тицам, скорость движения которых заметно пре-
вышает скорость движения свободной поверхности
0.6
v/vfs > 1.1 [7, 11, 46, 48, 50]. Средние размеры та-
ких частиц обычно не превышают нескольких мик-
0.4
рометров [46, 48]. О размерах более крупных фраг-
0.5
1.0
2.0
5.0
ментов разрушения, которые движутся вблизи сво-
2 1/2
(d )
,
мкм
бодной поверхности и могут давать [10,12] заметный
вклад в полную массу пылевого выброса, имеются
Рис. 10. Отношение поверхностной плотности пыления
к транспортной оптической толщине пылевого слоя как
только данные качественного характера [7,46].
функция среднеквадратичного размера частиц. Результа-
Распределение частиц по размерам [7, 11, 46,
ты расчета для Pb (нижняя кривая) и Sn (верхняя кри-
50] можно аппроксимировать распределением Юн-
вая) с лог-нормальным распределением частиц по разме-
ге [51]
рам (σ = 0.5). Длина волны излучения 1.55 мкм
n(d) = Ad,
где d — размер частицы, dmin < d < dmax, или
лог-нормальной зависимостью [47, 49]
(
)
сеяния частицей размером d, g — средний косинус
1
1
d
n(d, σ) =
exp
-
ln2
,
(14)
угла рассеяния, m = ρ0πd3/6 — средняя масса од-
2πσ2d
2σ2
dm
ной пылевой частицы, ρ0 — плотность материала.
где dm — медианный размер, σ — ширина распреде-
Интеграл от (15) по всей толщине пылевого слоя
ления.
позволяет связать его поверхностную плотность с
Степенную зависимость часто используют для
транспортной оптической толщиной
описания «хвоста» распределения при сравнитель-
но больших d [38,47]. Показатель α может варьиро-
σs(1 - g)
ваться в достаточно широких пределах [35,47,50] от
τtr =
dz σtr (z) =
ρa.
(16)
m
α = 1.5 [35] до α = 5.6 [50].
zfs
Лог-нормальное распределение с шириной σ =
= 0.2-1.0 хорошо описывает экспериментальные
Отношение ρatr зависит от характерного масшта-
данные по распределению частиц по размерам при
ба распределения частиц по размерам. Поскольку
v/vfs > 1.1 [47,49]. Средние размеры пылевых час-
регистрируемый в PDV-экспериментах сигнал про-
тиц могут варьироваться в интервале 1-10 мкм в
порционален сечению обратного рассеяния, которое
зависимости от условий эксперимента и характе-
для металлических частиц размером порядка или
ристик образца, в частности от угла гофрировки
больше длины волны пропорционально геометриче-
неровностей поверхности [49].
скому сечению πd2/4, отношение ρatr удобно пред-
В предположении, что расстояние между части-
ставить как функцию среднеквадратичного размера
цами превышает их размер и рассеяние зондирую-
частиц. Для лог-нормального распределения с ти-
щей волны на отдельных частицах происходит неза-
пичной для эксперимента шириной σ = 0.5 справед-
висимо, задавая распределение частиц по размерам,
лива оценка ρatr ≃ ρ0(d2)1/2.
можно связать оптические характеристики пылево-
Связь поверхностной плотности пыления с
го облака с его объемной и поверхностной плотно-
транспортной оптической толщиной τtr для об-
стями [30]. Рассеивающую способность пылевого об-
разцов Pb и Sn проиллюстрирована на рис. 10.
лака удобно характеризовать его транспортным ко-
Для расчета транспортного сечения рассеяния
эффициентом рассеяния [51-53], который выража-
использовались численный код [54] и оптические
ется через объемную плотность пыления соотноше-
постоянные веществ [55].
нием
Как следует из представленных на рис. 10 ре-
σs(1 - g)
зультатов, пылевые выбросы в экспериментах, опи-
σtr(z, tm) =
ρ(z, tm),
(15)
m
санных в разд. 2, характеризуются большой оп-
где σs(1 - g) — транспортное сечение рассеяния,
тической толщиной. Для размеров частиц d
усредненное по размерам частиц, σs — сечение рас-
1.5-3 мкм, которые следуют из оценки по накло-
415
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
v/vfs
тями — сильно затормозившиеся с начальной ско-
1.6
ростью v0 > vsh (vsh — скорость ударной волны в
воздухе [36]) и испытывающие слабое торможение
1.5
с начальной скоростью v < vsh. Поскольку каждая
1.4
точка на кривой v(z) отвечает определенному зна-
чению начальной скорости частиц, ей соответствует
1.3
определенное значение доли объемной плотности.
1.2
Наглядно это проиллюстрировано на рис. 12.
Для того чтобы далее рассчитать эволюцию объ-
1.1
емной плотности с учетом разброса частиц по раз-
1.0
мерам, показанные на рис. 12 профили объемной
плотности нужно просуммировать с весовыми мно-
1.0
1.1
1.2
1.3
1.4
1.5
1.6
жителями, отвечающими распределению размеров
z/zfs
частиц. Если полученный результат дополнитель-
но проинтегрировать по координате z, то полу-
Рис. 11. (В цвете онлайн) Семейство кривых v(z) для час-
чим распределение объемной плотности по скорос-
тиц Sn размером 2 мкм. Моменты времени tm = 0.5-4 мкс
ти, ρ(v, tm). Аналогично, если выполнить интегри-
(от верхней к нижней кривой), скорость свободной поверх-
рование по скорости при фиксированном значении
ности vfs = 1.23 км/с
z, приходим к распределению объемной плотнос-
ти ρ(z, tm).
ну PDV-спектрограмм опытов 1-3 (см. рис. 4), и
Временная эволюция транспортного коэффици-
результатов SPH-моделирования для поверхностной
ента рассеяния рассчитывается аналогично, если за-
плотности ρa, ожидаемые значения транспортной
менить в выражении для распределения объемной
оптической толщины лежат в диапазоне τtr = 5-15.
плотности массу отдельной частицы md = ρ0πd3/6
на транспортное сечение σs(1-g). Из-за перераспре-
деления частиц по скоростям и координатам в про-
3.3. Эволюция распределения объемной
цессе торможения равенство (15) нарушается, одна-
плотности в газовой среде
ко интегральное соотношение (16) остается в силе.
Торможение пылевых частиц в воздухе может
Результаты соответствующих расчетов для экс-
приводить к существенному изменению исходно-
перимента с образцом Sn (16 ГПа) проиллюстриро-
го распределения выброшенной массы по скоро-
ваны на рис. 13. Сравнение рис. 13а и 13б показыва-
сти движения v и координате z. Наиболее быстрые
ет, что те части распределений, которые определя-
пылевые частицы, согласно соотношению (1), тор-
ются вкладом частиц, движущихся медленнее удар-
мозятся на ранних стадиях разлета пылевого об-
ной волны в воздухе, слабо меняются со временем и
лака. Это хорошо видно на спектрограммах (см.
практически совпадают между собой. Вклад частиц,
рис. 3). Поэтому из-за процесса торможения проис-
первоначально движущихся быстрее ударной волны
ходит подавление высокоскоростной части распре-
в воздухе, с ростом времени уменьшается в распре-
деления объемной плотности по скорости.
делении по координатам медленнее, чем в распреде-
Для расчета эволюции распределения пылевых
лении по скоростям. Как следует из рис. 13а, с рос-
частиц по скоростям и координатам воспользуемся
том tm в распределении по скоростям становится за-
подходом, примененным к аналогичной проблеме в
метным вклад частиц, затормозившихся до v < vfs.
работе [56] (см. также Приложение B). Предполага-
ется, что торможение частиц в газовой среде подчи-
няется уравнению (1).
4. РАСЧЕТ ДОПЛЕРОВСКОГО СПЕКТРА
Для описания временной зависимости распреде-
ОБРАТНОРАССЕЯННЫХ ВОЛН
ления объемной плотности сначала для частиц за-
данного размера рассчитывается семейство траек-
Используемый в настоящей работе подход к ана-
торий, отвечающих различным значениям началь-
лизу спектров обратнорассеянных волн был пред-
ной скорости. На основе этих результатов на плос-
ложен в работе [28] и применялся затем в [29, 30]
кости vz формируется набор кривых v(z) для каж-
для восстановления параметров пылевых выбросов
дого момента времени tm (рис. 11). На заданной глу-
из экспериментальных PDV-данных. Он основан на
бине z оказываются частицы с различными скорос-
использовании связи между спектральной мощно-
416
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
а
10
1.6
v/vfs
а
1.4
1.2
1
1.0
10
0.1
1
0.1
0.01
1.0
1.2
z/zfs
1.4
1.0
1.1
1.2
1.3
1.4
1.5
1.6
1.6
1.6
v/vfs
v/vfs
1.4
б
1.2
1.0
10
б
10
1
1
0.1
1.0
1.2
0.1
z/zfs
1.4
1.6
1.6
0.01
v/vfs
1.4
в
1.2
1.0
1.1
1.2
1.3
1.4
1.5
1.6
1.7
1.0
z/zfs
10
1
Рис. 13. (В цвете онлайн) Эволюция распределений объем-
0.1
ной плотности ρ(v, tm) (а) и ρ(z, tm) (б) частиц Sn. Резуль-
1.0
таты расчетов для моментов времени tm = 0.5, 4, 8 мкс (от
1.2
правой к левой кривой). Начальное распределение массы
z/z
1.4
по скорости — расчет SPH-методом для образца Sn при
fs
1.6
16 ГПа, ϕ = 140 . Среднеквадратичный размер частиц
1.6
(d2)1/2 = 2 мкм
v/vfs
1.4
г
1.2
1.0
стью |J(ω)|2/Δt гетеродинных биений (Δt — вре-
менное «окно» фурье-преобразования) и временной
10
корреляционной функцией обратнорассеянных волн
1
[28, 30]. Как показано в работах [28, 30], спектраль-
0.1
/Δt гетеродинных биений с
ная мощность |J(ω)|2
1.0
точностью до общего множителя совпадает с допле-
1.2
ровским спектром |E(ω + ω0)|2/Δt обратнорассеян-
z/zfs
1.4
1.6
ных волн, где E(ω) — фурье-образ рассеянного поля,
— частота зондирующего излучения. Усреднен-
ω0
Рис. 12. Эволюция профиля объемной плотности частиц
ный по флуктуациям доплеровский спектр
заданного размера. Результаты расчета для частиц Sn раз-
мером 2 мкм для моментов времени tm = 0 (а), tm = 1 мкс
I(ω) = 〈|E(ω + ω0)|2/Δt〉
(б), tm = 3 мкс (в), tm = 5 мкс (г). Начальный профиль —
результат SPH-моделирования для ϕ = 140. Скорость
напрямую связан с корреляционной функцией об-
свободной поверхности vfs = 1.23 км/с
ратнорассеянного поля [51]:
I(ω) = dt eiωtI(r, n = -n0, t),
(17)
417
3
ЖЭТФ, вып. 3
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
где
пылевого облака зависимость этих коэффициентов
от z меняется со временем. Значение времени раз-
лета tm входит в эти коэффициенты как параметр.
I(r, n, t) = d3r exp(-ik0n · r +0t) ×
Зондируемый объем пылевого выброса ограни-
× 〈E (r+r/2, tm+t/2) E (r-r/2, tm-t/2)〉 ,
(18)
чен поперечными размерами поля зрения PDV-сис-
темы. Предполагается, что поперечные размеры об-
вектор r определяет положение детектора, n0 и n
ласти зондирования намного больше, чем среднее
направления распространения соответственно зон-
значение длины свободного пробега фотонов в обла-
дирующей и рассеянной волн, k0 = ω0/c — волно-
ке. В продольном направлении зондируемый объем
вое число. В предположении, что пылевые частицы
простирается от области наиболее быстрых частиц
расположены случайным образом и их число внутри
пылевого облака до свободной поверхности. Про-
зондируемого объема велико, зависящая от времени
дольный размер зондируемого объема меняется со
«интенсивность» I(r, n, t) удовлетворяет транспорт-
временем tm, но общее число пылевых частиц в этом
ному уравнению [51].
объеме, их суммарная масса и, соответственно, оп-
В отличие от работы [28], где расчеты приводи-
тическая толщина предполагаются неизменными.
лись простейшим двухпотоковым методом, в настоя-
щей работе при решении уравнения для I(r, n, t) ис-
В случае, когда пыление ударно-нагруженного
пользуется так называемое транспортное приближе-
образца происходит в вакууме, пылевое облако рас-
ние [52,53,57]. Транспортное приближение оправда-
ширяется в соответствии с автомодельным распре-
но, если дифференциальное сечение рассеяния мож-
делением (2). В этом приближении скорость частиц
но представить как суперпозицию изотропной ча-
однозначно связана с координатой z: v(z) = z/tm.
сти и дифракционного пика. Именно такая ситуа-
Входящий в интегральное слагаемое уравнения (19)
ция складывается при рассеянии на металлических
фактор рассеяния принимает вид
частицах. В транспортном приближении дифферен-
циальное сечение заменяется на сумму изотропного
〈σtr exp(-ik0(μ - μ)vt) =
сечения и δ-пика в направлении вперед. Интеграль-
(
)
z
ное по углам сечение изотропной части полагает-
= σtr(z)exp
-ik0(μ - μ)
t
(20)
tm
ся равным транспортному сечению рассеяния. Тогда
исходное транспортное уравнение сводится к урав-
При расширении пылевого облака в газовую сре-
нению для изотропного рассеяния и в рассматрива-
ду автомодельное приближение можно применять
емом случае плоской геометрии имеет вид
только в первые моменты времени после выхода
(
)
ударной волны на поверхность образца, пока тормо-
μ
+ σtr + κ I(z,μ,t) =
∂z
жением частиц в среде можно пренебречь. В усло-
1
виях торможения, когда нарушается однозначная
1
связь между скоростью v и координатой z (см.
=
〈σtr exp[-ik0(μ - μ)vt]〉I(z, μ, t),
(19)
2
рис. 11), входящий в уравнение (19) фактор рассея-
-1
ния можно представить в виде многогруппового раз-
где σtr и κ — соответственно транспортный коэффи-
ложения
циент рассеяния и коэффициент поглощения, угло-
вые скобки означают усреднение по скоростям пыле-
〈σtr exp[-ik0(μ - μ)vt] =
вых частиц. Координата z отсчитывается от началь-
ного положения свободной поверхности, μ — косинус
= σ(i)tr (z)exp[-ik0(μ - μ)vi(z)t],
(21)
угла между направлением распространения рассе-
i
янной волны и осью z. Граничное условие к этому
где суммирование проводится по группам частиц,
уравнению для падающей на среду из z = + в
для которых справедлива однозначная связь между
направлении μ = -1 плоской волны имеет вид
скоростью и координатой и зависимости vi(z) явля-
1
ются монотонными функциями z. Расчет входящих
I(z = +∞, μ, t) =
δ(μ + 1).
2π
в уравнение (21) величин основан на результатах
разд. 3.3 (см. также Приложение B).
Входящие в уравнение (19) коэффициенты σtr и
κ, характеризующие оптические свойства пылевого
Процедура решения уравнения (19) с фактором
облака, являются функциями z. Из-за расширения
рассеяния (21) изложена в Приложении C.
418
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
б
а
I, отн. ед.
I, отн. ед.
I,
отн. ед.
I, отн. ед.
8
8
8
8
6
6
6
6
4
4
4
4
2
2
2
2
0
0
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
I, отн. ед.
/ fs
I, отн. ед.
/ fs
I,
отн. ед.
/ fs
I, отн. ед.
/ fs
8
8
8
8
6
6
6
6
4
4
4
4
2
2
2
2
0
0
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
/ fs
/ fs
/ fs
I, отн. ед.
/ fs
I, отн. ед.
I, отн. ед.
I,
отн. ед.
8
8
8
8
6
6
6
6
4
4
4
4
2
2
2
2
0
0
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
I, отн. ед.
/ fs
I, отн. ед.
/ fs
I,
отн. ед.
/ fs
I, отн. ед.
/ fs
8
8
8
8
6
6
6
6
4
4
4
4
2
2
2
2
0
0
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
/ fs
/ fs
/ fs
/ fs
I, отн. ед.
I, отн. ед.
I,
отн. ед.
I, отн. ед.
8
8
8
8
6
6
6
6
4
4
4
4
2
2
2
2
0
0
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
/ fs
/ fs
/ fs
/ fs
Рис. 14. (В цвете онлайн) Эволюция экспериментально измеренного спектрального распределения 〈|J(ω)|2 (синие ли-
нии) и теоретически рассчитанного доплеровского спектра I(ω) (красные линии) в опытах с образцом Pb (13.5 ГПа, канал
1 (а) и канал 2 (б)) в интервале tm = 5.5-14.5 мкс. Момент времени выхода ударной волны на поверхность tm = 5 мкс.
Временное окно Δt = 125 нс, интервал усреднения Δtm = 0.5 мкс. Доплеровский сдвиг сигнала от свободной поверх-
ности ωfs = (2vfs/c)ω0. Расчеты приведены для значений транспортной оптической толщины τtr = 9 (а) и τtr = 5 (б),
среднеквадратичного размера частиц (d2)1/2 = 2.1 мкм и SPH-профиля объемной плотности с ϕ = 140
5. СРАВНЕНИЕ С ЭКСПЕРИМЕНТОМ.
результате получается набор спектров, отвечающих
ОБСУЖДЕНИЕ РЕЗУЛЬТАТОВ
различным моментам времени tm движения пылево-
го облака. Согласно работам [28,30], для подавления
Стандартная обработка экспериментальных дан-
флуктуаций в спектрах можно провести дополни-
ных PDV-измерений заключается в расчете фурье-
тельное усреднение спектральной мощности |J(ω)|2
спектров J(ω) гетеродинных биений фототока во
гетеродинных биений. Временная шкала разбивает-
временном «окне» Δt от tm - Δt/2 до tm + Δt/2. В
ся на квазистационарные интервалы Δtm. Значение
419
3*
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
I, отн. ед.
I,
отн. ед.
Δtm выбирается из условия, чтобы основные пара-
10
10
метры спектрального профиля |J(ω)|2 (ширина, по-
8
8
ложение максимума) оставались на таком времен-
6
6
ном отрезке практически неизменными. На каждом
4
4
квазистационарном интервале рассчитывается сред-
2
2
няя спектральная мощность
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
1
/ fs
/ fs
〈|J(ω)|2 =
|J(ω)|2n,
(22)
N
I, отн. ед.
I,
отн. ед.
n=1
10
10
8
8
где |J(ω)|n — амплитуда спектра в момент време-
6
6
ни tmn), N — число измеренных спектров на выбран-
4
4
ном временном отрезке Δtm. Результаты усреднения
2
2
оказываются достаточно устойчивыми относитель-
0
0
но вариации длительности окна Δt (расчеты прово-
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
дились для Δt = 0.05-0.25 мкс).
/ fs
/ fs
I, отн. ед.
I,
отн. ед.
В предположении о случайном расположении
10
10
пылевых частиц усредненная по флуктуациям вели-
8
8
чина 〈|J(ω)|2 должна быть пропорциональна опре-
6
6
деленному выше доплеровскому спектру обратного
4
4
рассеяния I(ω) (см. разд. 4, (17)). Результаты срав-
2
2
нения экспериментальных данных для 〈|J(ω)|2 с
0
0
рассчитанным теоретически спектром I(ω) показа-
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
ны на рис. 14-16. Эти спектральные распределения
/ fs
/ fs
I, отн. ед.
I,
отн. ед.
относятся к различным моментам времени разлета
10
10
пылевого облака.
8
8
Теоретический расчет доплеровских спектров
6
6
проводился по следующей схеме. По результатам
4
4
SPH-моделирования для объемной плотности вы-
2
2
брошенной массы (см. рис. 8) задавалось началь-
0
0
0.8
1.0
1.2
1.4
1.6
1.8
ное распределение пылевых частиц по скорости (т. е.
0.8
1.0
1.2
1.4
1.6
1.8
/ fs
/ fs
входящая в (2) функция F (v/vfs)). Затем по этому
распределению, в предположении лог-нормального
Рис. 15. (В цвете онлайн) Эволюция экспериментально из-
разброса частиц по размерам, рассчитывалась вре-
меренного спектрального распределения 〈|J(ω)|2 (синие
менная эволюция пространственной зависимости
линии) и теоретически рассчитанного доплеровского спек-
фактора рассеяния 〈σtr exp[-ik0(μ-μ)vt], т.е. вхо-
тра I(ω) (красные линии) в опыте с образцом Pb (23 ГПа)
дящие в (21) пространственные профили скорости
в интервале tm = 6.5-13.5 мкс. Момент времени выхо-
частиц и транспортного коэффициента рассеяния
да ударной волны на поверхность tm = 6 мкс. Временное
для различных моментов времени tm. Далее, по
окно Δt = 125 нс, интервал усреднения Δtm = 0.5 мкс.
определенной таким образом функции (21) численно
Расчеты приведены для транспортной оптической тол-
решалось транспортное уравнение и рассчитывал-
щины τtr
= 14, среднеквадратичного размера частиц
(d2)1/2
= 1.2 мкм и SPH-профиля объемной плотности
ся доплеровский спектр обратнорассеянного сигна-
с ϕ = 140
ла (Приложение C). Параметры, от которых зави-
сит наблюдаемый спектр, такие как транспортная
оптическая толщина пылевого слоя τtr и входящий
в SPH-моделирование средний угол гофрировки по-
нием относительно свободной поверхности [30]. Ско-
верхности образца ϕ, подбирались из условия наи-
рость свободной поверхности находилась из спек-
лучшего согласия рассчитанного и эксперименталь-
трального распределения при малых tm (см. спек-
но измеренного спектров в первые моменты времени
трограммы на рис. 3). Параметры распределения
(порядка 0.5 мкс) после выхода ударной волны на
частиц по размерам (среднеквадратичный размер
поверхность образца. Значение τtr определяется по
(d2)1/2 и ширина σ) подбирались из условия наилуч-
соотношению между шириной спектра и его смеще-
шего согласия рассчитанного и измеренного в экспе-
420
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
Таблица 3. Восстановленные параметры пылевых
Применимость приведенных выше теоретичес-
выбросов
ких расчетов доплеровских спектров имеет ограни-
. Это ограничение обусловлено
чение по времени tm
(d2)1/2,
ρa,
ρSPHa,
№ опыта τ
ϕ
двумя причинами.
tr
мкм мг/см2 мг/см2
Одна причина связана с пренебрежением в на-
1(1)
9
2.1
18
ших расчетах столкновительным механизмом тор-
13.8/9.5
140
можения в облаке пылевых частиц. Как следует из
1(2)
5
2.1
10
представленных на рис. 11 кривых v(z), с ростом
2
14
1.2
13
12.1/10.1 140
времени движения tm область толщин z, где од-
3
10
2.0
13
10.8/15.5 140
новременно находятся частицы с различными ско-
ростями (область неоднозначной зависимости v(z)),
расширяется. Одновременно увеличивается «опти-
ческая» толщина этой области относительно про-
рименте спектральных распределений для различ-
цесса столкновения пылевых частиц, движущихся
ных значений tm. На классе начальных распреде-
с различными скоростями, между собой. Когда со-
лений частиц по скорости, рассчитанных SPH-мето-
ответствующая оптическая толщина становится по-
дом, погрешность определения значений τtr и (d2)1/2
рядка единицы, пренебрегать вероятностью столк-
можно оценить как не превышающую 10-15 %.
новения частиц между собой уже нельзя. Процесс
Восстановленные таким образом значения пара-
столкновений приводит к выравниванию скоростей
метров пылевых выбросов приведены в табл. 3. Зна-
частиц в этой области и, соответственно, к эффек-
чение поверхностной плотности ρa найдено по опре-
тивному торможению более быстрых частиц и уско-
деленным из сравнения с экспериментом значениям
рению медленных. Время tm, при котором столк-
τtr и (d2)1/2 с помощью результатов, представлен-
новительный механизм торможения частиц в пыле-
ных на рис. 10. Для сопоставления в табл. 3 при-
вом слое становится заметным, зависит от оптичес-
ведены результаты расчета ρa SPH-методом (пер-
кой толщины облака пылевых частиц, движущихся
вое значение — при фиксированном λ, второе — при
быстрее ударной волны в воздухе. По оценкам для
фиксированном Rz).
первых двух опытов этот механизм становится важ-
ным при временах, превышающих 15 мкс (13.5 ГПа)
Некоторое расхождение между восстановленны-
и 5 мкс (23 ГПа) после выхода ударной волны на по-
ми из эксперимента значениями ρa и результатами
верхность образца.
SPH-моделирования можно в первую очередь свя-
зать с различиями между реальным профилем по-
Вторая причина обусловлена нарушением с рос-
верхности образцов и используемым в моделирова-
том tm приближения плоской ударной волны в
нии. На разброс параметров шероховатости поверх-
воздухе. В частности, в эксперименте с образцом
ности указывает, в частности, и различие данных,
Sn диаметр ударного воздействия составлял всего
полученных в опыте 1 для одного и того же образца
10 мм. Поэтому на временах, превышающих 3 мкс
с двух участков поверхности, разнесенных всего на
после выброса частиц, ударная волна в воздухе ста-
5 мм.
новится сферической, и ее скорость существенно па-
дает [58]. Расчет эволюции доплеровских спектров в
Для образца Sn оценка поверхностной плотности
этом предположении также показан на рис. 16. Как
пыления ρa, полученная из результатов наложения
следует из рис. 16, компактизация потока частиц
рассчитанных спектров на измеренные вPDV-экспе-
при скоростях, близких к скорости ударной волны в
рименте, несколько ниже значения ρa, определен-
воздухе, которая отражается на спектре в виде вы-
ного контактными методами [12]. Для образцов с
раженного пика при v/vfs 1.2, в этом случае не
близкими характеристиками поверхности в работе
наблюдается, и согласие с экспериментом прослежи-
[12] при давлении 27 ГПа были получены значения
вается до больших значений tm.
поверхностной плотности около 20 мг/см2. Имею-
щиеся различия в оценках ρa, возможно, обуслов-
В экспериментах с образцами Pb диаметр удар-
лены вкладом относительно крупных фрагментов,
ника составлял 30 мм, и приближение плоской удар-
движущихся вблизи свободной поверхности, а так-
ной волны в воздухе может нарушаться на суще-
же разницей в значениях давления в ударной волне
ственно больших временах (превышающих 15 мкс
и, соответственно, различным фазовым состоянием
(P = 13.5 ГПа) и 10 мкс (P = 23 ГПа) после выхода
выброшенных частиц вещества.
ударной волны на поверхность).
421
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
I, отн. ед.
I, отн. ед.
пыления (поверхностную плотность пыления, сред-
8
8
неквадратичный размер частиц, начальный про-
филь объемной плотности). Приведены результаты
6
6
экспериментов на образцах Sn и Pb заданной тол-
4
4
щины и качества обработки поверхности. Продемон-
2
2
стрирована зависимость процесса пыления от давле-
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
ния в ударной волне и других условий эксперимента.
/ fs
/ fs
Прямое численное моделирование процесса
I, отн. ед.
I, отн. ед.
выброса массы с поверхности ударно-нагруженных
8
8
образцов, выполненное бессеточным лагранжевым
6
6
SPH-методом, позволило определить поверхност-
4
4
ную плотность и начальный профиль объемной
2
2
плотности пыления. На основе этих результатов
0
0
рассчитана временная эволюция распределения
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
объемной плотности при расширении образовав-
/ fs
/ fs
I, отн. ед.
I, отн. ед.
шегося пылевого облака в газовую среду. Для
интерпретации усредненных спектральных данных
8
8
PDV-измерений использовано решение транспорт-
6
6
ного уравнения для второго момента рассеянного
4
4
волнового поля. Показано, что, подбирая значения
2
2
оптической толщины пылевого слоя, разброса
0
0
частиц по скоростям и размерам, удается добить-
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
ся хорошего согласия теоретических расчетов с
I, отн. ед.
/ fs
I, отн. ед.
/ fs
полученными из эксперимента спектральными про-
8
8
филями сигнала обратного рассеяния в широком
6
6
диапазоне времен tm. Такой подход дает возмож-
4
4
ность на основе данных доплеровских гетеродинных
измерений решать обратную задачу восстановления
2
2
основных характеристик возникающего при разру-
0
0
0.8
1.0
1.2
1.4
1.6
1.8
0.8
1.0
1.2
1.4
1.6
1.8
шении поверхности образца пылевого выброса.
/ fs
/ fs
Финансирование. Работа одного из авторов
Рис. 16. (В цвете онлайн)
Эволюция экспериментально
(В. В. Ж.) поддержана Российским научным фон-
измеренного спектрального распределения 〈|J(ω)|2 (си-
дом (проект № 19-19-00697).
ние линии) и теоретически рассчитанного доплеровско-
го спектра I(ω) (красные и черные линии для плоской
и сферической ударных волн, соответственно) в опыте
ПРИЛОЖЕНИЕ A
с образцом Sn (16 ГПа) в интервале tm = 10-17 мкс.
Момент времени выхода ударной волны на поверхность
Моделирование начальной стадии пыления
tm = 9.3 мкс. Временное окно Δt = 125 нс, интервал
контактным SPH-методом
усреднения Δtm = 0.5 мкс. Расчеты приведены для транс-
портной оптической толщины τtr = 10, среднеквадратич-
Метод сглаженных частиц
— SPH-метод
ного размера частиц (d2)1/2 = 2.0 мкм и SPH-профиля
представляет собой способ интерполяции некоторой
объемной плотности с ϕ = 140
функции f(r), значения которой заданы на дис-
кретном наборе из N точек пространства r1,...,rN .
В некоторой точке r пространства с использованием
6. ВЫВОДЫ
этих известных значений и сглаживающей функции
(ядра) W (r, h) можно построить аппроксимацию
В настоящей работе показано, что гетеродин-
значения функции
ные измерения доплеровского спектра сигнала об-
ратного рассеяния от ударно-нагруженных металли-
ческих образцов дают возможность количественно
f (r) =
fjΔVjW(|r - rj|, h),
(A.1)
определять основные физические характеристики
j=1
422
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
0.9
ство конечно-разностных SPH-аппроксимаций урав-
нений (5)-(7). В данной работе используется схема
0.8
с решением задачи Римана на контактной границе
между парой частиц [59], которая обеспечивает мо-
нотонность функций в окрестности разрыва. В этом
0.7
случае система уравнений (5)-(7) записывается в ви-
1 CPU: ячейка Вороного
де
0.6
Ячейки движутся
вместе с веществом
i
= 2ρi
mj (Ui - U∗ij) · ∇Wij,
(A.6)
0.5
dt
ρj
0
200
400
600
j=1
Количество итераций
dUi
2
Рис. 17. (В цвете онлайн) Отношение времени SPH-расче-
=-
mj P∗ij∇Wij,
(A.7)
dt
ρi
ρj
та к полному времени на один шаг моделирования (полез-
j=1
ная нагрузка) как функция количества итераций работы
алгоритма балансировки. Внизу показана пространствен-
dei
2
ная декомпозиция моделируемого образца на ячейки Во-
=
mj P∗ij(Ui - U∗ij) · ∇Wij,
(A.8)
dt
ρi
ρj
роного в ходе параллельного расчета. Частицы, принадле-
j=1
жащие одной ячейке, обрабатываются отдельным процес-
где U∗ij и P∗ij — значения контактной скорости и дав-
сором
ления для частиц i и j, которые получаются при
решении одномерной задачи Римана вдоль линии,
соединяющей их центры. Для расчета контактных
где ΔVj — соответствующий элемент объема, h
значений используется решение, приведенное в ра-
дистанция сглаживания. Cглаживающее ядро удов-
боте [60]. В качестве сглаживающей функции ис-
летворяет условиям
пользуется полином Wendland-C2 [61].
Для расчетов методом SPH с необходимой точно-
W (|r - R|, h)d3r = 1,
(A.2)
стью требуется большое число частиц (высокое про-
|r-R|<h
странственное разрешение). Используемый в этом
случае алгоритм параллельных вычислений должен
lim W (|r - R|, h) = δ(R),
(A.3)
быть оптимизирован по времени счета. В данной ра-
h→0
боте используется динамическое разбиение среды на
где δ(R) — дельта-функция. Удобство такого пред-
ячейки Вороного [62], за счет автоматического пере-
ставления заключается в том, что градиент функ-
распределения частиц между которыми достигается
ции f(r) выражается через производную сглажива-
баланс вычислительной нагрузки. Суть подхода по-
ющего ядра W (r, h), вычисленную для разных точек
яснена на рис. 17. Моделируемый образец состоит
пространства:
из большого числа частиц (несколько миллионов),
которые обрабатываются параллельно работающи-
ми процессорами. Каждый процессор обрабатывает
∇f(r) = fj ΔVj ∇W (|r - rj |, h),
(A.4)
только локальный набор данных о частицах (обычно
j=1
50-100 тысяч частиц на ядро процессора), которые
соответствуют одной ячейке Вороного в моделируе-
r-rj
мом пространстве. Диаграмма Вороного однознач-
∇W (|r - rj |, h) = W(|r - rj |, h)
(A.5)
|r - rj |
но определяется положениями центров ячеек, по-
этому изменение положения этих центров с целью
Таким образом, аппроксимация значения производ-
улучшить баланс вычислительной нагрузки приво-
ной функции в данной пространственной точке рас-
дит к изменению диаграммы и соответствующему
считывается с использованием ее ближайшего окру-
перераспределению частиц между рабочими процес-
жения на дистанции сглаживания путем простого
сорами. Таким образом параллельную декомпози-
суммирования.
цию удается быстро адаптировать к изменяющим-
ся условиям. Детали алгоритма подробно описаны в
Используя такое представление пространствен-
ной зависимости функции, можно построить множе-
работе [62].
423
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
(
)
Представление среды в виде набора SPH-частиц
Na
v(i)0
N (v, z, tm) =
F
×
имеет ряд преимуществ перед сеточными методами.
vfs
vfs
i
Эйлеровы сеточные схемы требуют отслеживания
(
)
dz
-1
контактных поверхностей и свободных границ, что
×
(v(i)0, tm)
δ v - v(v(i)0,tm) ,
(B.2)
dv0
является основным источником погрешности при
моделировании этими методами. Для лагранжевых
сеточных схем необходимы специальные процедуры
где суммирование выполняется по корням v(i)0 =
перестройки сетки в области формирования струи
= v(i)0(z,tm) уравнения z = z(v0,tm). Распределе-
при пылении, что приводит к значительной чис-
ние N(v, z, tm) на плоскости vz отлично от нуля на
ленной диффузии импульса и энергии в зоне пе-
кривой v = v(z) (см. рис. 11), которая получается
репланировки и, как следствие, к потере мелкомас-
исключением v0 из системы уравнений v = v(v0, tm)
штабных физических процессов, происходящих при
и z = z(v0,tm). Каждая точка на этой кривой отве-
струеобразовании. При моделировании бессеточным
чает определенному значению начальной скорости
SPH-методом не приходится преодолевать указан-
v0. Суммирование в (B.2) эквивалентно суммирова-
ные технические сложности, поэтому он оптимально
нию по группам частиц, отвечающих монотонным
пригоден для моделирования таких физических яв-
участкам в зависимости v(z).
лений, как струеобразование, вращательные и сдви-
Формулы (B.1), (B.2) описывают эволюцию рас-
говые течения сжимаемых сред, фрагментация рас-
пределения частиц заданного диаметра. Чтобы вы-
четной области на части вследствие эффектов раз-
числить распределение объемной плотности выбро-
рушения.
шенной массы, ρ(v, z, tm), с учетом разброса пыле-
вых частиц по размерам d, нужно умножить (B.1)
πd3/6
или (B.2) на массу отдельной частицы md = ρ0
ПРИЛОЖЕНИЕ B
и проинтегрировать по d с распределением n(d). По-
лученный результат можно выразить через поверх-
Распределение пылевых частиц при
ностную плотность пыления, воспользовавшись со-
торможении в газовой среде
отношением Na = ρa/md, где md — средняя масса
Из-за торможения пылевых частиц в воздухе ис-
частицы. В настоящей работе предполагается, что
ходное распределение выброшенной массы по ско-
распределение частиц по размерам, n(d), не зависит
рости и координате с течением времени изменяется.
от их начальной скорости (согласно работам [47,49],
Функция распределения частиц определенного раз-
это оправдано при v0/vfs > 1.1, т. е. для большей
мера d в заданный момент времени tm может быть
доли частиц, вовлеченных в процесс многократного
записана в виде
рассеяния зондирующего сигнала).
Na
( v0 )
Если интегрирование заменить дискретной сум-
N (v, z, tm) =
dv0 F
×
мой по размерам, то выражение для объемной плот-
vfs
vfs
ности ρ(v, z, tm) сведется к суммированию по груп-
× δ(z - z(v0,tm))δ (v - v(v0,tm)), (B.1)
пам частиц, каждая из которых отвечает монотон-
где Na — число частиц, вылетающих с единицы
ному участку зависимости v(z) для частиц заданно-
площади, (1/vfs) F (v0/vfs) — нормированное на-
го размера. В расчетах, результаты которых приве-
чальное распределение пылевых частиц по скорос-
дены в настоящей работе, мы использовали 16-груп-
ти (см. соотношение (2)), v(v0, tm) и z(v0, tm) =
повое приближение. Распределение n(d) моделиро-
tm
=
dt v(v0, t) — скорость и координата частицы
валось ансамблем частиц четырех различных раз-
0
с начальной скоростью v0 в момент времени tm.
меров, для каждого из которых на кривой v(z)
Зависимость v(v0, tm) рассчитывается с помощью
могут быть выделены четыре монотонных участка
численного интегрирования уравнения (1) c учетом
(см. рис. 11). Использовалась аппроксимация лог-
изменения закона торможения в области за фрон-
нормального распределения (14) с шириной σ = 0.5
том ударной волны в воздухе [36,37]. Распределение
набором из четырех дискретных значений размера
(B.1) получается из решения одномерного кинети-
d/dm = 0.8, 1.3, 1.8, 2.3 с соответствующими весо-
ческого уравнения, если в нем учесть силу тормо-
выми множителями 0.50, 0.32, 0.13, 0.05. На резуль-
жения и пренебречь столкновениями частиц между
таты расчета объемной плотности и транспортного
собой.
коэффициента рассеяния дополнительное увеличе-
Если выполнить интегрирование по v0, то выра-
ние числа дискретных значений d/dm практически
жение (B.1) можно представить в виде
не влияет.
424
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
ПРИЛОЖЕНИЕ C
Полная оптическая толщина пылевого слоя равна
τtot = ζ(z = vfstm). Входящая в уравнение (C.3) ве-
Решение транспортного уравнения при
личина Φ(0)i(ζ, t) соответствует парциальной плотно-
многогрупповом описании пылевого слоя
сти от источника в виде падающей на среду плоской
монохроматической волны,
Для решения транспортного уравнения (19) в
случае рассеяния волн на ансамбле частиц, который
Φ(0)i(ζ, t) = exp[-ζ - ik0vi(ζ)t] .
(C.6)
характеризуется неоднозначной зависимостью ско-
рости частиц от координаты z, воспользуемся мно-
Ядра hik(ζ, ζ, t) в системе интегральных уравнений
гогрупповым представлением фактора рассеяния
(C.3) определяются выражениями
1
〈σtr exp[-ik0(μ - μ)vt]
Λ σ(k)tr(ζ)
(ζ)
hik(ζ, ζ, t) =
η
×
2
σtr(ζ)
|μ|
μ
1
в виде разложения (21) на парциальные вклады. По-
(
)
ζ
сле подстановки (21) исходное транспортное уравне-
× exp
-
- ik0μ[vk(ζ) - vi(ζ)]t
,
(C.7)
ние (19) преобразуется к виду
μ
(
)
где предполагается, что альбедо однократного рас-
μ
+ σtr(z) + κ(z) I(z,μ,t) =
сеяния Λ = σtr /(σtr + κ) одинаково для всех групп
∂z
частиц и не зависит от координат. Значение Λ опре-
1
=
σ(k)tr(z)exp[-ik0μvk(z)tk(z, t), (C.1)
деляется только материалом образца и длиной вол-
4π
k
ны излучения. Расчеты по теории Ми [54] для длины
волны λ = 1.55 мкм, оптических постоянных из [55]
где суммирование проводится по группам, каждая
и размеров частиц 1-5 мкм дают для Λ значения
из которых соответствует монотонному участку за-
0.89 для Pb и 0.83 для Sn.
висимости скорости vk(z) частиц заданного диамет-
Система интегральных уравнений (C.3) являет-
ра от z и введена парциальная плотность
ся многогрупповым обобщением известного уравне-
ния Милна [53] на случай движущихся частиц. В
1
одногрупповом приближении (т. е. для случая одно-
Φk(z, t) = 2π dμ exp[ik0μvk(z)t]I(z, μ, t). (C.2)
значной функции v(z)) уравнение (C.1) рассматри-
-1
валось ранее [30]. В настоящей работе используется
16-групповое приближение (ансамбль частиц четы-
Входящие в (C.1) величины σtr(z), κ(z), σ(k)tr(z) и
рех различных размеров, для каждого из которых
vk(z) зависят от времени движения tm как от пара-
на кривой v(z) могут быть выделены четыре моно-
метра.
тонных участка (см. Приложение B)).
Преобразуя уравнение (C.1) к интегральному ви-
Численное интегрирование системы уравнений
ду и выполняя в нем интегрирование по угловой
(C.3) проводится на сетке значений ζn и μm, являю-
переменной μ, получим систему уравнений относи-
щихся корнями полиномов Лежандра. Для каждого
тельно парциальной плотности Φi:
значения времени t, входящего в доплеровские фа-
зовые множители в (C.6), (C.7), решается система
соответствующих алгебраических уравнений и, та-
Φi(ζ, t) =
hik(ζ, ζ, tk(ζ, t) +
ким образом, формируется массив значений плот-
k
0
ности Φk(ζn, t).
+ Φ(0)i(ζ, t), (C.3)
Интенсивность обратнорассеянного сигнала
I(μ = 1, ζ = 0, t) вычисляется по формуле
где вместо переменной z введена оптическая глуби-
на
Λ
σ(k)tr(ζ)
I(t) =
×
4π
σtr(ζ)
ζ = dz[σtr(z) + κ(z)].
(C.4)
k
0
z
× exp[-ζ - ik0vk(ζ)tk(ζ, t). (C.8)
Входящий в (C.4) транспортный коэффициент рас-
Доплеровский спектр сигнала обратного рассеяния,
сеяния σtr (z) определяется суммой парциальных
I(ω), вычисляется далее с помощью дискретного
вкладов:
σtr(z) =
σ(k)tr(z).
(C.5)
преобразования Фурье по t от выражения (C.8).
k
425
А. В. Андрияш, С. А. Дьячков, В. В. Жаховский и др.
ЖЭТФ, том 157, вып. 3, 2020
Результат (C.8) относится к случаю обратного
3.
W. T. Buttler, M. B. Zellner, R. T. Olson et al., J.
рассеяния волн от свободного слоя (т. е. предполага-
Appl. Phys. 101, 063547 (2007).
ется, что отражение волн от свободной поверхности
4.
M. B. Zellner, M. Grover, J. E. Hammerberg et al.,
отсутствует).
J. Appl. Phys. 102, 013522 (2007).
Обобщение (C.8) на случай отражающей свобод-
5.
M. B. Zellner, W. Vogan-McNeil, J. E. Hammerberg
ной поверхности может быть получено с помощью
et al., J. Appl. Phys. 103, 123502 (2008).
соответствующего граничного условия (см. [30]) и
представлено в виде разложения решения по после-
6.
M. B. Zellner and W. T. Buttler, Appl. Phys. Lett.
довательным отражениям от границы. Это разло-
93, 114102 (2008).
жение удается просуммировать для диффузно-отра-
7.
В. А. Огородников, А. Л. Михайлов, В. В. Бурцев
жающей поверхности [30]. Наиболее простой резуль-
и др., ЖЭТФ 136, 615 (2009).
тат получается в предположении, что коэффициент
8.
V. V. Igonin, G. B. Krasovsky, S. E. Kuratov et al.,
отражения поверхности равен единице и отраже-
Phys. Scripta 142, 014019 (2010).
ние носит зеркальный характер [28,63]. В этом слу-
9.
W. T. Buttler, D. M. Oró, D. L. Preston et al., J.
чае интенсивность обратного рассеяния можно вы-
Fluid Mech. 703, 60 (2012).
разить через решение для свободного слоя удвоен-
ной оптической толщины. Значение интенсивности
10.
Y. Chen, H. Hu, T. Tang et al., J. Appl. Phys. 111,
053509 (2012).
получается из выражения (C.8), если в него подста-
вить Φk(ζ, t) для слоя толщиной 2τtot, верхний пре-
11.
А. Л. Михайлов, В. А. Огородников, В. С. Сасик
дел интегрирования заменить на 2τtot, а вместо од-
и др., ЖЭТФ 145, 892 (2014).
ного экспоненциального множителя подставить сум-
12.
S. K. Monfared, D. M. Oró, M. Grover et al., J. Appl.
му
Phys. 116, 063504 (2014).
13.
W. T. Buttler, D. M. Oró, R. T. Olson et al., J. Appl.
exp[-ζ-ik0vk(ζ)t]+ exp[-2τtot+ζ-ik0(2vfs-vk(ζ))t].
Phys. 116, 103519 (2014).
14.
В. А. Огородников, А. Л. Михайлов, В. С. Сасик
Такой результат позволяет оценить максимальный
и др., ЖЭТФ 150, 411 (2016).
эффект отражения волн от поверхности на сигнал
15.
В. А. Огородников, А. Л. Михайлов, С. В. Ерунов
обратного рассеяния.
и др., ЖЭТФ 152, 1156 (2017).
Предположение о зеркальном отражении от по-
верхности [28, 63] оправдано по двум причинам. С
16.
M. V. Antipov, V. A. Arinin, A. B. Georgievskaya et
одной стороны, для многократно рассеянного излу-
al., J. Dynam. Behav. Mater. 3, 300 (2017).
чения характер отражения (зеркальный или диф-
17.
O. T. Strand, D. R. Goosman, C. Martinez et al.,
фузный) оказывается несущественным [52, 53]. С
Rev. Sci. Instr. 77, 083108 (2006).
другой стороны, коэффициент отражения излуче-
18.
B. J. Jensen, D. B. Holtkamp, P. A. Rigg et al., J.
ния от металлических образцов близок к единице
Appl. Phys. 101, 013523 (2007).
(для нормального падения r = 0.92 для Pb, r = 0.87
для Sn [55], усредненные по углам коэффициенты
19.
A. R. Valenzuela, G. Rodriguez, S. A. Clarke et al.,
соответственно равны 0.9 и 0.85).
Rev. Sci. Instr. 78, 013101 (2007).
Для реализуемых в экспериментах значений оп-
20.
E. A. Moro, M. E. Briggs, L. M. Hull et al., Appl.
тической толщины пылевых выбросов τtot 1 и
Opt. 52, 08661 (2013).
значений Λ, отвечающих частицам Pb и Sn, влияние
21.
W. J. Warren, E. A. Moro, M. E. Briggs et al., Appl.
отражения волн от свободной поверхности на спектр
Opt. 53, 4661 (2014).
сигнала обратного рассеяния оказывается незначи-
тельным (даже для канала 2 в опыте 1).
22.
М. В. Асташкин, В. К. Баранов, А. Б. Георгиевс-
кая и др., Письма в ЖЭТФ 99, 165 (2014).
23.
G. Prudhomme, P. Mercier, L. Berthe et al., J. Phys.:
ЛИТЕРАТУРА
Conf. Ser. 500, 142022 (2014).
1. В. А. Огородников, А. Г. Иванов, А. Л. Михайлов
24.
G. Prudhomme, P. Mercier, and L. Berthe, J. Phys.:
и др., Физика горения и взрыва 34, 103 (1998).
Conf. Ser. 500, 142027 (2014).
2. W. S. Vogan, W. W. Anderson, M. Grover et al., J.
25.
А. В. Федоров, А. Л. Михайлов, С. А. Финюшин и
Appl. Phys. 98, 113508 (2005).
др., ЖЭТФ 149, 792 (2016).
426
ЖЭТФ, том 157, вып. 3, 2020
Доплеровские гетеродинные измерения. . .
26.
А. В. Федоров, А. Л. Михайлов, Л. К. Антонюк и
45.
C. Roland, T. de Rességuier, A. Sollier et al., J.
др., Физика горения и взрыва 52, 115 (2016).
Dynam. Behav. Mater. 3, 156 (2017).
27.
А. В. Федоров, А. Л. Михайлов, С. А. Финюшин и
46.
D. S. Sorenson, P. Pazuchanics, R. P. Johnson et
др., Физика горения и взрыва 52, 122 (2016).
al., Technical Report LA-UR-14-24722, Los Alamos
National Laboratory (2014).
28.
А. В. Андрияш, М. В. Асташкин, В. К. Баранов и
47.
D. S. Sorenson, G. A. Capelle, M. Grover et al., J.
др., ЖЭТФ 149, 1121 (2016).
Dynam. Behav. Mater. 3, 233 (2017).
29.
J.-E. Franzkowiak, G. Prudhomme, P. Mercier et al.,
48.
S. K. Monfared, W. T. Buttler, D. K. Frayer et al.,
Rev. Sci. Instr. 89, 033901 (2018).
J. Appl. Phys. 117, 223105 (2015).
30.
A. V. Andriyash, M. V. Astashkin, V. K. Baranov et
49.
M. M. Schauer, W. T. Buttler, D. K. Frayer et al., J.
al., J. Appl. Phys. 123, 243102 (2018).
Dynam. Behav. Mater. 3, 217 (2017).
31.
А. В. Федоров, И. С. Гнутов, А. О. Яговкин,
50.
D. S. Sorenson, R. W. Minich, J. L. Romero et al., J.
ЖЭТФ 153, 92 (2018).
Appl. Phys. 92, 5830 (2002).
32.
J. J. Monaghan, J. Comp. Phys. 136, 298 (1997).
51.
А. Исимару, Распространение и рассеяние волн в
случайно-неоднородных средах, т. 1, 2, Мир, Моск-
33.
S. P. Marsh, LANL Shock Hugoniot Data, Univ. of
ва (1981).
California Press, Berkeley (1980).
52.
В. В. Соболев, Рассеяние света в атмосферах пла-
34.
O. Durand and L. Soulard, J. Appl. Phys.
114,
нет, Наука, Москва (1972).
194902 (2013).
53.
H. C. Van De Hulst, Multiple Light Scattering, Acade-
35.
O. Durand and L. Soulard, J. Appl. Phys.
117,
mic Press, New York (1980).
165903 (2015).
54.
P. Laven, Appl. Opt. 42, 436 (2003); http://www.
36.
Л. Д. Ландау, Е. М. Лифшиц, Гидродинамика,
philiplaven.com.
Наука, Москва (2000).
55.
E. D. Palik and G. Ghosh, Handbook of Optical Cons-
37.
L. D. Cloutman, Amer. J. Phys. 56, 643 (1988).
tants of Solids, Academic Press, San Diego (1998).
38.
A. B. Georgievskaya and V. A. Raevsky, AIP Conf.
56.
D. M. Oró, J. E. Hammerberg, W. T. Buttler et al.,
Proc. 1426, 1007 (2012).
AIP Conf. Proc. 1426, 1351 (2012).
39.
A. B. Georgievskaya and V. A. Raevsky, J. Dynam.
57.
L. A. Dombrovsky, Comp. Therm. Sci. 4, 297 (2012).
Behav. Mater. 3, 321 (2017).
40.
G. Dimonte, G. Terrones, F. J. Cherne et al., J. Appl.
58.
Л. И. Седов, Методы подобия и размерности в ме-
Phys. 113, 024905 (2013).
ханике, Наука, Москва (1977).
41.
J.-L. Shao, P. Wang and A.-M. He, J. Appl. Phys.
59.
A. N. Parshikov and S. A. Medin, J. Comp. Phys.
116 (2014).
180, 353 (2002).
42.
S. Dyachkov, A. Parshikov, and V. Zhakhovsky, AIP
60.
J. K. Dukowicz, J. Comp. Phys. 61, 119 (1985).
Conf. Proc. 1793, 100024 (2017).
61.
H. Wendland, Adv. Comp. Math. 4, 389 (1995).
43.
Я. Б. Зельдович, Ю. П. Райзер, Физика ударных
волн и высокотемпературных гидродинамических
62.
M. S. Egorova, S. A. Dyachkov, A. N. Parshikov et
явлений, Физматлит, Москва (2008).
al., Comp. Phys. Comm. 234, 112 (2019).
44.
J. R. Asay and M. Shahinpoor, High-Pressure Shock
63.
J.-E. Franzkowiak, P. Mercier, G. Prudhomme et al.,
Compression of Solids, Springer, New York (1993).
Appl. Opt. 57, 2766 (2018).
427