CC BY
48
РАСЧЕТ СИЛЫ ДАВЛЕНИЯ СВЕТА НА КРУГЛЫЙ ДИЭЛЕКТРИЧЕСКИИ ЦИЛИНДР С ИСПОЛЬЗОВАНИЕМ БЫСТРОГО ИТЕРАТИВНОГО АЛГОРИТМА И НА ОСНОВАНИИ АНАЛИТИЧЕСКОГО РЕШЕНИЯ
Котляр В.В., Налимов А.Г, Личманов М.А.
Институт систем обработки изображений РАН Самарский государственный аэрокосмический университет
Аннотация
Проведено сравнение аналитического решения дифракции непараксиального гауссова пучка на круглом цилиндре с быстрым итеративным алгоритмом расчета. Рассчитаны силы давления, оценена эффективность оптического захвата круглого диэлектрического цилиндра с помощью гауссова пучка с острой фокусировкой в двумерном случае.
Введение
В [1, 2] рассмотрена дифракция непараксиального гауссова пучка на круглом диэлектрическом цилиндре. При этом центр круга может не совпадать с центром перетяжки цилиндрического гауссова пучка. Приведены аналитические формулы, описывающие напряженность электрического и магнитного полей дифракции в виде рядов по цилиндрическим функциям Бесселя и Ханкеля. В [3-7] получены общие выражения для силы давления света на трехмерный диэлектрический объект, основанные на электромагнитном тензоре Максвелла [3]. Сила давления электромагнитной волны на диэлектрический объект выражается через интеграл по поверхности, охватывающей это тело, от квадратичной комбинации проекций векторов напряженности электрического и магнитного полей. Поэтому чтобы рассчитать силу давления света требуется решить задачу дифракции света на объекте.
В [8, 9] был предложен быстрый итеративный метод расчета дифракции двумерной электромагнитной волны с ТЕ- и ТМ-поляризациями на диэлектрическом бесконечно-протяженном цилиндре с произвольной формой сечения. В двухмерной постановке мощность падающего на объект света рассчитывается как линейная плотность, то есть мощность цилиндрического пучка света на единицу длины цилиндра. Также рассчитывается и сила давления света на цилиндр, как сила на единицу длины цилиндра.
В данной работе проведено численное сравнение двух методов расчета дифракции ТЕ-волны на круглом диэлектрическом цилиндре - аналитического [1, 2] и итеративного [8, 9]. Получены условия, при которых оба метода отличаются не более чем на 1%. Кроме этого показано, что при увеличении смещения центра круглого цилиндра из центра перетяжки гауссова пучка растет число ненулевых членов ряда, которые должны быть учтены в аналитическом расчете, чтобы сохранить точность решения.
Также в работе показана возможность оптического захвата цилиндрическим непараксиальным гауссовым пучком с радиусом перетяжки, равным половине длины волны, круглого диэлектрического цилиндра с диаметром, равным длине волны, и диэлектрической проницаемостью, равной 2.
1. Расчет поля дифракции и силы давления света на диэлектрический цилиндр
Для расчета сил, действующих на объект со стороны света, необходимо рассчитать все компоненты электрического и магнитного полей дифракции лазерного излучения на этом объекте. В частности для двумерного случая с ТЕ-поляризацией необходимо рассчитать электрическую компоненту Ех и магнитные компоненты НУ, Н2.
у 2
Далее рассмотрим быстрый итеративный алгоритм расчета поля дифракции и расчет аналитическим методом с разложением по цилиндрическим волнам. Рассчитанные поля дифракции затем подставим в расчет силы давления света на двумерный объект.
1.1 Итеративный алгоритм
В [8] описан итеративный метод расчета дифракции произвольной электромагнитной волны, удовлетворяющей уравнению Гельмгольца, на цилиндрических диэлектрических микрообъектах с произвольной формой сечения (ТЕ-поляризация). Ядром итеративного алгоритма является уравнение Фредгольма второго рода:
Еи+1 (2, У) = уЕ0, (2, У) + (1 - у)Еп (2, У) + +УРДО Еп & П)Н01) (ц/ (2-О2 + (У-П)2) , (1)
где (2, у) е X ; Еп+1 и Еп - амплитуда электрического поля внутри области X на (п+1)-ом и п-ом шагах итераций; у - постоянная релаксации алгоритма, которая регулирует скорость его сходимости, _ 1к2 (е-1) 1 2п
р =-, к = — - волновое число света, е -
4 X
диэлектрическая проницаемость цилиндра, (§, п) е X - декартовы координаты внутри объекта,
Н0(1)(х) - функции Ханкеля 1-го рода нулевого порядка. В качестве начального приближения можно брать падающее поле:
Е,( х, у) = Ех ((х, у). (2)
Данный алгоритм сходится при условии кёе< 4п, где ё - диаметр сечения цилиндра. Для моделирования оптического захвата цилиндриче-
ских микрообъектов можно воспользоваться аналитическим решением дифракции гауссова пучка.
1.2. Аналитическое решение задачи дифракции
Рис. 1. Схема падения гауссова пучка с фокусом в точке (- 2о, Уо) на круглый цилиндр с центром в точке (0;0) Следуя [2], рассмотрим дифракцию двумерного непараксиального гауссова пучка на круглом однородном цилиндре (рис. 1). Для случая ТЕ-поляризации, когда (Ех, Ну, Н*) - отличны от нуля, напряженность электрического поля для непараксиального гауссова пучка можно записать в виде:
Ео х (г, ф) = Ео x (кг у
(3)
С =
Ю0л/п Г -—] ехр
к2®2д2
1 - д2 - гкдуо - т атс8т д
dд.
(4)
Напряженность электрического поля вне цилиндра равна сумме рассеянного и падающего полей:
Г
Ех = Ео X (кг У- +
Г
+Ео X Г«„СНЧкг) е'пф (5)
где
а„ = (к1 ^ (к,К) ¿п (К - ып Цф) Гп (кК)) / / ( (к1 К) НП1'> (кК) - Ып к К) Н'пт (кК)) , (6)
где к1 = —^ё^ - волновое число в среде с диэлек-
X
трической проницаемостью 81.
Связь между проекциями Н , Н* и Ех следует из уравнений Максвелла:
Н =-1 дЕк
у к д*
н2 = 1 дЕ±
к ду
(7)
1.3 Сила давления света на объект
Сила, действующая на объект со стороны падающего поля, рассчитанная в общем случае с помощью электромагнитного тензора Максвелла [9] в 2Б случае, может быть вычислена интегрированием по контуру £ вокруг объекта:
Ех = 0,
¥У = 18о | {2 Ну |2-81 |Ех|2-| Н*|2
+Яе (НуН*) dy},
^ = 2ео|{2[IН* Г-81 Ех|2-|Ну
+Яе (Н*Н*) dz},
dz +
dy +
(9)
где 8о - диэлектрическая постоянная,
8о = 8.85.1о-12^-
о Н • м2
(в системе СИ).
2. Численное моделирование Исследуем работоспособность итеративного и аналитического методов расчета дифракции электромагнитной волны на микрообъектах. Далее, на основании результата вычислений поля дифракции электромагнитной волны на микрообъекте, рассчитаем силу давления света на этот микрообъект.
2.1. Сравнение итеративного и аналитического расчетов.
а) б)
Рис. 2. а) Амплитуда Ех поля дифракции непараксиального гауссова пучка на цилиндре, б) модуль проекции вектора Умова-Пойнтинга на ось распространения света 2.
Параметры: длина волны X =1 мкм, диаметр круглого цилиндра ¿1=1 мкм, диэлектрическая проницаемость цилиндра 8 =2, среды 81 =1, все поле 1,5x1,5 мкм, цилиндр помещен в центр перетяжки с радиусом ст =0,5 мкм
На рис. 2 представлено поле дифракции, рассчитанное по формуле (5). СКО между данным решением и расчетом картины дифракции при помощи алгоритма, описанного в [8], составляет о,о8% при параметрах: длина волны X = 1 мкм, диаметр круглого цилиндра ¿=1 мкм, диэлектрическая проницаемость цилиндра 8 =2, среды 81 =1, все поле 1,5х1,5 мкм, 256х256 отсчетов, объект помещен в центр перетяжки с радиусом ст =о,5 мкм, количество коэффициентов ряда Л/=8. Вектор Умова-Пойнтинга рассчитывается по формуле:
£ = 8оС2 [Е X Н ] (1о)
а его проекция на ось распространения света 2 будет равна:
£* =8оС2ЕхНу.
(11)
2
77 = — '-«
На рис. 3 приведены те же поля для цилиндра, смещенного на 0,25 мкм (вверх) по оси Y из перетяжки. По картине модуля проекции вектора Умова-Пойнтинга на ось распространения света хорошо заметно, что энергия в пучке отклоняется в сторону смещения цилиндра.
а) б)
Рис. 3. а) амплитуда Ex поля дифракции непараксиального гауссова пучка на цилиндре, б) модуль проекции вектора Умова-Пойнтинга на ось распространения света Z. Параметры те же, что и на рис. 2, цилиндр смещен на L=0,25 мкм по оси У
Точность расчета поля дифракции аналитическим методом зависит от количества коэффициентов ряда Сп (4). Зависимость невязки 5 между амплитудами электрических полей, полученных при помощи итеративного алгоритма и по формуле (5), от количества коэффициентов показана на рис. 4, рассчитываемая по формуле:
X (Ех (п, т) - Еп (п, т))2
(12)
Количество коэффициентее;
Рис. 4. Зависимость невязки 5 между амплитудами поля Ех, рассчитанного итеративным и аналитическим методами, от количества коэффициента Сп ряда (5)
Количество коэффициентов ряда Сп изменяется
от -3..3 до -9..9, остальные параметры моделирования те же, что и для рис. 2.
Диапазон параметров, при которых поле дифракции рассчитывается корректно у аналитического метода шире, чем у итеративного алгоритма. Однако аналитический метод дает погрешность при расчете дифракции гауссова пучка на цилиндре, находящегося в теневой области. На рис. 5 представлен график зависимости невязки 5 между итера-
тивным и аналитическим методами расчета дифракции от смещения круглого цилиндра вдоль вертикальной оси У (рис. 1) по линии, проходящей точно через перетяжку (7.=0).
а)
б)
Смешение 1_, мкм
Рис. 5. Зависимость невязки 5 между амплитудами
электрических компонент электромагнитной волны, рассчитанных итеративным и аналитическим методами.
В аналитическом методе бралось: а) N=7 коэффициентов; б) N=9 коэффициентов
Как можно увидеть из рис. 5, невязка сильно возрастает при смещении цилиндра из фокуса в теневую область. Однако это можно исправить взятием большего количества коэффициентов Сп. Такая необходимость взятия большего количества коэффициентов Сп возникает из-за смещения значимых для расчета больших по модулю коэффициентов Сп от 0-го номера в сторону больших номеров, что проиллюстрировано на рис. 6.
На рис. 6а, б изображена зависимость модуля коэффициентов Сп от номера п. Как видно из графиков, при взятии того же количества коэффициентов и смещении цилиндра в вычислении принимают участие не все значимые коэффициенты, что дает увеличение погрешности расчетов.
Зависимость невязки 5 между картинами дифракции, рассчитанными итеративным и аналитическим методами, от смещения вдоль оси 2 дает иную картину. Так, на рис. 7 можно увидеть, что невязка незначительно возрастает при смещении цилиндра из фокуса вдоль оси распространения света на большее смещение, чем вдоль оси У. Параметры вычислительного эксперимента те же, что и для рис. 2.
Это объясняется другой зависимостью коэффициентов Сп от смещения, отображенной на рис. 8.
а)
б)
Рис. 6. Зависимость модуля коэффициента Сп от номера п: а) цилиндр совпадает с фокусом гауссова луча, б) цилиндр смещен из фокуса гауссова пучка на Ь=1 мкм по оси У
Рис. 7. Зависимость СКО между амплитудой Ех поля дифракции, рассчитанного итеративным и аналитическим методами, от смещения Ь цилиндра вдоль оси 2 через фокус при У=0
а коэффициента (
Рис. 8. Зависимость модуля коэффициентов Сп от
номера п. Перетяжка смещена от цилиндра по оси 2 на расстояние Ь=-1,56 мкм
Как видно из рис. 8, ширина полосы значимых коэффициентов Сп расширяется при смещении пе-
ретяжки из центра координат (и от центра цилиндра), и, чтобы это компенсировать, необходимо брать тем больше коэффициентов Cn, чем дальше смещена перетяжка гауссова пучка. В данном примере (рис. 8) взятие коэффициентов Cn от -9 до 9 дает невязку 5 между амплитудами Ex, рассчитанных аналитическим и итеративным методами, порядка 1% (рис. 7).
Взятие большего количества коэффициентов Cn сопряжено с увеличением вычислительной сложности и, следовательно, с увеличением временных затрат на расчет поля дифракции. Так, например, расчет поля дифракции размером 5x5 мкм (256x256 отсчетов) при параметрах: длина волны X = 1 мкм, диэлектрическая проницаемость среды е1=1, диэлектрическая проницаемость цилиндра е =2, диаметр цилиндра d=1 мкм (перетяжка расположена в центре координат) аналитическим методом при взятии коэффициентов Cn от -7 до 7 требует 40 секунд для компьютера с ЦП Celeron 1000, СКО между амплитудами Ex, рассчитанными аналитическим и итеративным методами составляет 0,038%. Расчет поля дифракции в том же случае, но при расположении перетяжки в координатах (Z=-1,56 мкм, Y=0 мкм) и взятии коэффициентов Cn от -13 до 13 требует 19 минут 53 секунды, невязка 5 между амплитудами Ex, рассчитанными аналитическим и итеративным методами составляет 0,077%.
2.2 Расчет силы давления
На рис. 9. приведен пример графиков зависимо -сти проекций силы, действующей на круглый цилиндр со стороны гауссова пучка, вдоль осей Y, Z от смещения цилиндра по тем же осям, из которых видна возможность «захвата» объекта лазерным излучением в оптическую ловушку.
а
Рис 9. Зависимость проекций силы, действующей на цилиндр, от его смещения вдоль той же оси: а) вдоль оси 2; б) вдоль оси У
Параметры моделирования: все поле размером 5x5 мкм, 256x256 отсчетов, длина волны X =1 мкм, диэлектрическая проницаемость среды 81=1, диэлектрическая проницаемость цилиндра 8 =1,3, диаметр цилиндра ¿=1,5 мкм, радиус перетяжки в фокусе ст =о,5 мкм. Линия смещения объекта на обоих графиках проходит через центр перетяжки. Из
рис. 9 видно, что сила, действующая на микроци-
и
линдр, равна по порядку величины 10-10 — при
м
мощности лазера гауссова пучка 100
мВт
Из рис. 9 видно, что при смещении в любом направлении из фокуса гауссова пучка сила действует на частицу в направлении, противоположном смещению. Проекция силы вдоль оси 2 складывается из градиентной и рассеивающей сил [10]. Градиентная сила всегда направлена в максимум интенсивности излучения, а рассеивающая - вдоль оси распространения света. Поэтому центр захвата частицы не совпадает с геометрическим центром перетяжки гауссова пучка, и возможность захвата определяется балансом между рассеивающей и градиентной составляющими проекциями силы вдоль оси 2.
На рис. 10 представлен график зависимости эффективности захвата Q от смещения вдоль осей У, 2, рассчитываемой по формуле:
Q = р (13)
пР
где с - скорость света в вакууме, ^ - сила, действующая на объект, п =
те - показатель преломления объекта, Р - мощность излучения лазера.
а)
б)
Рис. 10. Графики зависимости эффективности захвата вдоль оси от смещения вдоль той же оси: а) ось распространения света 2; б) ось У
В случае силы, направленной вдоль рассматриваемой оси, Q>0, и Q<0 в обратном случае, что показывает возможность оптического «захвата» частицы. Физический смысл эффективности захвата -это часть излучения, пошедшего на создание силы, действующей на частицу.
Заключение В работе получены следующие результаты:
- проведено сравнение точности расчета дифракции гауссова пучка на диэлектрических микрообъектах аналитическим и итеративным методами (рис. 4, 5, 7)
- Рассмотрено влияние количества взятых для расчета коэффициентов Cn на точность расчета поля дифракции при смещении объекта из центра перетяжки непараксиального гауссова пучка вдоль осей Y (рис. 5) и Z (рис. 7)
- оценена эффективность захвата диэлектрического цилиндра отдельно для смещения вдоль осей Y и Z (рис. 10) по формуле (10).
Благодарности
Работа поддержана российско-американской
программой «Фундаментальные исследования и
высшее образование» (BRHE), грант CRDF REC-
SA-014-02 и президентским грантом НШ-
1007.2003.01.
Литература
1. Zimmerman E., Dandliner R., Souli N. Scattering of an off-axis Gaussian beam by a dielectric cylinder compared with a rigorous electromagnetic approach. // J. Opt. Soc. Am. A, 1995, v. 12, p. 398-403.
2. Wu Z., Guo L.. Electromagnetic scattering from a multi-layerd cylinder arbitrarily located in a Gaussian beam, a new recursive algorithms. // Progress in electromagnetics research, 1998, PIER, v. 18, p. 317-333.
3. Ландау Л. Д. , Лифшиц Е.М. . Краткий курс теоретической физики. Механика. Электродинамика // Книга 1, М., Наука, 1969.
4. Lock J.A.. Calculation of the radiation trapping force for laser tweezers by use of generalized Lorenc-Mie theory I. Localized model description of an on-axis tightly focused laser beam with spherical aberration. // Appl. Opt., 2004, V. 43, p. 2532-2544.
5. Ganic D. , Gan X. , Gu M.. Exact radiation trapping force calculation based on vectorial diffraction theory. // Opt. Express, 2004, v. 12, no. 12, p. 2670-2675.
6. Nieminen T.A. , Heckenberg N.R. , Rubinstein-Dunlop H. . Computational modeling of optical tweezers. // Proceeings of SPIE, 2004, v. 5514, p. 514-523.
7. Pobre R. , Saloma C. . Radiation forces on nonlinear mi-crosphere by a tightly focused Gaussian beam. // Appl. Opt., 2002, v. 41, no. 36, p. 7694-7701.
8. Котляр В.В., Налимов А.Г., Скиданов Р.В. Быстрый метод расчета дифракции на цилиндрических диэлектрических микрообъектах // Компьютерная оптика, 2004, №25, с. 24-28.
9. Котляр В.В., Налимов А.Г., Скиданов Р.В. Расчет вектора Умова-Пойнтинга и силы давления электромагнитной волны на однородный диэлектрический цилиндр // Известия СНЦ РАН, 2004, т. 6, № 2.
10. Rohrbach A., Stelzer E.H.K. . Optical trapping of a dielectric particles in arbitrary fields. // J. Opr. Soc. Am. A, 2001, v. 18, p. 839-853.
м
