УДК 622.179
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ И РАСЧЕТ РАСПРОСТРАНЕНИЯ НАПРАВЛЕННЫХ ИМПУЛЬСОВ В ОДНОРОДНЫХ ПОГЛОЩАЮЩИХ СРЕДАХ
© 2011 Б.В. Скворцов, И.А. Лезин, А.В. Солнцева
Самарский государственный аэрокосмический университет имени академика С.П.Королева (национальный исследовательский университет)
Поступила в редакцию 26.04.2011
Импульсное зондирование применяется для исследования физико-химических свойств сред и объектов. Импульсные сигналы, взаимодействующие с контролируемой средой или объектом неисчерпаемы по своим информационным и технологическим возможностям. Любая среда или объект, в которых распространяется и с которыми взаимодействует импульс, характеризуется комплексным волновым вектором, компоненты которого определяются их свойствами. В работе с общих позиций, математически описывается энергетическое поле, созданное излучающим импульсом при распространении его в однородной поглощающей среде. Дается трехмерное математическое описание процесса распространения направленных импульсных сигналов, распространяющихся в поглощающих средах. Изменение амплитуды и спектрального состава зондирующего импульса, прошедшего через исследуемую среду дает диагностическую или прогнозирующую информацию о среде. Разработан алгоритм и приведены примеры расчета распространения электромагнитного импульса в пассивной поглощающей среде, иллюстрирующие искажение формы импульса, пришедшего в заданную точку пространства.
Ключевые слова: направленный импульс, волновой вектор, поглощающая среда, распространение, математическое описание, алгоритм расчета.
Импульсное зондирование известно и применятся для исследования физико-химических свойств сред и объектов [1]. Импульсные сигналы, взаимодействующие с контролируемой средой или объектом неисчерпаемы по своим информационным и технологическим возможностям.
Направленный зондирующий импульс любой физической природы и произвольной формы, действующий в какой- либо точке пространства в общем случае является вектором (рис. 1).
Импульс может быть произвольной формы и описываться описывается выражением:
р(О,t) = ipx(0,t) + jpy(0,t) + Kpz(0,t), (1)
где под точкой "О" будем понимать координаты точки излучения (х0, y, z). Отметим, что мы намеренно не совмещаем точку излучения с началом координат, так как среда может облучаться синхронно из нескольких точек, что может привести в дальнейшем к неточности описания. Импульс исчерпывающим образом описывается комплексной спектральной амплитудой S(jrn), которая определяется как преобразование Фу-
Скворцов Борис Владимирович, доктор технических наук, профессор, научный руководитель НИЛ «Аналитические приборы и системы». E-mail: aps@ssau.ru. Лезин Илья Александрович, кандидат технических наук, ассистент кафедры информационных систем и технологий. E-mail: aps@ssau.ru.
Солнцева Александра Валерьевна, инженер НИЛ «Аналитические приборы и системы». E-mail: aps@ssau.ru.
рье от исходного сигнала и представляет собой набор синусоидальных сигналов разной амплитуды и фазы.
S( О, ja ) = F[p(O, t)] = J p( О, t )e - Jatdt =
—TO
= iSx(O,ja) + ~jSv(O,ja) + KS JO, ja)
(2)
Изменение амплитуды и спектрального состава зондирующего импульса, прошедшего через исследуемую среду, дает информацию о среде или объекте. Распространяясь в среде, импульс можем в общем случае изменять форму и направление, что определяется свойствами среды. Рассмотрим распространение импульса в однородной среде, в которой направление импульса не изменяется.
В работе [5], на основе анализа работ [2,3], получено аналитическое выражение, моделирующее распространение импульса в однородной поглощающей среде, которое имеет вид:
p(r,t) = 1 J S (O, j®y[t—k(ja>r ]a
2n
— ] Jp(O,j®y[)a—k(ja)r]Tda 2 П
(3)
где г = i (х - x0)+ }(у - у0)+ к(2 - z0 ) - радиус-вектор заданной точки пространства относительно точки излучения,
—то
а) б)
Рис. 1. Иллюстрация к постановке задачи распространения импульса в неподвижной среде
к 0'®) = 1кх 0'® )+ 'ку ('®) + К кг 0® ) = = к'(®)-'к''(®)
- волновой вектор, (4)
кх ') = к'х (®) + ]к"х{с®),
ку с/®)=ку (®)+к(®), _ кг(/®)=к;(®)+ 'к;(®) к'(®) = 7к'(®) + 'к'(®) + Кк'(®) , к"(® ) = Тк"(® ) + 'к"(® ) + Кк"(® ) .
При этом к' (®) - определяется фазовой скоростью распространения энергетической волны в среде [2]:
= ®Уф (®)
(5)
к' (®,Я ) = ■
®
Уф(®У \ф(®)2 ;
к" (®) = а(®) - коэффициент поглощения среды.
Здесь и далее следует различать ' и к - как орты декартовой системы координат, к - как волновой вектор, ' = ^/-Г ;
Г = \г\ = ^1 (X - х0)2 + (у - у0)2 + (г - ;0)2 , (6)
к (]®)-Г =(х - х0 )кх 0®) +
+ (У - Уо)ку 0®)+(; - ;о)К 0®)
(7)
ной, проявляющейся, например, в виде нагревания среды, и реактивной, определяемой электромагнитными или квантово-механическими процессами в зоне следования импульса.
Учитывая (4) можно записать:
®УеЛ* -т)®-к' '
Чтс1а
(8)
1 да да _
РМ = Iр(о,тУк"с
2П-да-да
Отметим, что в общем случае для неоднородных сред волновой вектор является функцией координат, что приводит к тому, что импульс может двигаться по криволинейной траектории и менять свое направление. В однородной среде импульс не может изменить своего направления, поэтому вектора рСо, г), к (у®), к' (®), к" (®), г, р(г, г) колле-ниарны, то есть совпадают до направлению. Из коллениарности следует, что к'(® ) ■ г = к'(®) г , к"(®) ■ г = \к "(®)| г
Тогда для однородных сред выражение (8) можно переписать в виде:
р(г,г) = — да (9)
- скалярное произведение волнового вектора на радиус вектор заданной точки пространства, т- формальный параметр интегрирования.
Отметим, что в классической теории волновых процессов в не поглощающих однородных средах волновой вектор определяет число длин волн, умещающихся на единице длины [1/м]. В поглощающих средах волновой вектор также имеет размерность [1/м], однако здесь его физическая трактовка больше согласуется с комплексным коэффициентом поглощения пакета энергии распространяющегося в среде, которая по аналогии с электричеством может быть актив-
1 2п
Выражение (9) связывает параметры импульса, появившегося в точке "г"пространства с параметрами импульса, запущенного в точке "О". Оно определяет каждую координатную составляющую импульса:
1
Рх,(у),(; )М = — I I Рх,(у),(; ){О,Т)Х
2-да -да
-¡к'Щге'[-г®-\Р(® Х^®
(10)
х е 1 1 е
Сигнал, вычисленный по формулам (3), (9), (10) в общем случае является комплексным:
р(г,г) = Яе(г,г) + '1т(г,г). (11)
Как указано в [2, 3], физический смысл следует придавать действительной части выражения (13). Однако при анализе энергетических характеристик импульса не следует отказываться от мнимой составляющей выражения (11).
Прошедший через однородную среду импульс, в общем случае изменяет свою форму, скорость и в точке приема несет в себе информацию о свойствах среды или об объекте. При этом важно знать амплитудный спектр принятого импульса, который определяется как прямое преобразование Фурье:
S(r, j—) = F[p(r, t)]. (12)
Следует отметить, что, как правило, все приемники импульсных сигналов, являются коорди-натно-ориентированными, поэтому важно вычислять координатные составляющие прошедшего через среду импульса и соответствующие им амплитудные спектры Sx (r,jm), Sy (r,jm), Sz (r,jm).
Приведенные выражения, определяющие в общем виде форму и местоположение зондирующего импульса, прошедшего через среду, являются математической основой для исследования сред и объектов методом импульсного зондирования. Задавая параметры среды, в виде волнового вектора, определяющего коэффициент поглощения и фазовую скорость распространения, с учетом направления излучающего импульса можно численно моделировать прохождение импульса в пространстве или объекте, наблюдать искажение его формы, определять место и время прибытия в какую-либо точку пространства.
Выражение (10) является скалярным и определяет каждую координатную проекцию импульса, пришедшего в заданную точку среды как функцию времени и координат. Формулы (3), (9) и (10) предполагают, что импульс может появиться только в той точке пространства, которая лежит на направлении распространения излученного импульса.
Для простых однородных сред, без частот-нойзависимости_(диссипации) волнового вектора k(m) = k' - jk" выражение (9) примет вид:
-\к"\г
где ¡и0 =1,256637-10-6 Гн/м, е0 =8,85416 10-12 Ф/ м, а0 , ц0, е0 удельная электрическая проводимость, относительные магнитная и диэлектрическая проницаемости среды соответственно,
к ' =
а2 + s2s20m2 + ss0— )
2
коэффициент распространения;
к' =
jj—
(-Jа2 + s2s02)m2 -ss0—)
2
(15)
(16)
коэффициент затухания.
Фазовая скорость определяется по формуле:
уФ =—= Ф г
2m
¡2 2 2 2 (1/) jj0(s0srn+^a +ss0m )'у '
Определим вид функции, которая описывала бы форму сигнала в некоторой точке {х, у, . Для простоты сделаем допущение, что зондирующий сигнал имеет только одну координатную составляющую р(0, ^) = грх (О, ^), направленную только вдоль оси х. При расчете воспользуемся формулой (10).
Для проверки полученных выражений рассмотрим простейший случай, когда проводимость среды а = 0 и волновой вектор содержит только действительную часть:
k(j—) = —^jjo) ,
(18)
с - скорость света в вакууме.
Подставив данное выражение в формулу (10), мы приходим к следующей форме записи:
1 ш ш j\ (t-T—-x—JJs I
Px (r,t) = — JJ Px (O,T) ^ c } dTdm =
^ —rf\—rf\
1 ш ш / Ч -j—\ T-t+—-[jE
- k r ш ш r, . ,-,, j 1 I ,
p(r,t)= e- J J p(Ot) [-T)m-\kj)r 'dTdm .(13) = T J J Px OT
/77" 2 ТС —m —m
dTdm =
2п
о -|Р|-Г
Здесь сомножитель е 1 1 определяет изменение амплитуды сигнала по координатам, все остальное описывает форму импульса во времени в заданной точке г.
Рассмотрим пример расчета распространения электромагнитного импульса в однородной среде, выпущенного из начала координат 0(0,0,0). Для электромагнитных сигналов в неограниченной однородной среде волновой вектор является скаляром, зависит от частоты и определяется выражениями [1, 4]:
= J Px (0,t)S\t-\t - ^^VJSj |dT =
(19)
= Px\ 0,t-xVJS| = Px(0,t--x-),
Уф
k(—) = y]mj—)i0 [s—S— - ja——] = = k'—)- jk'—) '(14)
где УФ = С- фазовая скорость распространения сигнала.
То есть в данном случае сигнал приходит в искомую точку с задержкой, определяемой фазовой скоростью распространения электромагнитных импульсов в рассматриваемой среде, но без затухания и искажения формы, что хорошо известно и вполне очевидно.
-ш
Рассмотрим общий случай, когда волновой Учитывая, что к'х (0)= 0 и к"х (0)= 0, выра-
вектор имеет действительную и мнимую состав- жение (25) приводится к виду: ляющие и описывается выражениями (14)-(16). Функция (10) примет следующий вид:
1 да да
Px М = —¡ í Px {0,r)>
— rr-, — rv-,
х е'((-т®кх (® к ®)х )с1ж1®
Зададим зондирующий сигнал в виде прямо угольного импульса:
Px (O,t) =
0,t < 0, A,0 < t < At, 0,t > At.
(21)
Подставим это выражение в (10) и после преобразований получим:
1 да At p-( )=^ííAe
¿((t-T)c^-(k'x (w)-jk- (w))x
drdw =
л да -k- (w).
Are xy ' п - w
-да
-k"x (w)
a r e ~ . i mAt j ( к / \ wAtj , — I-sinI-IcosI tw-kx(w)x--Idw +
п w 12 I l 2 I
(22)
wAt
jA да e 1 w>" . (wAt \.( wAt j,
— I-sinI-IsinI tw-kx(w)x--Idw.
п w 1 2 I l 2 I
, () A да e"k'(w)x (wAt
fx k )=-í -Щ —
П -да w l 2 _
x cos|(tw - k'x (w)x - wAjdw
%
{(,w) =
:>-k''(w)x
(wAt
-sm\ -
w
x cosí tw -k'x(w)x-
l 2 wAt I
~2~ J
(24)
%
(t,0 )= lim
e
- k''(w).
lx . (wAtj Sin\ - I x
w^Q w l 2
wAt j
x cos\ tw - k'x(w)x -
<px (t,0 )= lim
w^Q
(e0 wAt
— sin-cos0
w2
At 2.(26)
Таким образом, приближенно интеграл (23) (20) вычисляется по формуле:
fx (((t,0)+l.(%x (t,iAw)+
П i=i
+ %(t,-iAw))-%x((,nAw)+%x((,-nAwm\Áw (27)
2
Из условий допустимой погрешности интегрирования найдем значения n и Aw .
Величина шага интегрирования Aw определяет погрешность численного вычисления интегралов периодических тригонометрических функций. Например, для обеспечения погрешности вычисления A < 1% тригонометрической функции с периодом Г достаточно, чтобы шаг интегрирования был равен Aw = Т/20.
Период функции sin (wAt¡2):
т = 4п
sin = At '
Для определения искомой формы результирующего сигнала используется действительная часть выражения (22):
Период функции cos\ tw - k'x(w)x -
(28) wAt I
T J
определяется через приближенное значение параметра k x (w):
(23) k'x(w) =
LLow
2 , 2 2 2, o +s s0w +ssqw
Решение задачи возможно только методами численного интегрирования. Приближенное вычисление значения интеграла в выражении (23) осуществляется по методу трапеций. Обозначим подынтегральную функцию как:
Liliqwo O
l1-^-0-,w << —
2 ss0
(29)
w
л/^Д
0ss0 ,w >>
o
ss0
Таким образом, оценить период можно для обоих случаев:
Оценка значения функции (24) в точке ®= 0 определяется через предел:
Tcos (w) ¡
2п
t-AL - x \LLQ0
2 V 2w 2п
,w<<
o
ssr
At
t - — - xy LL
0ss0
,w >>
O.(30)
ss.
Из формулы (29) видно, что шаг интегрирова-(2 5) ния в общем случае является функцией от частоты:
х
2
х
Лс(с) =
тгп
п[т,Тсо* (с)]
20
(31)
Таким образом, выражение (27) преобразуется к новому виду, что означает переход к интегрированию с переменным шагом:
(
п
/х (,п)~— РРх (г,0)+1(Рх ((-1 + Лс(с,-1 ))-
г=1
+ Рх(, с— _ Лсо(со1_1)))О)
(32)
Следующим действием является определение п или границ интегрирования. Для этого необходимо вычислить значение п, после которого вычислением оставшейся части интеграла (23) можно пренебречь. При вычислении значения функции (32) значение п можно не задавать заранее, а определить его по ходу выполнения вычислений.
При достаточно больших значениях с составляющая волнового вектора к"х(сС практически не меняется:
к хСС^
ИИ0а\
2 . 2 2 2 а + е гпю -ее с
а^ц/
Л /ее,
.(33)
Обозначим Лст приращение частоты, соответствующее периоду Тгп -Тс03 (с) и, подставив их, рассмотрим следующий интеграл:
1х() = - I П с
А с+лст е Н ее0
( сЛг
-£/п| - |Х
с
V 2
(
Х С08
V
I- ЛгЛ
с| г - цц0ее0 - —
(34)
dс
При бесконечном увеличении значения с интеграл (34) стремится к 0, поэтому ряд (32) сходится при п ^ да и может быть вычислен с заданным требованием точности. Если значения функции (32) практически не отличаются при п и п + Лст /Лс , то интеграл (23) можно считать вычисленным.
Общий алгоритм вычисления формы результирующего импульсного сигнала, прошедшего через среду при зондировании её сигналом с заданными характеристиками имеет следующий вид.
Шаг 1. Задание границ отображения графика г1 и гг, определение количества точек т на графике, вычисление шага приращения
Лг = (гг _ г1 )/(т _ 1).
Шаг 2. Начальная точка графика г = г1.
Шаг 3. По формуле (32) вычисляется приближенное значение результирующего сигнала в точке г .
Шаг3.1. Задается начальное значение п = 0.
Шаг 3.2. По формулам (28) и (30) определяются периоды, вычисляется приращение Лс по формуле (31).
Шаг 33. Рассчитывается значение функции /х (г,п) по формуле (32).
Шаг 3.4. Увеличивается приближение п = п +1, соответственно с = с + Лс .
Шаг 3.5. Если
\/х (г,п/2)_ /х (,п)
меньше
/х Ы
некоего наперед заданного е, то приближенное вычисление точки считается завершенным. В противном случае возвращаемся к шагу 3.2.
Шаг 4. Аналогично шагу 3 вычисляется значение мнимой части в точке г.
Шаг 5. На графиках отмечаются вещественная и мнимая точки, а также величина модуля сигнала в точке г .
Шаг 6. г = г + Лг. Если г > гг, то алгоритм завершается, иначе переходим к шагу 3.
Аналогично можно вычислить мнимую часть выражения (22). Спектр исходного сигнала рассчитывается по формуле прямого преобразования Фурье:
тш
Ксх (с)= I Рх (0,г) - =
да
Л _ ш , 2 А . сЛг = А|е }с^г =-¡¡т-
(35)
сЛг
2
Учитывая, что результирующий сигнал не имеет аналитического выражения, его спектр вычисляется по табличным значениям:
т г1+] / \
Sx(г,}с)=Ъ I ((г,г)+}рхтМ>) =
г=0 г,
гп "1+1 гп 1+1
= 1 р& (г,г)| е-}с}^г + }Iрхт (г,г)| е_}сгЛ
(36)
г=0
г=0
где р^е(г,г) и рхтг (г,г) - табличные значения действительной и мнимой части результирующего сигнала, соответственно.
На рис. 2 приведены результаты расчета распространения прямоугольного электромагнитного импульса длительностью Лг = 10_7 сек, амплитуды А = 1, через среду с проводимостью а = 10 ~17 Сим, относительной диэлектрической проницаемостью е = 2, магнитной проницаемостью ц = 1 для точек контроля, находящихся на расстоянии 1 м, 50 м и 100 м от точки излучения. Графики иллюстрируют время прихода импульса в заданную точку, затухание и искажение его формы в зависимости от расстояния. При этом действительная часть импульса соответствует
2
его физическому содержанию, которое наблюдается на экране осциллографа. Мнимая часть позволяет с высокой точностью фиксировать временные параметры импульса (момент прихода и длительность), так как на фронтах импульса его мнимая составляющая имеет острый, явно выраженный экстремум. Причем минимум на переднем фронте, максимум - на заднем. Это дает возможность идентифицировать "тело" импульса при периодических импульсных облучениях. Модуль позволяет с высокой точностью определить энергетические характеристики, вычислить полные, активные и реактивные потери, обусловленные движением импульса в среде.
На рис. 3 показаны результаты расчета прохождения электромагнитным импульсом расстояния 100 м при различных значениях диэлектри-
ческой проницаемости среды ( е = 1 ,е = 2 ,е = 4 ) и постоянных проводимости (а = 10~17 Сим) и магнитной проницаемости / = 1.
Анализ графиков показывает, что диэлектрическая проницаемость среды значительно влияет на скорость движения импульса. С ростом диэлектрической проницаемости время прибытия импульса в заданную точку увеличивается, затухание возрастает, что полностью соответствует аналитическим выражениям (16), (17) и подтверждает правильность разработанной методики расчетов.
На рис. 4 показаны результаты расчета прохождения электромагнитным импульсом расстояния 100 м при различных значениях проводимости среды а = 10~10 , а = 10~5, а = 10~4 Сим и постоянных диэлектрической е = 2 и магнитной /= 1 проницаемостей. Анализ графиков
Рис. 2. Составляющие импульса при разных расстояниях до точки контроля
Рис. 3. Составляющие импульса на расстоянии 100 м при изменения е
Рис. 4. Составляющие импульса на расстоянии 100 м при изменения а
показывает, что проводимость среды сильно влияет на затухание сигнала, искажает его форму, но при выбранных значениях практически не повлияло на время прилета импульса, что вполне объяснимо, так как групповая скорость мало зависит от проводимости.
СПИСОК ЛИТЕРАТУРЫ
1. Исследование объектов с помощью пикосекундных импульсов / Г.В.Глебович, А.В.Андриянов, Ю.В.Вве-
денский. М.: Радио и связь, 1984. 256 с.
2. Гинзбург В.Л. Распространение электромагнитных волн в плазме. М.: Наука, 1967. 684 с.
3. Вайнштейн Л.А. Распространение импульсов // Успехи физических наук. 1976. Т.118. Вып. 2. С.339-369.
4. Никольский В.В. Электродинамика и распространение радиоволн. М.: Наука, 1978. 444 с.
5. Скворцов Б.В. Математическое моделирование распространения направленных импульсных сигналов в поглощающих средах // Авиакосмическое приборостроение. 2010. №12. С.28-32.
MATHEMATICAL MODELING AND CALCULATION OF THE EXTENDED DIRECTION PULSE IN HOMOGENEOUS ABSORBING ENVIRONMENTS
© 2011 B.V. Skvortsov, I.A. Lezin, A.V. Solntseva
Samara State Aerospace University
Impulsive probing used to study physical and chemical properties of environments and objects. Pulsed signals, interacting with the controlled environment or subject, endless to information and technological capabilities. Any environment or object, in which extends and with which pulse interacted, is characterized by a complex wave vector, whose components are defined by their properties. Working with common positions mathematically described the energy field, created by the emitted pulse as it propagates in the homogeneous absorbing environments. Given a three-dimensional mathematical description of the distribution process directed pulse signals, propagating in absorbing environments. Change of amplitude and spectral composition of the probe pulse, transmitted through the environment, gives a diagnostic or predictive information about the environment. Developed an algorithm and examples of calculating the propagation of an electromagnetic pulse in passive absorbing environment, illustrating the distortion of the pulse, arriving at a given point in the area. Key words: directed pulse, wave vector, absorbing environment, dissemination, mathematical description, algorithm for calculating.
Boris Skvortsov, Doctor of Technics, Professor, supervisor at the Research Laboratory «Analytical Devices and Systems». E-mail: aps@ssau.ru.
Ilya Lezin, Candidate of Technics, Assistant Lecturer at the Information Systems and Technology Department. E-mail: aps@ssau.ru.
Alexandra Solntseva, Engineer at the Research Laboratory «Analytical Devices and Systems». E-mail: aps@ssau.ru.