Например, Бобцов

Установление закономерностей и моделирование диффузионного режима хаотической генерации в сильно рассеивающих средах

УДК 535.374
УСТАНОВЛЕНИЕ ЗАКОНОМЕРНОСТЕЙ И МОДЕЛИРОВАНИЕ ДИФФУЗИОННОГО РЕЖИМА ХАОТИЧЕСКОЙ ГЕНЕРАЦИИ В  СИЛЬНО РАССЕИВАЮЩИХ СРЕДАХ

© 2012 г.

В. П. Ящук, канд. физ.-мат. наук; М. В. Журавский; О. А. Пригодюк Киевский национальный университет имени Тараса Шевченко, г. Киев, Украина Е-mail: yavasil@ukr.net

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

Ключевые слова: хаотическая генерация, диффузионный режим, моделирование, м­ етод Монте-Карло, сильно рассеивающая активная среда.

Коды OCIS: 290.4210, 140.3460.

Поступила в редакцию 03.05.2012.

Введение
Многократное упругое рассеяние способно существенно увеличивать поглощение и вынужденное излучение в среде, заполняющей пространство между рассеивающими частицами, благодаря увеличению пути, пройденному излучением в этой среде. Последнее обстоятельство приводит к так называемой хаотической генерации1 (ХГ) – возникновению излучения с  узким спектром и некоторыми другими характеристиками, свойственными лазерному излучению [1, 2]. При увеличении накачки спектр излучения сужается до минимального значения его ширины, которое в дальнейшем практически не изменяется, а интенсивность излучения в максимуме спектра проявляет характерную для лазера S-образную зависимость. Но, в отличие от обычного лазера, переход от режима спонтанного излучения к  генерации происходит не скачкообразно, а в конечном
1 Англ. random lasing.

диапазоне изменения интенсивности накачки (усиления). Модам хаотического лазера (ХЛ) соответствуют не резонансы резонатора, а различные траектории фотонов в среде, вследствие чего излучение ХЛ является ненаправленным, а спектр – сплошным, который не содержит узких линий, характерных для мод резонатора.
Возникающие экспериментальные закономерности обусловлены интегральным проявлением тесно переплетающихся эффектов вынужденного и спонтанного излучения, поглощения, упругого многократного рассеяния, а также просветления среды и отражения на границах образца. Поэтому выяснение особенностей формирования излучения ХГ является достаточно сложной задачей, которая требует комплексного подхода, включающего как экспериментальные исследования, так и моделирование. Существующие работы по моделированию ХГ направлены преимущественно на выяснение общих черт формирования излучения ХГ или носят иллюстративный характер [3–7]. В  настоящей работе делается попытка экспе-

30 “Оптический журнал”, 79, 9, 2012

риментально и с помощью моделирования изучить влияние параметров активной рассеивающей среды на характеристики излучения хаотической генерации, а также динамику развития ХГ. Условия, при которых наблюдается ХГ, соответствуют диффузионному режиму распространения света2 в многократно рассеивающей активной среде, при котором фазовые соотношения между рассеянными волнами можно не учитывать. В этом случае наиболее удобным для моделирования является метод Монте-Карло (ММК) [4, 7–12].

1. Экспериментальные закономерности хаотической генерации в суспензиях диэлектрических частиц

На рис.  1,  2 приведены эксперименталь-

ные характеристики излучения ХГ в много-

кратно рассеивающей среде, представляющей

собой высококонцентрированную суспензию

диэлектрических частиц в поливинилацетат-

ной (ПВА) матрице с растворенным в ней ро-

дамином 6Ж. Оптическая накачка образцов

осуществлялась второй гармоникой ИАГ-Nd3+-

лазера (l = 532 нм, tимп = 15 нс). Данные рис. 1 соответствуют суспензии высокопреломляю-

щих частиц (боразон (BN) – показатель преломления n  =  2,1, средний диаметр d–  =  2,5  мкм,

концентрация c = 2,5×109 см–3) в концентрированном растворе красителя (Ск = 3 мМ/л). Данные рис.  2 соответствуют суспензии низкопре-

ломляющих d–  =  1  мкм,

частиц (аэросил c  =  3×1011  см–3

(SiO2)  – n  =  1,46, (1), 1,7×1011 (1′),

5×1010 (1′′)) в менее концентрированном раство-

ре красителя (Ск  =  1  мМ/л). Длина свободного пробега фотона, рассчитанная как lc  =  (sc)–1 (s  – площадь поперечного сечения частицы),

составляет 77  мкм в суспензии частиц боразо-

на (СБ) и 4, 7,5 и 24  мкм в суспензиях аэро-

сила (СА) соответственно приведенным выше

концентрациям, что при использованной тол-

щине образцов h  =  1,5  мм обеспечивает много-

кратное рассеяние. Таким образом, кратность

рассеяния, пропорциональная отношению

h/lc, больше в СА и возрастает с увеличением его концентрации.

На представленных для СБ зависимостях

ширины спектра излучения Dl(Iн) и спектраль-

2 Определяющим для диффузионного режима является диффузионный характер распространения света, при котором в результате многократного рассеяния полностью разрушается когерентность рассеянных волн и поэтому когерентные эффекты несущественны.

r, отн. ед.

1 (а)
1 2

0
540 560 580 l, нм 600

Dl, нм
40

1

(б) rmax, E, усл. ед.
2 3

20 I

III

II

0 IE-3

0,01

0,1 1
Iн, МВт/мм2

Рис. 1. а – спектры излучения родамина 6Ж в  многократно рассеивающей суспензии частиц боразона в ПВА в режиме спонтанного излучения (1) и ХГ (2); б  – экспериментальные зависимости ширины спектра (1), спектральной плотности энергии в максимуме спектра (2) и интегральной по спектру энергии (3) излучения родамина 6Ж в многократно рассеивающей суспензии частиц боразона в ПВА от интенсивности накачки.

ной плотности энергии rмакс  =  (dE/dl)lмакс  = =  f(Iн) в максимуме этого спектра от интенсивности излучения накачки Iн четко выделяются две точки перегиба (рис.  1б). Они делят эти кривые на три участка, соответствующих трем режимам излучения образца: I  – спонтанное излучение (люминесценция), II  – усиленное спонтанное излучение (УСИ), III  – хаотическая генерация. Переход между этими режимами достаточно плавный, поэтому пороговая интенсивность накачки для возникновения УСИ и ХГ определяется условно по точкам перегиба этих кривых. Рисунок  1а показывает различие спектров излучения в  крайних режимах – спонтанного излучения и ХГ.

“Оптический журнал”, 79, 9, 2012

31

Dl, отн. ед.
1

0 0,01

rmax, усл. ед.
1,0

1

(а) rmax, E, усл. ед. 2 3
1′′
1′ 1
0,1 1
Iн, МВт/мм2 (б) Dl, нм
2
12

0,3 8
50 200 350 500 650
h, мкм
Рис.  2. а  – экспериментальные зависимости ширины (1, 1′, 1′′), спектральной плотности энергии в максимуме (2) и интегральной по спектру энергии (3) излучения родамина 6Ж в многократно рассеивающей суспензии частиц аэросила в ПВА от интенсивности накачки при концентрации частиц 3×1011  (1), 1,7×1011  (1′), 5,3×1010  (1′′)  см–3; б  – экспериментальные зависимости спектральной плотности энергии в максимуме спектра (1) и ширины спектра (2) излучения родамина 6Ж в  многократно рассеивающей суспензии частиц аэросила в ПВА от толщины образца.
Энергия, излученная во всем спектре в фиксированный телесный угол (кривая 3), изменяется линейно. Это доказывает, что наблюдающийся нелинейный рост rмакс(Ін) в  диапазоне II обусловлен исключительно ­перераспределением энергии в спектре при п­ ереходе к ХГ.
Кривые Dl(Iн) и rмакс(Ін) образцов СА (рис.  2а) значительно сдвинуты в сторону большей интенсивности накачки относительно аналогичных кривых СБ (рис.  1б). Поэтому в аналогичных диапазонах изменения Iн в СА четко проявляется только первый порог Iп1 (возник-

новения УСИ), а второй порог Iп2 (возникновения ХГ) проявляется только при самых больших концентрациях рассеивающих частиц.
Увеличение порогов Iп1, Iп2 в СА, несмотря на меньшую величину lc в этих суспензиях, указывает на то, что условия возникновения УСИ и ХГ определяются не только длиной свободного пробега и кратностью рассеяния, но и показателем преломления частиц, который влияет на их отражательную способность. Для частиц с диаметром d  >>  l эта способность может быть приближенно оценена коэффициентом отражения Френеля r = (n – nм)2/(n + nм)2, где n и nм  – показатели преломления рассеивающих частиц и матрицы соответственно. С  учетом сказанного, условия возникновения ХГ более объективно отражает величина lr = (rsc)–1, которая учитывает вклад в коэффициент рассеяния только отраженных фотонов, в отличие от величины lc, которая описывает вклад в рассеяние всех фотонов, попавших на частицы.
Отражение, резко меняющее направление упавшего на частицу излучения, приводит к большему удлинению траекторий и, соответственно, более длительному удержанию излучения в активной среде, чем преломление или дифракция на частицах. Исходя из этих соображений, следует ожидать, что уменьшение величины lr способствует понижению порога возникновения ХГ. Оценки подтверждают этот вывод, поскольку соотношение lr, СА  >>  lr, СБ (2,3 и 0,3  см соответственно) коррелирует с соотношением между величинами пороговой интенсивности в СА и СБ: Iп1,  СА  >>  Іп2,  СБ (0,6 и 0,04 МВт/мм2 соответственно).
По своему смыслу длина lr близка к средней транспортной длине lт  =  lc(1  –  〈cоsq〉)–1, соответствующей расстоянию, на котором пучок света теряет свою направленность [13]. Входящий в  эту формулу средний косинус угла рассеяния 〈cоsq〉 определяется дифференциальным сечением рассеяния частиц. Для шарообразных диэлектрических частиц он может быть вычислен по формулам Ми. Однако для частиц неправильной формы эти формулы могут дать лишь оценочный результат, поэтому для оценок использование lr (вместо lт) является вполне оправданным.
На формирование ХГ существенно влияет отражение излучения на границах образца. Как видно из рис.  2б, зависимость rмакс(h) является немонотонной: спектральная плотность энергии в максимуме спектра сначала растет

32 “Оптический журнал”, 79, 9, 2012

при уменьшении толщины образца до некоторой величины, а затем резко падает. Возрастание rмакс сопровождается сужением спектра излучения ХГ. Эти закономерности могут быть объяснены постепенной заменой менее эффективной положительной обратной связи за счет рассеяния в объеме образца на более эффективную связь, обусловленную отражением от тыльной (по отношению к входной грани, на которую падает накачка) границы образца [14]. Коэффициент этого отражения соответствует среднему коэффициенту отражения рассеянного излучения на границе образца с воздухом, который для ПВА составляет 57% [15]. Это отражение, возвращая излучение ХГ назад, в  активную часть образца, увеличивает его эффективное усиление, что и приводит к сужению спектра и росту rмакс. При уменьшении толщины образца меньше толщины активной (возбужденной) области возрастание rмакс сменяется его падением, поскольку падает объем излучающей области образца.
2. Методика моделирования
Для моделирования ХГ методом Монте-Карло была использована модель распространения излучения в активной хаотической среде, разработанная на основе описанной в [9] модели. Среда состоит из одинаковых сферических частиц (рассеивающих центров), хаотически распределенных в среде со спектрально зависимыми коэффициентами усиления k(l, x, y, z) и  поглощения a(l, x, y, z) (рис.  3а). Диаметр ч­ астиц d принимался бо′льшим, чем длина волны излучения l, коэффициент отражения (рассеяния) частиц был принят равным r  =  1, а  угловая зависимость вероятности рассеяния на частице задавалась как P(q)  ~  (1–cоsq). Коэффициенты k и a зависят от пространственных координат в соответствии с заданным распределением населенности энергетических уровней активных молекул, обуславливающих усиление и поглощение в среде. Предполагалось, что таковыми являются молекулы родамина 6Ж.
Коэффициент отражения и индикатриса рассеяния были выбраны исходя из следующих соображений. При включении в индикатрису рассеяния отраженного и преломленного компонентов коэффициент отражения непоглощающей диэлектрической частицы может быть принят за единицу, а вид индикатрисы определяется суммарной угловой зависимостью

коэффициентов преломления и отражения для больших частиц (с  диаметром d  >  l) или теорией Ми для частиц малого диаметра (d  ~  l). Однако при многократном рассеянии, как показали результаты моделирования, статистические характеристики рассеянного излучения слабо зависят от вида угловой зависимости вероятности рассеяния P(q). Поэтому, для упрощения расчетов, использовалась зависимость P(q)  ~  (1  –  cоsq), которая соответствует индикатрисе отражения сферической частицы с диаметром, бо′льшим длины волны падающего излучения, и коэффициентом отражения ­равным единице.
В соответствии с идеологией ММК, для моделирования начального (затравочного для развития ХГ) спонтанного излучения в активной рассеивающей среде генерировался начальный ансамбль фотонов, распределенных в  пространстве и спектре в соответствии с пространственным распределением возбужденных молекул и контуром люминесценции красителя. Во времени фотоны ансамбля составляли прямоугольный импульс длительностью tи. Количество фотонов в ансамбле определялось требуемой точностью вычислений и составляло 2×103–5×105. Процесс распространения гене-

ha Накачка

da dм

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

“Оптический журнал”, 79, 9, 2012

33

рированных фотонов в среде рассматривался, с  корпускулярной точки зрения, как марковская цепь случайных столкновений фотонов с  рассеивающими частицами. Направление фотонов после столкновения определялось как случайная величина с угловым распределением, соответствующим индикатрисе рассеяния. При расчете точки следующего рассеяния ­предполагалось, что наиболее вероятные положения рассеивающих частиц по отношению к  точке последнего рассеяния находятся на сфере с  радиусом, равным среднему расстоянию между частицами x–  =  c–1/3, где с  – количество рассеивающих частиц в единице объема (рис.  3б). Вероятность отклонения рассеивающих частиц от этой сферы описывалась нормальным распределением Гаусса G(Dx–) с дисперсией, которая принята равной удвоенному диаметру частиц 2d. Зависимость вероятности рассеяния от концентрации частиц определялась “прозрачностью” сферы – долей ее поверхности, занятой проекциями частиц. На границе образца фотон отражался с вероятностью r, соответствующей среднему коэффициенту отражения рассеянного излучения на границе ПВА с воздухом [15].
На отрезке пути между двумя последовательными столкновениями вес фотона p изменялся вследствие вынужденного излучения и поглощения в молекулах красителя. Для учета обратного влияния излучения на плотность инверсии населенности уровней красителя (эффект насыщения) образец разбивался на 105–106 элементарных объемов (в зависимости от требуемой точности). Изменение веса фотона при прохождении такого элементарного объема приводило к такому же по величине, но противоположному по знаку изменению плотности инверсии в этом объеме. Начальное распределение инверсии в среде задавалось законом Б­ угера с учетом просветления красителя в слое, прилегающем к входной грани образца:

D(h)

=

ïïïïìîíDDììààêêññ

, h < hï exp(-a(h

-



)),

h ³ hï ,

(1)

где hп – толщина просветленного слоя, определяемая из уравнения

Iíàñ = Iíexp(-aðhï ),

(2)

Iнас – интенсивность насыщения красителя в образце (примерно 5  МВт см–2 [3]), a  =  aр  +  aпогл – полный коэффициент ослабления, включающий коэффициенты рассеяния aр и поглощения aпогл, Iн  – интенсивность на-

качки на входной границе образца, Dмакс  =  c  – плотность инверсии в просветленной части образца.
Таким образом рассчитывались веса всех фотонов ансамбля, приобретаемых ими в процессе усиления и поглощения. Спектр и мощность излучения ХГ в зависимости от интенсивности накачки, концентрации рассеивающих частиц и красителя, размеров образца, поперечных размеров активной области (пучка накачки) определялись суммированием весов фотонов каждой энергии (длины волны излучения), рассчитанных при этих параметрах.
Для определения динамики излучения рассчитывалось пространственное распределение фотонов накачки в образце после каждого акта рассеяния фотона. При усреднении по результатам многих (104) испытаний это распределение соответствовало моментам времени, кратным среднему времени свободного пролета фотона tc = lc/v (v – скорость света в среде).
Процедура вычислений включала следующие шаги. Вначале генерировался ансамбль начальных затравочных фотонов спонтанного излучения (распределенных в пространстве и  спектре аналогично предыдущему), общее число которых соответствовало плотности фотонов накачки в образце в начальный момент времени. Затем рассчитывались положения и  веса этих фотонов на момент первого tp1 рассеяния, и генерировался новый ансамбль затравочных фотонов в соответствии с плотностью фотонов накачки в этот момент времени, после чего аналогичные расчеты проводились для фотонов обоих ансамблей. С  момента времени tр2 расчеты проводились для трех ансамблей и т. д. Динамика изменения спектра и энергии излучения ХГ определялась путем усреднения данных всех испытаний, полученных после каждого столкновения фотонов ХГ с частицами. При достаточно большом количестве испытаний полученные результаты соответствовали пространственному и спектральному распределению фотонов ХГ в моменты времени кратные tc.
3. Результаты расчетов и их обсуждение
На рис. 4 приведены расчетные зависимости Dl(Iн) и rмакс(Ін) для суспензии частиц диаметром 3  мкм в растворе родамина 6Ж. В  двойном логарифмическом масштабе (рис.  4а, вставка) характер этих зависимостей хорошо коррелирует с экспериментальными кривыми

34 “Оптический журнал”, 79, 9, 2012

(рис.  1, 2). На расчетных кривых также выделяются три области, отличающиеся разным характером зависимости от накачки, которые относятся к различным режимам излучения: спонтанного излучения, УСИ и ХГ. Эти области характеризуются соответственно не зависящим от накачки широким спектром и линейным ростом rмакс; быстрым сужением спектра, сопровождающимся нелинейным ростом rмакс; узким (примерно в шесть раз у′же, чем в спонтанном режиме) спектром, слабо зависящим от Iн. Такой характер эти зависимости проявляют лишь при достаточно высоких концентрациях рассеивающих частиц (превышающих с = 5×108 см–3), при которых в среде превалирует многократное рассеяние.
Повышение концентрации частиц увеличивает в области УСИ темп роста спектральной плотности энергии в максимуме спектра излучения (рис.  4а) и темп сужения спектра (рис.  4б) с увеличением интенсивности накачки. Это сопровождается также снижением порога ХГ. Эти изменения коррелируют с увеличением кратности рассеяния q и связанным с этим ростом длины пробега излучения в среде [10, 11].
Согласно расчетам кратность рассеяния является приблизительно степенной функцией от

концентрации рассеивающих частиц q  ~  c1,9, поэтому темп снижения порога с увеличением c возрастает.
Другим важным параметром, влияющим на порог ХГ, является диаметр пучка накачки, при уменьшении которого ниже определенного значения наблюдается резкое увеличение этого порога. Эту зависимость можно условно назвать поперечным размерным эффектом. Он  может быть объяснен наличием в многократно рассеивающей среде некоторой области, в которой преимущественно локализуются траектории излучения ХГ [9, 10]. Однако адекватнее считать, что те траектории, в которых излучение взаимосвязано благодаря наличию общих участков, могут быть сгруппированы в  пучки и сопоставлены модам ХЛ, которые также локализуются в ограниченной области пространства. В  такой интерпретации упомянутое резкое возрастание пороговой интенсивности накачки возникает тогда, когда диаметр активной области dа (определяемый диаметром пучка накачки) становится меньше диаметра dм области локализации мод ХГ и часть объема моды оказывается за пределами возбужденной области, где преобладает поглощение (рис. 3а). Это приводит к падению эффективного усиления мод, для компенсации которого необходи-

rmax, усл. ед.
rmax, усл. ед. 4,0×109

(а)
Dl, отн. ед. 4 3 2 16

4

4

Dl, отн. ед.
1,0

(б)

2,0×109 1Е7 0,1

2
5 0
1 Iн, МВт/мм2

3 2
1

0,5

1
2 3 4

0,1

1 Iн, МВт/мм2

0,1

1 Iн, МВт/мм2

Рис.  4. а  – расчетные зависимости спектральной плотности энергии излучения rmax в модельной многократно рассеивающей среде от интенсивности накачки при различной концентрации рассеивающих частиц: 109 (1), 5×109 (2), 1010 (3), 3×1010 (4) см–3. На вставке: сравнение этих зависимостей в двойном логарифмическом масштабе с расчетной зависимостью ширины спектра излучения (5) от интенсивности накачки при максимальной концентрации частиц 3×1010 см–3; б – участки расчетных зависимостей ширины спектра излучения от интенсивности накачки при концентрации рассеивающих частиц 109 (1), 5×109 (2), 1010 (3), 3×1010 (4) см–3, нормированные в точке с минимальной накачкой.

“Оптический журнал”, 79, 9, 2012

35

мо поднять интенсивность накачки, что проявляется как увеличение порога ХГ.
Размерный эффект должен проявляться и в продольном направлении (вдоль пучка накачки) при изменении толщины образца. Однако его проявление может маскироваться отражением от границ образца, которое, возвращая излучение назад, увеличивает положительную обратную связь и эффективное усиление, испытываемое этим излучением в образце. Поэтому при утоньшении образца, когда его тыльная поверхность приближается к границе области локализации мод ХГ, вместо ожидаемого ухудшения условий ХГ (в  частности увеличения пороговой интенсивности накачки), может наблюдаться их улучшение, что и проявляется экспериментально (рис. 2б).
Представленные на рис. 5 расчетные данные подтверждают эти выводы. Видно, что с уменьшением толщины образца вплоть до 200  мкм количество фотонов Nи (пропорциональное энергии генерированного излучения) в образце возрастает, а спектр ХГ сужается. Резкое падение количества фотонов при уменьшении h ниже этого значения не может быть связано с проявлением пространственного эффекта, поскольку отражение от тыльной грани образца препятствует выходу излучения моды за пределы активной области, как в поперечном эффекте. Поэтому сокращение толщины образца до h    hа), либо из-за уменьшения излучающего объема образца (h tп, где tп  =  h/v  – время пролета фотоном образца (v  – скорость света в среде), люминесценция является более инерционным процессом, чем

“Оптический журнал”, 79, 9, 2012

37

влияние рассеяния на накопление излучения в образце. Поэтому динамика развития излучения ХГ, обусловленная многократным рассеянием и вынужденным излучением, в реальных условиях является скрытой. В  частности, при длительных импульсах накачки (tн  >>  tз) это приведет к зависимости Nи(t), практически повторяющей импульс накачки Nн(t), а при коротких (tн