Функция Грина для случайно-неоднородной среды

Главным образом интерес к вопросу распространения волн в случайно-неоднородных средах (какой является, например, атмосфера) можно объяснить бурным развитием спутниковых технологий. В этом случае становится важной задача расчета характеристик (например, амплитуды) волны, прошедшей через среду, и установления их связей с параметром неоднородности среды. Важную роль здесь играет функция Грина для случайно-неоднородной среды, зная которую можно определить эти характеристики. Рассматривается прохождение света через среду с флуктуирующей диэлектрической проницаемостью.

Волновое уравнение в случайно неоднородной среде. Функция Грина

Рассеяние электромагнитных волн в такой среде определяется системой уравнений Максвелла. Основные отличительные черты рассеяния можно рассматривать для упрощенной модели: скалярного поля $u=u(\mathbf {r} ,\;t)$. Этим скалярным полем заменяются напряженности электрического и магнитного полей; тогда $u=u(\mathbf {r} ,\;t)$ удовлетворяет волновому уравнению:

${\frac {\varepsilon (\mathbf {r} ,\;t)}{c^{2}}}{\frac {\partial ^{2}u(\mathbf {r} ,\;t)}{\partial t^{2}}}-\Delta u(\mathbf {r} ,\;t)=0,$

$\varepsilon (\mathbf {r} ,\;t)=\varepsilon _{0}+\delta \varepsilon (\mathbf {r} ,\;t),$

где $c$ — скорость света в вакууме, $\varepsilon _{0}$ — среднее значение диэлектрической проницаемости, $\delta \varepsilon (\mathbf {r} ,\;t)$ — флуктуации диэлектрической проницаемости. Обратим внимание, что среднее значение диэлектрической проницаемости $\varepsilon _{0}$ предполагается не зависящим от координат и времени, то есть в среднем система однородна и изотропна. Также с хорошей точностью в первом приближении можно считать, что и неусреднённая диэлектрическая проницаемость $\varepsilon (\mathbf {r} ,\;t)=\varepsilon (\mathbf {r} )$ не зависит от времени. Это объясняется тем, что характерные времена, отвечающие за молекулярные процессы в системе, на несколько порядков больше характерных времён электромагнитного поля, поэтому среда как бы «не успевает среагировать» на изменение поля.

Волновое уравнение с такой диэлектрической проницаемостью на самом деле является примером стохастического уравнения, так как в нём присутствует случайный параметр ${\delta \varepsilon (\mathbf {r} )}$. Этот параметр входит в уравнение с помощью умножения, то есть мультипликативно, а не с помощью сложения (аддитивно), как в известном уравнении для броуновского движения.

Описывая рассеяние, интересны характеристики поля $u=u(\mathbf {r} ,\;t)$, усреднённые по флуктуациям диэлектрической проницаемости. Эти характеристики: среднее значение поля $\langle u(\mathbf {r} ,\;t)\rangle$ и интенсивность $I=I(\mathbf {r} ,\;t)$, которую определяет средний квадрат поля (усреднение так же ведётся по флуктуациям диэлектрической проницаемости) $\langle {u(\mathbf {r} ,\;t)}^{2}\rangle$. Статистику флуктуаций считаем заданной, а также учитываем, что усреднённое отклонение от среднего значения диэлектрической проницаемости равно нулю:

$\langle \delta \varepsilon (\mathbf {r} )\rangle =0.$

Начальное однородное волновое уравнение всегда имеет решение в виде $u(\mathbf {r} ,\;t)=0$. Это очевидное тривиальное решение. Легко показать, что при отсутствии флуктуаций ненулевым решением является плоская монохроматическая волна вида:

$u(\mathbf {r} ,\;t)=u_{0}\exp[i\mathbf {k_{0}} \mathbf {r} -iwt].$

Подставим это выражение в волновое уравнение. Получаем:

$-{\frac {w^{2}\varepsilon _{0}}{c^{2}}}u(\mathbf {r} ,\;t)+k_{0}^{2}u(\mathbf {r} ,\;t)=0.$

Отсюда ясно, что предложенное решение будет удовлетворять уравнению, если частота плоской волны $w$ и волновой вектор $k_{0}$ связаны дисперсионным соотношением:

$k_{0}^{2}={\frac {w^{2}\varepsilon _{0}}{c^{2}}}.$

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

Определим функцию Грина $G(\mathbf {r} ,\;t)$. Пусть эта функция является решением начального волнового уравнения, в правую часть которого добавлен расположенный в начале координат монохроматический источник (частота источника $w$). Полагаем, что источник «адиабатически включился в бесконечно далёком прошлом», для этого дополним правую часть множителем $\exp[\alpha t]$, где $\alpha$ — малая положительная величина. В окончательных выражениях будем устремлять её к нулю. Итого функция Грина удовлетворяет уравнению:

${\frac {\varepsilon (\mathbf {r} )}{c^{2}}}{\frac {\partial ^{2}G(\mathbf {r} ,\;t)}{\partial t^{2}}}-\Delta G(\mathbf {r} ,\;t)=e^{-iwt+\alpha t}\delta (\mathbf {r} ).$

Удобно искать решение этого уравнения в виде $G(\mathbf {r} ,\;t)=e^{-iwt+\alpha t}G(\mathbf {r} )$. Подставляя это выражение в уравнение для функции Грина, получаем:

${\frac {\varepsilon (\mathbf {r} )G(\mathbf {r} )}{c^{2}}}{\frac {\partial ^{2}e^{-iwt+\alpha t}}{\partial t^{2}}}-e^{-iwt+\alpha t}\Delta G(\mathbf {r} )=e^{-iwt+\alpha t}\delta (\mathbf {r} ).$

От двойного дифференцирования экспоненты по времени появится множитель $-(w+i\alpha )^{2}$, тогда получаем уравнение для функции $G(\mathbf {r} )$:

${\frac {-\varepsilon (\mathbf {r} )(w+i\alpha )^{2}}{c^{2}}}G(\mathbf {r} )-\Delta G(\mathbf {r} )=\delta (\mathbf {r} ).$

Нужно решить это уравнение для некоторой диэлектрической проницаемости $\varepsilon (\mathbf {r} ),$ а затем усреднить это решение по всевозможным отклонениям $\delta \varepsilon (\mathbf {r} )$. Но оказывается, что нет возможности получить решение этого уравнения при произвольной диэлектрической проницаемости, поэтому решение ищется с использованием теории возмущений, полагая отклонение $\delta \varepsilon (\mathbf {r} )$ малой величиной.

Функция Грина для среды без флуктуаций диэлектрической проницаемости

Для начала необходимо найти функцию Грина $G(\mathbf {r} ,\;t)$, отвечающую волновому уравнению без отклонений диэлектрической проницаемости, то есть $\varepsilon (\mathbf {r} )=\varepsilon _{0}$:

${\frac {\varepsilon _{0}}{c^{2}}}{\frac {\partial ^{2}G_{0}(\mathbf {r} ,\;t)}{\partial t^{2}}}-\Delta G_{0}(\mathbf {r} ,\;t)=e^{-iwt+\alpha t}\delta (\mathbf {r} ).$ (1)

Снова ищем решение в виде $G_{0}(\mathbf {r} ,\;t)=e^{-iwt+\alpha t}G_{0}(\mathbf {r} )$. Тогда $G_{0}(\mathbf {r} )$ удовлетворяет уравнению:

$-k_{0}^{2}G_{0}(\mathbf {r} )-\Delta G_{0}(\mathbf {r} )=\delta (\mathbf {r} ),$ (2)

где величиной $k_{0}={\frac {{\sqrt {\varepsilon _{0}}}(w+i\alpha )}{c}}$. Видно, что у $k_{0}$ присутствует мнимая положительная часть, далее нам это понадобится. Уравнение $(2)$ удобно решать с помощью преобразования Фурье вида:

$F(\mathbf {k} )=\int \limits _{-\infty }^{\infty }F(\mathbf {r} )e^{-i\mathbf {k} \mathbf {r} }\,d\mathbf {r} ,$ (3)

$F(\mathbf {r} )={\frac {1}{(2\pi )^{3}}}\int \limits _{-\infty }^{\infty }F(\mathbf {k} )e^{i\mathbf {k} \mathbf {r} }\,d\mathbf {k} ,$ (4)

Выражение $(3)$ — прямое Фурье-преобразование, $F(\mathbf {k} )$ — Фурье-образ функции $F(\mathbf {r} )$, выражение $(4)$ — обратное Фурье-преобразование. Образ функции Грина $G_{0}(\mathbf {r} )$ будем обозначать через $G_{0}(\mathbf {k} )$. Применяя Фурье-преобразования к уравнению $(2)$ и учитывая, что $\delta$-функция является Фурье-образом единицы, получаем:

$(k^{2}-k_{0}^{2})G_{0}(\mathbf {k} )=1,$ (5)

$G_{0}(\mathbf {k} )={\frac {1}{k^{2}-k_{0}^{2}}}.$ (6)

Чтобы получить функцию $G_{0}(\mathbf {r} )$, делаем обратное Фурье-преобразование $G_{0}(\mathbf {k} )$:

$G_{0}(\mathbf {r} )={\frac {1}{(2\pi )^{3}}}\int \limits _{-\infty }^{\infty }{\frac {e^{i\mathbf {k} \mathbf {r} }}{k^{2}-k_{0}^{2}}}\,d\mathbf {k} .$ (7)

Будем вычислять этот интеграл в сферической системе координат, выбрав полярную ось вдоль вектора $\mathbf {r}$ (под полярной осью понимается ось, от которой отсчитывается угол $\theta$):

$G_{0}(\mathbf {r} )={\frac {1}{(2\pi )^{3}}}\int \limits _{-\infty }^{\infty }{\frac {e^{i\mathbf {k} \mathbf {r} }}{k^{2}-k_{0}^{2}}}\,d\mathbf {k} ={\frac {1}{(2\pi )^{3}}}\int \limits _{0}^{2\pi }\,d\varphi \int \limits _{0}^{\pi }\sin \theta \,d\theta \int \limits _{0}^{\infty }k^{2}{\frac {e^{ikr\cos \theta }}{k^{2}-k_{0}^{2}}}\,dk=$

$={\frac {2\pi }{8\pi ^{3}}}\int \limits _{\pi }^{0}\,d\cos \theta \int \limits _{0}^{\infty }k^{2}{\frac {e^{ikr\cos \theta }}{k^{2}-k_{0}^{2}}}\,dk={\frac {1}{4\pi ^{2}}}\int \limits _{0}^{\infty }{\frac {k^{2}}{k^{2}-k_{0}^{2}}}{\Bigl [}{\frac {e^{ikr\cos \theta }}{ikr\cos \theta }}{\Bigr ]}_{\pi }^{0}\,dk=$

$={\frac {1}{4\pi ^{2}}}\int \limits _{0}^{\infty }{\frac {k^{2}}{k^{2}-k_{0}^{2}}}{\Bigl (}{\frac {e^{ikr}}{ikr}}-{\frac {e^{-ikr}}{ikr}}{\Bigr )}\,dk={\frac {1}{4\pi ^{2}ir}}{\frac {1}{2}}\int \limits _{-\infty }^{\infty }{\frac {k(e^{ikr}-e^{-ikr})}{k^{2}-k_{0}^{2}}}\,dk=$

$={\frac {1}{8\pi ^{2}ir}}\int \limits _{-\infty }^{\infty }{\frac {ke^{ikr}}{k^{2}-k_{0}^{2}}}\,dk-{\frac {1}{8\pi ^{2}ir}}\int \limits _{-\infty }^{\infty }{\frac {ke^{-ikr}}{k^{2}-k_{0}^{2}}}\,dk={\frac {2\pi ie^{ik_{0}r}}{8\pi ^{2}ir2}}+{\frac {2\pi ie^{ik_{0}r}}{8{\pi }^{2}ir2}}={\frac {e^{i{k_{0}}r}}{4\pi r}}.$

Для вычисления интеграла по сферическим координатам мы воспользовались четностью функции ${\frac {k(e^{ikr}-e^{-ikr})}{k^{2}-k_{0}^{2}}}$, а также последние интегралы брались по вычетам. Для первого слагаемого контур интегрирования замыкался сверху, в этой полуплоскости затухает $e^{ikr}$, тогда вычет берётся в $k=k_{0}={\frac {{\sqrt {\varepsilon _{0}}}(w+i\alpha )}{c}}$. Для второго слагаемого замыкали контур в нижней полуплоскости, и тогда срабатывает вычет в точке $k=-k_{0}$, при этом необходимо не забыть, что контур обходится по часовой стрелке, тогда как в теореме по вычетам используется обход против часовой стрелки. Направление обхода можно легко изменить, добавив во втором слагаемом множитель $(-1)$.

Итоговое выражение для функции Грина будет:

$G_{0}(\mathbf {r} ,\;t)={\frac {e^{ik_{0}r-iwt+\alpha t}}{4\pi r}}.$

Это расходящаяся сферическая волна. Амплитуда этой волны убывает как ${\frac {1}{r}}$ по мере удаления от источника.

Функция Грина с учётом флуктуаций диэлектрической проницаемости

Перепишем уравнение

$-{\frac {\varepsilon (\mathbf {r} )(w+i\alpha )^{2}}{c^{2}}}G(\mathbf {r} )-\Delta G(\mathbf {r} )=\delta (\mathbf {r} )$

в виде

$(-k_{0}^{2}-\Delta )G(\mathbf {r} )={\frac {\delta \varepsilon (\mathbf {r} )}{\varepsilon _{0}}}k_{0}^{2}G(\mathbf {r} )+\delta (\mathbf {r} ).$

Для использования теории возмущения, в которой мы будем считать $\delta \varepsilon (\mathbf {r} )$ малой величиной, удобнее перейти к интегральному аналогу предыдущего уравнения:

$G(\mathbf {r} )=G_{0}(\mathbf {r} )+{\frac {{k_{0}}^{2}}{\varepsilon _{0}}}\int G_{0}(\mathbf {r} -\mathbf {r_{1}} )\delta \varepsilon (\mathbf {r_{1}} )G(\mathbf {r_{1}} )\,dr_{1}.$

Тогда можно легко написать итерационное решение в виде ряда:

$G(\mathbf {r} )=G_{0}(\mathbf {r} )+{\frac {k_{0}^{2}}{\varepsilon _{0}}}\int G_{0}(\mathbf {r} -\mathbf {r_{1}} )\delta \varepsilon (\mathbf {r_{1}} )G_{0}(\mathbf {r_{1}} )\,d\mathbf {r_{1}} +$

$+{\frac {k_{0}^{4}}{\varepsilon _{0}^{2}}}\int G_{0}(\mathbf {r} -\mathbf {r_{1}} )\delta \varepsilon (\mathbf {r_{1}} )\int G_{0}(\mathbf {r_{1}} -\mathbf {r_{2}} )\delta \varepsilon (\mathbf {r_{2}} )G_{0}(\mathbf {r_{2}} )\,d\mathbf {r_{1}} \,d\mathbf {r_{2}} +\ldots$

Величина $G(\mathbf {r} )$ — случайная величина. В дальнейшем её необходимо усреднить по всевозможным флуктуациям диэлектрической проницаемости. Это представляет собой следующий трудоёмкий шаг.

Литература

• Рытов С. М., Кравцов Ю. А., Татарский В. И. Введение в статистическую радиофизику. — Ч. 2. Случайные поля. — М.: Наука, 1978.
• Исимару А. Распространение и рассеяние волн в случайно-неоднородных средах. — Т. 1, 2. — М.: Мир, 1981.