Численное решение нелинейного уравнения Шредингера в радиально-симметричном случае Текст научной статьи по специальности «Физика»

ЧИСЛЕННОЕ РЕШЕНИЕ НЕЛИНЕЙНОГО УРАВНЕНИЯ ШРЕДИНГЕРА В РАДИАЛЬНО-СИММЕТРИЧНОМ СЛУЧАЕ

А.А. Дегтярев, А.Е Деркач.

Самарский государственный аэрокосмический университет

Аннотация

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

Введение

Как известно [1], нелинейное уравнение Шредингера является частным случаем волнового уравнения в параболическом приближении и с учетом эффекта самовоздействия. Эффект самовоздействия проявляется при распространении оптического излучения в средах с кубичной нелинейностью (поляризация пропорциональна напряженности электрического поля в третьей степени).

С учетом этого эффекта уравнение Шредингера в радиально-симметричном случае запишется как [2]

где Лr =

dU i -+ —

dz Ikn, 0 < r < R,

dr2

-Л rU —

kn2

0

2n

0

U2

U = 0,

(1)

0 < z < L,

Л

1

+—

r

д_

dr

оператор Лапласа в ради-

ально-симметричном случае,

и - напряженность электрического поля, к - волновой вектор, п0 - показатель преломления среды, пнл - изменение показателя преломления под действием поля распространяющейся волны.

Уравнение (1) необходимо дополнить граничным и начальным условиями

К =0 =Ф(г)

Иг=я = о

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

Без учета нелинейных эффектов уравнение Шредингера принимает вид

(2)

dU + _J_

dz 2kn

■Л rU = 0

(3)

0

В дальнейшем уравнение (1) будем называть нелинейным уравнением Шредингера, а уравнение (3) - линейным уравнением.

1. Расчетная конечно-разностная схема

с итерационным уточнением Для численного решения системы (1)-(2) построим консервативную разностную схему [2]. Определим сетку

zk = k ■ hy, k = 0,K -1, hy = L /K

гг = (г + 0,5) • Нг, г = 0,1 -1, Нг = К /(I + 0,5) где , Нг - шаги по переменной 7, описывающей направление распространения волны, и по радиальной переменной соответственно.

Оператор Лапласа по переменной г аппроксимируем следующим разностным оператором [3]:

1 / ч А (г + Нг) - А(г) Л гА(г) =-(г + 0.5НГ))-г--— -

rh

hr

- (г - 0.5Иг )А(Г) - А(Г - Иг),

гкг Иг

причем Лг =ЛгА + о(кГ: / г) для достаточно гладкой функции А. Задаче (1)-(2) поставим в соответствие следующую двухслойную нелинейную разностную схему:

a(z + hz, r)- a(z, r)

kn2

2n0

+ Л rA( r)-

2kn0

-\A(z, r)|2 A(z, r) = 0

(4)

а(0, г ) = р(г) а(, К )= 0

А(т,, г) = 0.5(а^ + Н2, г) + а(, г)).

Система разностных уравнений (4) является нелинейной. Для нахождения неизвестной функции а(р + Нг, г) можно использовать итерационный метод последовательных приближений в сочетании с разностной прогонкой по г:

a (z + hz, r)- a(z, r) i s+},

—^-^—+-Лr A (z, r)-

2kn0

hz

kn„,

A( z, r)

A(z, r) = 0,

2n(

a(z + hz, R )= 0, A(z, r)= 0.5| a(z + hz, r)+ a(z, r)

(5)

i(z + hz, r)= u(z, r) s = 0,1,..

Теоретические исследования сходимости схем вида (5) приведены в работе [2].

r

h

z

2

Необходимо также отметить, что так как данная схема является консервативной, она - аналог интегральной формы закона сохранения электромагнитной энергии в среде.

2. Вычислительные эксперименты для среды с постоянным показателем преломления

Рассмотрим сначала аналитическое решение линейного уравнения, то есть пнл = 0. В этом случае решение задачи (1) запишем в виде

и(г)=£Ст ехр

( 2 Л -1 ^ а

V ^2 у

J о

Я

где а = —

2кп0

Ст =№гУт [^г^Г ,

Я

а /ит - нули функции Бесселя нулевого порядка.

Для простоты выберем начальное условие в виде функции Бесселя нулевого порядка

г ) = J о

Я

тогда, используя свойство ортонормированности Бесселевых функций, получаем

Со = 1, Ст = 0, т = 1,2..., в результате решение системы (1) - (2) будет иметь вид

и ( г ) =

ехр

( 2 Л ■ ^

-1-а

V Я2 у

J о

Я

(6)

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

Пучок, описываемый формулой (6), является стационарным, так как и (г,= и (г,0)2, то есть распределение интенсивности в плоскости ОХУ не зависит от расстояния г. Это распределение изображено на рис. 1.

Рис. 1. Распределение интенсивности волны в линейной среде.

На рис. 2-3 приведены результаты численного решения линейного уравнения Шредингера на различных расстояниях от входа в среду.1

В расчетах были использованы следующие значения параметров: Я=50 мкм;

Х=0,63 мкм, причем Х=2л/к; П0=1.

Рис. 2. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 1м, пг = 1000, пг = 100).

Рис.3. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 10м, пг = 1000, пг = 100).

Как видно из приведенных графиков результат численного решения уравнения Шредингера полностью совпадает с результатом аналитического решения. Это подтверждает правильность выбора консервативной разностной схемы (5) для вычислительных экспериментов с уравнением Шредингера.

3.2. Численное решение нелинейного уравнения Шредингера

Самофокусировка пучка

На рис.4-6 приведены результаты численного решения нелинейного уравнения Шредингера со

следующими параметрами:

Я=50 мкм;

Х=0,63 мкм, причем Х=2л/к;

В дальнейшем использованы следующие обозначения: пг - количество интервалов разбиений по радиусу; пг - количество интервалов разбиений по оси 2. Белым цветом изображается график решения линейного уравнения (3), а серым цветом - график численного решения уравнения Шредингера (1)-(2).

Г

т=0

1

г

г

n0 = 1;

пнл 2=0,001, причем пнл >0.

В работе [1] показано, что при пнл >0 происхо-пучка на расстоянии

дит фокусировка

f = R

п

2

2п2шф(0)

от входа в нелинейную среду.

Самодефокусировка пучка На рис. 7-10 приведены результаты численного решения нелинейного уравнения Шредингера со следующими параметрами: R = 50 мкм;

X = 0,63 мкм, причем X=2n/k; п0=1;

пнл 2=0,001, причем пнл<0.

Рис. 4. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 500мкм, пг = 1000, пг = 100).

Подставляя экспериментальные значения параметров, получаем, что / « 1100 мкм. Вычислительные эксперименты (рис. 3-5) показывают, что фокусировка пучка происходит на расстоянии порядка 800-900 мкм. Таким образом, результаты численного решения достаточно близки к теоретическим оценкам, что еще раз подтверждает практическую пригодность выбранной разностной схемы.

Рис. 5. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 700мкм, пг = 1000, пг = 100).

Рис. 6. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 750мкм, пг = 1000, пг = 100).

Рис. 7. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 500мкм, пг = 1000, пг = 100).

Рис. 8. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 1мм, пг = 1000, пг = 100).

Рис. 9. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 1,5мм, пг = 1000, пг = 100).

В работе [1] показано, что при пнл <0 происходит самодефокусировка пучка в нелинейной среде, то есть смещение энергии пучка в периферийную зону. Приведенные выше графики численного моделирования демонстрируют этот факт.

Рис. 10. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 2,5мм, п2 = 1000, пг = 100).

3. Вычислительные эксперименты для среды с показателем преломления, распределенным по параболическому закону Пусть показатель преломления зависит от радиуса следующим образом

(

2 2 П0 = П

1 - 2Д| -

а

2 Л

(7)

причем пнл = 0.

Будем искать решение уравнения в виде

И ( г ) = А(г )в(). (8)

Подставляя (8) в уравнение (1), получаем фак-торизованное уравнение

Ш дВ 1

B dz n0A

Л rA = -f,

где f - произвольная константа.

Г f

Таким образом, B(z) = C exp| i — z

V 2k

(9)

(10)

Представим A(r) в виде произведения двух функций

A(r ) = X (r )exp

Г 2 Л

Г

(11)

Осуществим еще одну замену переменных

2а

о0

Тогда получим

д2 X

а

да'2

„2 Г

h Л дХ + (1 -а)-—т +

да

2

4

Л

(16)

о

X = 0

о У

Решением уравнения (16) являются многочлены Лагерра

X = Ls (а' ) = Ls причем

П 2 1

s =-<о0--

8 0 2

Г 2r 2 Л

(17)

(18)

Таким образом, с учетом (8), (10) и (11) получаем искомое решение

U(z,r) = Cexp| i2kz L

Г 2r 2 Л

x exp

2Л

r

(19)

3.1. Результаты численного решения линейного уравнения Шредингера

Для вычислительных экспериментов выберем, например, моду Гаусса-Лагерра с номером s=10 Тогда из (18) и (14) получаем, что

f =

422 Д

о0 =-

21Д

(20)

(21)

Для того чтобы обеспечить выполнение граничного условия, выберем аргумент многочлена Лагерра, равный десятому нулю функции.

2R2

(22)

После подстановки (11) и (10) в (8) получаем

д2 X 4r dX 1 dX +

дг2

<0 дг r дг

4

<0

X+

4r2

(12)

+—X + >0 X = 0 <0

Осуществим замену переменных а = r . Тогда (12) преобразуется к виду

д2 X да 2

1-

2а

2

X

да

+

Г а аД П Л о 4а2

0

X = 0.

fn1

1

X+

j0

1 ДП

Пусть — =

< 04 4a 2

(13)

где Ш10=29,9206970123 - десятый нуль многочлена Гаусса-Лагерра.

На рис. 11 изображено распределение интенсивности в пучке, описываемом формулой (19).

Рис. 11. Распределение интенсивности волны в линейной среде.

а =

1

+

2

na

2

a

= m

10

V <0 у

<

+

4

+

На рис. 12-13 приведены результаты численного моделирования линейного уравнения со следующими параметрами: Я = 50 мкм;

X = 0,63 мкм, причем Х=2п/к; п1 = 1.

Рис.12. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 1м, пг = 1000, пг = 100).

Рис. 13. Распределение интенсивности волны на расстоянии Ь от входа в среду (Ь = 10м, пг = 1000, пг = 100).

Как видно из приведенных графиков, результат численного решения уравнения Шредингера полностью совпадает с результатом аналитического решения.

Заключение Итак, в результате проведения вычислительных экспериментов получено подтверждение сходимости разностной схемы (5), использованной для решения нелинейного уравнения Шредингера. Процесс итерационного уточнения обнаружил быструю сходимость: для уточнения решения до величин порядка 10-4 требуются 2-3 итерации.

Проведенные исследования подтвердили высокую точность численного решения задачи. Так для модели среды, не учитывающей эффект самовоздействия, численное решение практически совпадает с аналитическим, рассчитанным по формуле (6).

Для среды с "параболическим профилем" показатели преломления моды Гаусса-Лагерра, рассчитанные с помощью схемы (5), практически не претерпевают изменений в процессе распространения волны на расстояние порядка десятков метров, что соответствует теории [4].

Для модели среды, учитывающей эффект самовоздействия, численный расчет позволяет наблюдать как явление самофокусировки излучения, так и его самодефокусировки, что спрогнозировано теоретическими оценками, например в работе [1].

Литература

1. Виноградова М.Б., Руденко О.В., Сухоруков А.П.

Теория волн // М., Наука. 1979.

2. Карамзин Ю.Н., Сухоруков А.П., Трофимов В. А. Математическое моделирование в нелинейной оптике // М., Изд-во Моск. ун-та. 1989.

3. Самарский А.А. Теория разностных схем // М., Наука. 1977.

4. Адамс М. Введение в теорию оптических волноводов // М., Мир. 1984.