8. Thacher Н.С., Jr. Complete Elliptic Integrals // Communications of the ACM. 1963. Apr. 6. P. 163.
9. Stewart D.E., Ley к Z. Meschach: matrix computations in С // CMA Proceedings. 32. Canberra: Australian National University, Australia, 1993.
Поступила в редакцию 01.10.04
УДК 517.958:535.14
А. Г. Волков, В. А. Трофимов
О РОЛИ СПЕКТРАЛЬНОГО ИНВАРИАНТА
ПРИ КОМПЬЮТЕРНОМ МОДЕЛИРОВАНИИ НЕЛИНЕЙНОГО
РАСПРОСТРАНЕНИЯ ФЕМТОСЕКУНДНЫХ ИМПУЛЬСОВ1
(кафедра вычислительных методов факультета ВМиК)
1. Введение. Как известно [1, 2], анализ распространения фемтосекундного лазерного импульса основан на так называемом комбинированном (обобщенном) нелинейном уравнении Шредингера (КНУШ), которое отличается от традиционного нелинейного уравнения Шредингера (НУШ) наличием производной по времени от нелинейного отклика среды. Ее присутствие предъявляет новые требования к построению разностных схем для данного класса уравнений, которые обусловлены новыми (по сравнению с традиционными НУШ) свойствами КНУШ. Так, производная от нелинейного отклика приводит к зависимости скорости распространения волны от ее интенсивности. Результатом этого является возможность формирования оптических ударных волн [1, 2].
Другие свойства КНУШ удалось выявить, используя подход [3, 4], развиваемый авторами для подобных задач. Так, построенные инварианты показали присутствие "особой" спектральной гармоники в начальном условии для уравнения КНУШ: ее амплитуда изменяется по определенному закону, т.е. в процессе распространения светового импульса имеет место спектральный инвариант. Очевидно, данное свойство необходимо учитывать при построении разностных схем для уравнения КНУШ. Именно влияние точности сохранения спектрального инварианта на получаемое при компьютерном моделировании решение исследуется в настоящей работе. При этом рассмотрение проводится на примере нескольких схем, в том числе и консервативной схемы.
2. Основные уравнения. Процесс распространения фемтосекундного импульса в среде с кубичной нелинейностью с учетом производной по времени от нелинейного отклика среды (дисперсии нелинейности) при отсутствии влияния дифракции оптического излучения (что реализуется либо при его распространении в оптическом волокне, либо при длине среды много меньшей дифракционной длины светового пучка) и с учетом дисперсии второго порядка описывается следующим безразмерным комбинированным нелинейным уравнением Шредингера (КНУШ):
дА д2А д
— + 10— + 1а\А\2А + а1—(\А\2А) = 01 ¿>0, 0 < I < 1и (1)
£ = 0,5(А ехр(шЬ — 1кг)) + к. с. с начальными и граничными условиями
А|г=0 = А0(;), (2)
Здесь £ — напряженность электрического поля; и и к соответственно безразмерная частота и волновое число светового импульса; А(г,£) — нормированная на максимальное значение на входе в нелинейную
1 Работа выполнена при частичной финансовой поддержке РФФИ (грант № 02-01-727).
среду комплексная амплитуда импульса, распространяющегося вдоль координаты z, которая измеряется в единицах характерной длины среды; t — нормированное на длительность импульса время в сопровождающей его системе координат; D — коэффициент, равный отношению длины среды к дисперсионной длине; а — отношение начальной мощности импульса к характерной мощности самовоздействия; 7 — параметр, обратно пропорциональный произведению длительности импульса на его частоту; Lt — безразмерный временной интервал, в пределах которого анализируется процесс распространения импульса; к.с. — комплексно-сопряженное слагаемое.
Для уравнения (1) справедливы следующие инварианты [3]:
Li Lt
Isp{z) = Isp{ISP(z) = j Ae^ dt, IA(z) = j |A|2 dt = const, (3)
о 0
Lt t
/Г i(v-t)
EA* dt = const, E{z,t) = / A(z,r])e di],
о 0
(4)
0
где 1зр{г) — спектральный инвариант; г) — безразмерная энергия импульса; Д(^), ^2(2) — инварианты, которые характеризуют фазовые эффекты при распространении светового импульса.
3. Разностные схемы. Введем в области = (0, Ьг) X (О, Ь^) сетку и = X
= {гт = т/г, т = 0,1,..., Л^, /г = = {¿„ = пт, п = 0,1,..., Л^, т = Lt/Nt}. (5)
Определим сеточную функцию А наши введем также стандартные безындексные обозначения [6]:
0,5 „ 0,5
А = А(гт,гп), А = А(гт+1,гп), А±1 = А(гт, гп±1), А=0,5(А + А), А = 0,5(А±1 + А±1),
= 0,5( А +|А|2),
±1
= 0,5(
A±i
+ |А±1|2), Аии =
ип+1 - 2ип + и П — 1
n=l,...,Nt-l. (6)
Тогда для задачи (1), (2) с учетом введенных обозначений (5), (6) запишем следующие разностные схемы:
0,5
(А - A)/h + iDAB А + ia
0,5
А + «7 F = О,
F =
F =
4-1
+1
0,5
-1
0,5
А
-1
/2т,
2 А+1 + |А+112 А+1- Л'А- \А\2 А ) /2т,
F =
i+i
41
+ \А+1\ А+1 -
А—1 - |A_i | A_i /4т
с начальными и граничными условиями
A(0,i„) = A0(tn), га = 0,..., Ntl A(zm, 0) = A(zm, Nt) = 0, m = О,
(7)
(8)
(9) (10)
(П)
Следует подчеркнуть, что схема (7) с аппроксимацией производной по времени от нелинейного отклика в виде (8) консервативна [5] по записанным выше инвариантам. Схема (7), (10) наиболее близка к консервативной схеме, но она является условно консервативной: ее сеточное решение при достаточно малых шагах совпадает с решением, полученным на основе консервативной схемы. Наихудшей среди записанных схем является схема (7), (9). Она также условно консервативна на ограниченной трассе распространения светового импульса. Однако при любом шаге по времени г существует значение
шага по координате z такое, что при h ^ h* сеточное решение может неограниченно возрастать. Таким образом, для данной схемы имеет место ограничение шага h снизу, что ограничивает точность получаемого решения. Нетрудно видеть, что схема (7) с аппроксимацией производной по времени от нелинейного отклика в виде (8) или (10) второго порядка, а схема (7) с (9) является схемой первого порядка по t из-за аппроксимации производной от нелинейного отклика.
4. Результаты численных экспериментов. Как показал выполненный линейный анализ модуляционной неустойчивости распространения фемтосекундного импульса, которое описывается уравнением (1), существует "особая" гармоника иъ = I/7 в спектральном пространстве, на которой в линейном приближении реализуется конвективная неустойчивость решения [7]: возмущение на этой частоте возрастает вдоль координаты z. Но именно для этой частоты справедлив спектральный инвариант. Следовательно, при адекватном описании свойств исходной дифференциальной задачи сеточное решение не должно неограниченно возрастать на этой частоте. Поэтому прежде всего рассмотрим эволюцию малых возмущений на частоте иъ = I/7 гауссова начального импульса
A0(t) = ехр(-т]2) (1 + 0,1 cos(wBJ7)), т] = 1 - Lt/2. (12)
Для определенности амплитуда возмущения выбрана равной 0,1, тем не менее это не ограничивает общности рассмотрения, так как она определяет лишь длину трассы распространения, на которой возмущение возрастает (убывает). При этом изменение амплитуды возмущения в среде будем оценивать как отношение ее спектральной интенсивности в толще среды к ее значению на входе в среду:
2 00
Ф)= , ,|А"в(г)1м2' [ A(z,t)coS(uBt)dt. (13)
\AuB(z = 0)1 J
—00
Из проведенных расчетов для различных шагов г, h следует, что сохранение спектрального инварианта существенно зависит от частоты возмущения. Если возмущение имеет частоту, отличную от I/7 более чем на 20%, то все рассматриваемые схемы сохраняют спектральный инвариант даже для достаточно больших шагов г и h. Отклонение же частоты возмущения лишь на 10% от значения I/7 на порядок улучшает сохранение спектрального инварианта по сравнению со случаем наличия возмущения на частоте I/7 при неизменных шагах сеток. В поведении схем наблюдаются существенные различия, если возмущение имеет частоту, близкую к I/7. Так, при наличии возмущения на этом интервале частот иъ ~ I/7 спектральный инвариант не изменяется для схем (7), (8) или (10) при достаточно малых шагах г и /г, которые зависят от значения параметра 7. Для схемы же (7) с (9) сохранение инварианта имеет место только при превышении шага h некоторого h* (г). Рассмотрим этот случай более подробно на примере распространения импульса с параметрами D = 1, а = 10, 7 = 0,05, = I/7.
Результаты проведенных расчетов содержатся в таблице и на рис. 1, 2. В таблице представлена эволюция инвариантов вдоль трассы распространения в зависимости от соотношения шагов сеток. Для демонстрации общих закономерностей приведены расчеты при существенно различающихся шагах по времени. Хорошо видно, что при фиксированном шаге по времени с уменьшением шага по координате z условно консервативная схема (7), (9) дает бесконечный рост спектрального инварианта, энергии импульса и его максимальной интенсивности. При этом чем меньше шаг сетки по z, тем меньше трасса распространения, на которой сеточное решение бесконечно возрастает. В качестве подтверждения этого приведен рис. 1. Важно подчеркнуть, что рост максимальной интенсивности имеет место в том же сечении среды, в котором возрастает инкремент возмущения. Следовательно, причиной неограниченного роста интенсивности и инвариантов является рост возмущения, т.е. несохранение спектрального инварианта. Это демонстрирует важность его учета при моделировании распространения фемтосекундного импульса. Имеются отличия и в эволюции инвариантов Ii(z), /2(2:). Так, в проведенных расчетах они сохраняются для консервативной схемы (7), (8) всегда, а для условно-консервативных схем (7), (9) и (10) сохраняются лишь при определенном соотношении шагов сетки.
Из таблицы также следует, что при уменьшении шагов (особенно по времени) консервативная разностная схема для рассматриваемых параметров взаимодействия светового импульса со средой сохраняет спектральный инвариант с высокой точностью. Отсутствие его сохранения консервативной
2Параметры расчетов: D = 1, Lt = 100, а = 10,7 = 0,05, шв = 1/7-
Эволюция инвариантов разностных схем вдоль продольной координаты г в зависимости от соотношения шагов сеток2
т /г г Ы*)
Схема (7), (8), (10) Схема (7), (9)
0,1 0,025 0,0 0,0 0,0 1,26
0,25 0,1 0,11 1,35
0,5 0,16 0,21 1,61
0,75 0,17 оо оо
1,0 0,11
0,1 10"4 0,0 0,0 0,0 1,26
0,25 0,14 0,19 1,36
0,5 0,16 оо оо
0,75 0,0
1,0 0,15
0,01 0,025 0,0 0,0 0,0 1,26
0,25 0,17 0,17 1,27
0,5 0,05 0,05 1,28
0,75 0,16 0,17 1,29
1,0 0,1 0,1 1,3
0,01 10"3 0,0 0,0 0,0 1,26
0,25 0,13 0,13 1,27
0,5 0,18 0,18 1,28
0,75 0,12 0,12 1,29
1,0 0,02 0,02 1,3
0,01 10"4 0,0 0,0 0,0 1,26
0,25 0,03 0,03 1,27
0,5 0,06 0,06 1,28
0,75 0,09 оо оо
1,0 0,11
0,005 10"4 0,0 0,0 0,0 1,26
0,1 0,003 0,003 1,26
0,3 0,01 0,01 1,27
0,5 0,017 0,017 1,27
0,7 0,024 0,024 1,27
0,9 0,031 0,031 1,28
1,0 0,034 0,034 1,28
0,005 10"5 0,0 0,0 0,0 1,26
0,1 0,003 0,003 1,26
0,3 0,009 оо оо
0,5 0,014
0,7 0,02
0,9 0,026
1,0 0,029
схемой при грубых шагах по времени обусловлено, во-первых, применяемым для разрешения нелинейной разностной схемы итерационным процессом. Вторая причина связана с плохой аппроксимацией решения на грубой сетке по времени, которая ограничивает воспроизводимый спектральный состав светового импульса, обогащаемого из-за нелинейности задачи.
При смещении частоты возмущения иъ в область неустойчивости задачи по линейному приближению (например иъ = 5, у = 1) численное решение требует уменьшения шага по к для схем (7), (8) и (10). При этом в отличие от условно-консервативной схемы (7), (10) по консервативной схеме (7), (8) вычисления производятся для более крупного шага по пространству при фиксированном шаге по времени (см. рис. 2,а) (например, при к = 0,0125, г = 0,01 первая схема с начальным распределением (12) дает верный результат, а схема (7), (10) — неограниченный рост решения). Заметим, что из-за грубой сетки г) изменяется примерно на 10-15%. Другие же инварианты для консервативной разностной схемы сохраняются с высокой точностью. Следует подчеркнуть, что полученное сеточное решение сильно "изрезано" на спадающем участке импульса. Тем не менее максимальная интенсивность практически не превосходит значения, вычисленного для более мелкой сетки, и усредненное распределение совпадает с ним. Важно также подчеркнуть, что при отсутствии возмущения для данных шагов обе схемы дают ограниченное решение, но консервативная схема сохраняет спектральный инвариант на 3 порядка лучше, чем схема (7), (10). Условно-консервативная схема дает изменение спектрального
Рис. 1. Эволюция вдоль 2 максимальной интенсивности (а, в, д) и инкремента возмущения (б, г, е) первоначально гауссова импульса, распространяющегося в нелинейной среде для параметров И = 1, а = 10, у = 0,05, шв = 1/7 для шагов сеток (г; К) = (10~ ; 10" ) (а, б), (0,02; 10" ) (в, г), (5- 10" ; 10" ) (с?, е) при расчете по консервативной (7), (8) (пунктир) и условно-консервативной (7), (10) (сплошная кривая) схемам
инварианта на 3%. Инвариант энергии для консервативной схемы также не изменяется.
При приближении к частоте иъ = 1/7 (рис. 2,6) условие на шаг по пространству к становится более жестким: шаг необходимо уменьшать по сравнению с предыдущим случаем для достижения ограниченного решения. Заметим, что при отсутствии возмущения ограничения на шаг снижаются, а схемы (7), (8) и (7), (10) при тех же шагах сеток дают решение, аналогичное показанному на рис. 2, б,
И(м) I2
Рис. 2. Распределения интенсивности импульса в сечении г = 1 первоначально гауссова импульса, распространяющегося в нелинейной среде для параметров £)=1,а=10,7 = 1 при наличии возмущения на частоте шв = 5 (а), 1 (б), рассчитанные с шагами (т;/г) = (10_2;10_3) (сплошная линия), (Ю-2; 0,0125) (пунктир, а),
(Ю-2; 0,01) (пунктир, б)
но амплитуда осцилляций интенсивности в два раза меньше.
5. Выводы. Приведенные в работе результаты демонстрируют принципиальную роль, которую играет спектральный инвариант при компьютерном моделировании распространения фемтосекунд-ного импульса в оптическом волокне на основе уравнения КНУШ. Отсутствие консервативности разностной схемы по этому инварианту может привести к неограниченному росту сеточного решения рассматриваемой задачи. При этом уменьшение шагов сеток для некоторых условно-консервативных схем не устраняет, а, наоборот, усиливает проявление неустойчивости на частоте возмущения I/7.
Спектральный инвариант описывает изменение в среде спектральной гармоники именно на этой частоте и показывает, что она не должна возрастать в процессе распространения светового импульса. Тем самым он является неотъемлемой частью математической модели, применяемой для адекватного описания распространения фемтосекундного импульса в рамках КНУШ.
Введенное в начальное распределение возмущение может появиться в процессе расчетов непосредственно из-за ошибок округления, в результате чего при определенных шагах сеток будут иметь место описанные выше закономерности на более длинных трассах распространения.
СПИСОК ЛИТЕРАТУРЫ
1. Ахманов С.А., Выслоух В. А., Чирки н А. С. Оптика фемтосекундных лазерных импульсов. М.: Наука, 1988.
2. Агравал Г. Нелинейная волоконная оптика. М.: Мир, 1996.
3. Трофимов В. А. К вопросу об инвариантах нелинейного распространения фемтосекундных импульсов // Изв. вузов. Радиофизика. 1992. 35. № 6-7. С. 618-621.
4. Трофимов В.А. О новом подходе к моделированию нелинейного распространения сверхкоротких лазерных импульсов // ЖВМиМФ. 1998. 38. № 5. С. 835-839.
5. Варенцова С.А., Волков А.Г., Трофимов В.А. Консервативная разностная схема для задачи распространения фемтосекундного лазерного импульса в кубично нелинейной среде // ЖВМиМФ. 2003. 43. № 11. С. 1709-1721.
6. Самарский A.A. Теория разностных схем. М.: Наука, 1989.
7. Выслоух А. В., Иванова И. С., Магницкий С. А., Трофимов В. А. Модуляционная неустойчивость световых пучков и импульсов при их распространении в поглощающих средах // Оптика и спектроскопия. 2000. 88. № 3. С. 456-464.
Поступила в редакцию 24.12.03
УДК 519.6
М. С. Никольский
О ДОСТАТОЧНОСТИ ПРИНЦИПА МАКСИМУМА ПОНТРЯГИНА
В НЕКОТОРЫХ ОПТИМИЗАЦИОННЫХ ЗАДАЧАХ1
(кафедра оптимального управления факультета ВМиК)
В математической теории оптимального управления одним из важных и интересных вопросов является вопрос о достаточности принципа максимума Понтрягина (см., например, [1-12]). В настоящей статье мы рассматриваем этот вопрос для задач с фиксированным временем и с различными концевыми условиями.
1 Работа выполнена при поддержке РФФИ (проекты 02-01-00334, 03-01-00474).