Математическое моделирование процессов филаментации в средах с кубической нелинейностью

( Mathematical Modeling of Filamentation Processes in the Medium with a Cubic Nonlinearity
Preprint, Inst. Appl. Math., the Russian Academy of Science)

Балашов А.Д., Пергамент А.Х.
(A.D.Balashov, A.Kh.Pergament)

ИПМ им. М.В.Келдыша РАН

Москва, 2004

Аннотация

В работе рассмотрена проблема распространения интенсивного лазерного излучения в среде с кубической нелинейностью. Исследован вопрос о выборе математических моделей и адекватных алгоритмов, в частности, о выборе оптимальных сеток и разностных схем для математического моделирования процессов филаментации. На основе математического моделирования параксиального приближения и развитой Беспаловым и Талановым теории неустойчивости плоской волны в нелинейной среде получены результаты как для отдельных филамент (filament), так и для их взаимодействия. Результаты моделирования различных режимов филаментации показали, что начальная стадия неустойчивости с большой степенью точности описывается моделью Беспалова–Таланова. Показано, что сближение филамент на фоне плоской волны может по-разному влиять на развитие каждой из них в зависимости от их параметров.

Abstract

In the thesis the problem of intensive laser beam propagation in the medium with a cubic nonlinearity is considered. The choice of mathematical models and adequate algorithms is investigated, partially, the choice of optimal grids and difference schemes for the mathematical modeling of filamentation processes. On the basis of paraxial approximation and the Bespalov–Talanov theory of the plane wave instability the results have been received both for filaments and their interaction. The mathematical modeling of the filamentation regimes has shown that the Bespalov–Talanov theory describes the initial stage of the instability with high accuracy. It is established that the approach of filaments on the plane wave background may influence the wave evolution in a different way depending on their parameters.


Введение

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

При этом необходимо решить несколько различных по характеру задач. В начале нужно выбрать и обосновать математическую модель, описывающую конкретные режимы распространения излучения в усилительных каскадах мощных лазерных систем. Одной из широко распространенных моделей является приближение Фока-Леонтовича, более известное в научной литературе как параксиальное приближение, которое представляет собой параболическое нелинейное уравнение Шредингера (НУШ) для излучения в нелинейной среде. На основе модели НУШ можно исследовать эффекты схлопывания пучков в процессе самофокусировки [8], характер поля вблизи нелинейного фокуса [6,8], взаимодействие пучков.

В процессе первоначального исследования мелкомасштабной самофокусировки Ч.Х.Таунсом [1] было установлено наличие критической мощности для аксиально-симметричных пучков. Впоследствии В.И.Беспалов и В.И.Таланов «на кончике пера» открыли поперечную неустойчивость плоской волны в нелинейной среде [2].

Цель настоящей работы – разработка численных методов, позволяющих осуществить математическое моделирование распространения мощного лазерного излучения в нелинейных средах, в том числе явления мелкомасштабной самофокусировки или, как сейчас принято называть, филаментации, включая описание нелинейной стадии процесса развития неустойчивости. Этот алгоритм должен, в частности, адекватно воспроизводить значение критической мощности и значения инкрементов неустойчивости, что, как известно, является одной из самых сложных задач вычислительного анализа. Кроме того, на нелинейной стадии для достаточно больших значений величины, характеризующей отношение амплитуды сигнала к начальной амплитуде, критерием корректности численного метода может служить точность аппроксимации асимнтотического автомодельного решения вблизи точки коллапса [8].

В данной работе приведены результаты исследования решения НУШ и его гидродинамической аналогии, проведен анализ численных результатов неустойчивости Беспалова–Таланова как в одномерном, так и в двумерном случаях. Этот анализ дает более подробные знания о влияния шума на развитие пучка в среде с кубической нелинейностью. Именно наличие шума в пучке на входе в нелинейную среду (начальное условие для НУШ) остается пока главным объяснением развала цилиндрической симметрии пучка и появления нитей мелкомасштабной самофокусировки. Потребность в таком детальном рассмотрении возникла из-за появившихся в последнее время других объяснений этого явления, например, в работе [9], где предлагается альтернативное детерминированное объяснение, основанное на векторных эффектах.

Также исследовано влияние друг на друга при взаимном сближении двух филамент, имеющих гауссовское распределение с разными параметрами. Полученные результаты распространения филамент в нелинейной среде предполагается использовать при расчетах и анализе численных результатов распространения светового пучка с заданными случайным образом амплитудно-фазовыми возмущениями.

 

1.    Физическая постановка задачи. НУШ

 

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

в отклике можно выделить слагаемые на частоте внешнего поля

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

, где  и .

Минуя подробный вывод из уравнения Геймгольца по предложенному М.А.Леонтовичем методу разделения искомого решения на быстроменяющуюся функцию, служащую геометро-оптическим «костяком» решения, и медленно меняющуюся «функцию ослабления», приведенный в [8], запишем для поперечных компонент стационарного электромагнитного поля в изотропных средах с кубической нелинейностью уравнение:

 

,

 

где – оператор Лапласа в поперечной лучу плоскости (x,y),  – волновое число,  – линейная диэлектрическая проницаемость, которая связана с нелинейной диэлектрической проницаемостью уравнением

.

 

Перейдя к безразмерным амплитудам циркулярно поляризованных полей:

,

перепишем уравнение в виде:

 

,

где. Эта система в случае круговой ( или  равно нулю) или линейной () поляризации пучков соответствующей заменой переменных , ,  или приводится к стандартному виду НУШ:

.

 

2.    Исследование неустойчивости Беспалова–Таланова

 

Впервые теория образования нитевидной структуры светового пучка в результате самофокусировки была дана Беспаловым и Талановым в 1966 г. [2] Было показано, что в нелинейном диэлектрике амплитудно-фазовые возмущения плоской электромагнитной волны приводят к ее распаду на отдельные пучки, имеющие разные длины самофокусировки в зависимости от масштаба первоначального возмущения. Для нахождения данной зависимости исследуем решение НУШ типа плоской волны, т.е.

,

где  – амплитуда невозмущенной волны (, ).

После подстановки линеаризуем полученные уравнения, разделим действительную и мнимую части ():

;      .

 

Для определения устойчивости решения подставим уравнения возмущения типа

(произвольное возмущение представимо в виде суперпозиции таких полей). Получим следующую зависимость

.

Далее, в зависимости от значений , можно определить характер развития возмущений (рис. 2.1).

 

Рис. 2.1. Зависимость характера возмущений от значений

 

 

При численном решении задач важно соблюдать сохранение интеграла

,

чтобы обеспечить консервативность схемы.

 

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

 

Рассмотрим одномерный случай. Начальное возмущение зададим в виде

, т.е. Ψ0=20

 

Расчеты проводились для разных значений  и полностью подтвердили изложенную выше теорию. На рис. 2.2 представлено поведение решения в точке  для разных режимов. Если величина волнового вектора  достаточно велика, т.е , то решение устойчиво и периодически изменяется по x и z (1-я кривая на рис. 2.2, рис. 2.3.а). Значение  - это граница устойчивости, и решение для таких  не меняется по z (2-я кривая на рис. 2.2, рис. 2.3.б). Если , то решение экспоненциально растет в направлении z. В начальной стадии оно при этом остается гармоническим по x. Далее, с ростом z решение перестает быть гармоническим по x (рис. 2.3-г и 2.3.д). Таким образом, как и следовало ожидать, только начальные этапы развития неустойчивости соответствуют экспоненциальному росту (кривые 3, 4, 5 на рис. 2.2). Случай наискорейшего развития неустойчивости, при котором инкремент неустойчивых возмущений максимален (), соответствует кривой 4 на рис. 2.2.

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.2.

Поведение амплитуды колебаний  в точке с координатой x=0 в зависимости от z. Цифрами обозначены соответствующие графики:

1 - Рис. 2.3.a; 2 - Рис. 2.3.б; 3 - Рис. 2.3.в; 4 - Рис. 2.3.г; 5 - Рис. 2.3.д


На рис. 2.3.a-д показано распределение возмущений в начальный момент времени t=z=0 (сплошная линия) и в последующие моменты времени (пунктирная линия), которые отмечены на рис. 2.2 на соответствующих графиках.

 

Frame 001
Created with Tecplot 9.0-2-3

Frame 001
Created with Tecplot 9.0-2-3

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.3.a.

, т.е. возмущения устойчивы. Происходят периодические колебания

 

Рис. 2.3.б.

, соответствует границе между двумя состояниями, возмущения никак не должны развиваться

 

Рис. 2.3.в.

т.е. возмущения начинают развиваться,

увеличивая разброс амплитуд на фронте волны

 

Frame 001
Created with Tecplot 9.0-2-3

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.3.г. , соответствует максимальному инкременту развития возмущений.

Минимальное время развития неустойчивости t=0.016

Рис. 2.3.д. ,

т.е. возмущения развиваются не так быстро

t=0.021 > t=0.016 (ср. рис. 2.3.г)

 

Описанная выше картина развития периодических возмущений в целом сохраняется и при рассмотрении в двумерном случае начального возмущений вида:

 

, т.е. Ψ0=20, .

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.4.а, , т.е. устойчивость.

z=0 – черный, z=0.002 – серый.

Периодическое движение графика по z с периодом z=0.004

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.4.б. , т.е. граничный случай.

При увеличении z график почти не меняется

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.4.в. , т.е. неустойчивость.

z=0 – черный, z=0.01 – серый.

Далее по z нарастание неустойчивости продолжается

Из приведенных результатов видно, что неустойчивыми, а значит изменяющими сигнал, являются шумы с более низкими частотами. Уменьшить их влияние можно, увеличив их частоту, например, задавая на входе в нелинейную среду сигнал и регулярный высокочастотный шум.

 

Картина становиться более ясной, если полученные на каждом шаге по z значения Ψ изображать на комплексной плоскости (ReΨ, ImΨ). Тогда значения рассматриваемого квадрата на плоскости (x,y) в сечении z=0 перейдет в горизонтально расположенный отрезок на плоскости (ReΨ, ImΨ), длина которого равна . Далее на рисунках 2.5.а-в приведены графики, соответствующие предыдущим результатам 2.4.а-в.

 

Frame 001
Created with Tecplot 9.0-2-3

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.5.а. ,

т.е. устойчивость

Рис. 2.5.б. ,

т.е. граница устойчивости


 

В случае устойчивости  (рис. 2.5.а) этот отрезок при увеличении z, сохраняя свою длину, будет вращаться против часовой стрелки вокруг (0,0), одновременно с этим вращаясь вокруг собственного центра против часовой стрелки.

·        Повороту вокруг собственного центра на угол 2π соответствует один период колебаний величины  в плоскости (x,y). Частота колебаний определяется величиной .При  и z=0.004 значение угла равно .

·        Поворот вокруг (0,0) против часовой стрелки соответствует сделанному предположению , т.е. угол . Величина  при z=0.004 принимает значение , что соответствует периоду колебаний на рис. 2.5.а. При  (рис.3.б) поворот вокруг собственной оси отсутствует, h=0.

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.5.в. Чисто мнимые начальные условия.

, т.е. неустойчивость

 

В области неустойчивости сохраняется поворот вокруг (0,0), но при этом изменяется длина, а потом и кривизна «отрезка». Если задать чисто мнимые начальные условия, то картина на комплексной плоскости (ReΨ, ImΨ) будет отличаться только поворотом на угол  вокруг центра (0,0) от картины, получаемой при чисто вещественных начальных условиях (рис. 2.5.в). При этом .

Таким образом, в устойчивом режиме при сохранении характера решения (гармонические колебания) сдвиг по фазе между действительной и мнимой частью функции Y(x,y,z) меняется от 0 до 2p в течение одного периода колебаний по z. Колебания в плоскости (x,y) становятся ангармоническими, понятие сдвига по фазе между ReΨ и ImΨ теряет смысл, т.к. эта величина различна в разных точках плоскости (x,y).

 

2.2.        Описание численных методов. Сравнение результатов

 

Были рассмотрены различные модели исследуемого процесса и методы их решения.

1. Квазиодномерная задача (цилиндрически симметричный случай) решалась в эйлеровых и лагранжевых координатах. При введении лагранжевых координат использовалась гидродинамическая аналогия, предложенная в [6]. Далее приведен вывод соответствующих уравнений и ход решения задачи.

Комплекснозначная функция Ψ представляется в виде . Далее, переписав исходное уравнение, выделим мнимую и действительную части и в этих двух уравнениях перейдем к переменным t=z, v= (для этого от одного из уравнений необходимо взять частную производную ). Получим уравнения гидродинамики с дифференциальной связью давления p с плотностью ρ:

;

 

;

где

;     .

Введем энергетическую координату , тогда уравнения перепишутся в виде:

;     ;     ;

где

;     .

 

Краевые условия для исходного уравнения переходят в

 

v (0,t)=0, v(I1,t)=0;

 при , ,

 

а начальные данные – в

v(m,0)=v0(m), ,

где .

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

;     .

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

 

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

Text Box: Рис. 2.6. 
Нерегулярная сетка
2. Двумерная задача решалась на прямоугольной и полярной сетках, а также на нерегулярной сетке, совмещающей в себе оба вида сеток (рис. 2.6), о которой см. ниже. В качестве численного метода брался метод Ньютона, шаг по z выбирался таким, чтобы схемная вязкость была маленькой, а именно изменение интеграла мощности I1 не превосходило бы 0.1%. Граничные условия в центре возмущения – равенство нулю потока; на границе области – равенство нулю самой функции Ψ или потока через границу.

 

В качестве сравнения численных результатов на различных сетках приведем зависимости величины интеграла I1(z) и инкремента  на оси возмущения с параметрами A=1, b=4 (рис. 2.7.а-б). Цифрами на графике изображены следующие численные решения:

 

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

 

2.     Квазиодномерная задача (цилиндрический случай). Решение в эйлеровых координатах. При решении этой и последующих описанных задач ширина области бралась в  раза больше характерного размера изменения амплитуды входного луча. Отрезок [0,18] изменения координаты r равномерно разбивался на 400 ячеек. Шаг по времени уменьшался обратно пропорционально квадрату модуля интенсивности на оси луча.

 

3.     Двумерная задача. Сетка полярных координат (r,φ) разбивалась на 16 секторов по φ и на 400 параллелей по r. Шаг по времени уменьшался обратно пропорционально возрастанию модуля интенсивности на оси луча. В случае, если шаг по времени будет аналогичен случаю 2, кривые 2 и 3 на рис.2.7.б будут совпадать.

 

4.     Двумерная задача. Ортогональная сетка (rx, ry) разбивалась на  квадратов. Шаг по времени уменьшался обратно пропорционально возрастанию модуля интенсивности на оси луча.

 

Рис. 2.7.а. Сравнение инкремента

Рис. 2.7.б. Сравнение мощности

 

Особое внимание при исследовании НУШ уделяется описанию поведения решения в окрестности коллапса. Поэтому необходимо разработать метод, эффективный на нелинейной стадии развития неустойчивости. Выбор ортогональной сетки с описанными включениями нерегулярной сетки при решении задачи продиктован целью расчета развития возмущения при задании начального шума. Каждое «шумовое» возмущение будет задаваться отдельно. В последующем исследование эволюции случайных возмущений будет осуществляться с помощью параллельных алгоритмов, основанных на использовании нерегулярной сетке. Необходимость решения такой численной задачи указана в [5]: «При постановке численных расчетов возникает следующая трудность. Поперечные размеры пучка значительно превышают размеры возникающих в сечении пучка неоднородностей. Трехмерная задача распространения света с учетом дифракции, нелинейных эффектов и нестационарности усиления является весьма громоздкой».

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

, ,

с применением интегрального тождества

.

Тогда для произвольной ячейки (обозначим ее индексом 0, как центральную среди рассматриваемых соседних ячеек, обозначаемых k(0)):

,

где  – поток из рассматриваемой центральной ячейки в соседнюю k-ую ячейку,  – среднее значение  на границе 0k,  – вектор нормали к этой границе,  – вектор, направленный по нормали к границе рассматриваемых ячеек и по величине равный длине этой границы.

 

Text Box: Рис. 2.8. Нерегулярная сетка. Подробное описаниеРассмотрим в качестве примера верхнюю шестиугольную ячейку рассматриваемого разбиения на нерегулярной сетке. Для нее из вышеприведенных формул следует:

 

где

 – объем ячейки,

, , , .


2.3.        Применение анализа неустойчивости Беспалова–Таланова к одиночному возмущению

 

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

.

В результате численно подтвердился тот факт, что при возрастании начальной ширины возмущения до определенного момента инкремент неустойчивости растет. При этом решение до определенного момента времени (t=z) представляет собой уединенный максимум: процессы последующего распада максимума не рассматриваются. Первоначальная полуширина пучка есть b. Если b мало, то это означает, что пучок обладает широким спектром и если доля низкочастотных составляющих (которые, согласно теории Беспалова–Таланова, экспоненциально возрастают), мала, то такой пучок на линейной стадии неустойчивости не будет возрастать. Для таких начальных условий сначала произойдет уширение пучка, сопровождающееся падением амплитуды. Действительно (рис. 2.9.б), полуширина пучков с b=0.05, 0.08, 0.11 первоначально растет. Одновременно для тех же b падает амплитуда поля (рис. 2.9.а).

 

Frame 001
Created with Tecplot 9.0-2-3

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.9.а. Изменение по z амплитуды поля  на оси пучка   при Ψ0=10, А=1 и различных     значениях b

Рис. 2.9.б. Изменение по z полуширины пучка l(z)

 при Ψ0=10,

А=1 и различных значениях b


Для исследования уединенного максимума напишем преобразование Фурье для функции :

 

, где .

 

Т.е.

есть разложение первоначального возмущения в суперпозицию рассмотренных волн. Чем меньше значение b, тем шире спектр сигнала и тем меньшая доля его приходится на область малых k (здесь ), т.е. для малых b сигнал будет расти медленнее, чем для больших. Но в целом такое начальное возмущение неустойчиво, т.к. всегда существуют k, удовлетворяющие условию Беспалова–Таланова.

Если начальная ширина пучка достаточно велика, то в силу того, что доля низкочастотных растущих компонент велика, амплитуда такого пучка возрастает, а его полуширина падает. Для гауссовых пучков спектральное распределение также гауссово, при этом большая часть спектра сосредоточена в области . Для того, чтобы эта часть спектра удовлетворяла условиям Беспалова–Таланова, требуется , т.е. . В данном случае b>0.14, т.к. . Расчеты показывают, что вплоть до амплитуда сначала падает. Одновременно с падением амплитуды растет полуширина пучка. В дальнейшем на нелинейной стадии такие пучки будут себя вести так же, как более широкие, но с большим запозданием по z.

По мере развития неустойчивости подавляются те части спектра, где  и происходит уширение пучка, т.е. сужение спектра к малым k (рис. 2.9.а-б для ). Т.е. высокочастотные составляющие при уширении пучка, удовлетворяющие условию устойчивости, затухают, тем самым увеличивая преобладание низкочастотных (неустойчивых) составляющих (это происходит из-за сохранения первого интеграла – суммарной мощности пучка).

Однако значение имеет соотношение спектра не только с границей неустойчивости, но и с областью максимального инкремента. При дальнейшем увеличении b спектр продолжит сужаться, тем самым все большая часть k будет меньше , что соответствует более медленному развитию неустойчивости (рис. 2.9.а-б для b=0.2 и b=0.4). В этом случае устойчивые составляющие отсутствуют, вследствие чего уширения пучка не происходит. Данный анализ справедлив лишь на начальной, «линейной» стадии развития неустойчивости.

 

Аналогичным образом можно описать эволюцию пучка при отсутствии фона. При этом аналогом собственного фона для гауссовского возмущения будет являться величина

.

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

 

Frame 001
Created with Tecplot 9.0-2-3

Frame 001
Created with Tecplot 9.0-2-3

Рис. 2.10.а.

Изменение амплитуды  поля

на оси ограниченного пучка (t=z)

при А=1 и различных значениях b.

Рис. 2.10.б.

Изменение полуширины l(t) (t=z) пучка

при А=1 и различных значениях b.

 


Полученный результат согласуется также с общей теорией, развитой Таунсом [1], где утверждается, что необходимым условием для неустойчивости возмущения является превышение мощностью пучка определенного значения, так называемой критической мощности самофокусировки , которая достигается в рассмотренном примере при b=2.72.

 

3.    Два возмущения. Сближение и взаимное влияние.

 

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

Рассмотрим два возмущения на ненулевом фоне Y0=10. L – расстояние между центрами возмущений (L=1, L=0.6).

 

Случай 1 – разные амплитуды (А1=1, А2=1.3.) при одинаковой ширине (b=0.3). Рассматриваются значения z=0 и z=0.57.

1.1. L=1 (рис.3.1-б)

Взаимное влияние почти не наблюдается. При рассмотрении уединенного максимума установлено наличие кольцевых структур (рис. 3.2), которые представляют собой волны на фоне, затухающие по мере удаления от возмущения. На представленном рисунке уровень фона равен 10. Кольцевая структура здесь представлена углублением (т.е. фон развития возмущения понижен) и небольшим возвышением над уровнем 10 примерно на 0.03. Для двух максимумов также наблюдаются подобные структуры (рис. 3.3). В рассмотренном случае заметно отличие только в области, расположенной ниже уровня Ψ0=10, углубление заметно больше в области между возмущениями, чем на соответствующем одном возмущении с параметрами A=1, b=0.3 (рис. 3.1-а) или A=1.3, b=0.3 (рис. 3.1-г).

1.2. L=0.6 (рис. 3.1-в)

Взаимное влияние существенно. Амплитуды обоих возмущений существенно уменьшились по сравнению с амплитудой таких же возмущений на том же срезе, рассмотренных независимо (рис. 3.1-а,г). По сравнению с предыдущим (L=1) случаем углубление кольцевой структуры несколько изменилось. Теперь вместо одного минимума, как в предыдущем случае, их стало два, они существенно глубже (рис. 3.4), но между самими возмущениями величина Y не упала. Именно это явилось причиной более медленного развития амплитуд обоих взаимодействующих возмущений.

 

 

Frame 001
Created with Tecplot 9.0-2-3

    а         |                б              |              в            |            г         

  

Рис. 3.1. Взаимное влияние двух возмущений

с одинаковой шириной b=0.3 друг на друга.

а – одно возмущение A=1;

б – два возмущения A=1 и A=1.3 на расстоянии L=1;

в – два возмущения A=1 и A=1.3 на расстоянии L=0.6;

г – одно возмущение A=1.3

 

 

   Frame 001
Created with Tecplot 9.0-2-3

Рис. 3.2. Вид сверху. Возмущение A=1, b=0.3. Соответствует рис. 3.1-а

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 3.3. Вид сверху.

Два возмущения

A=1, b=0.3 и A=1.3, b=0.3

на расстоянии L=1. Соответствует рис. 3.1-б

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 3.4. Вид сверху.

Два возмущения

A=1, b=0.3 и A=1.3, b=0.3

на расстоянии L=1. Соответствует рис. 3.1-в

 

Случай 2 – разные ширины (b1=0.3, b2=0.5) при одинаковой амплитуде (A=1). Рассматриваются значения z=0 и z=0.54.

2.1. L=1 (рис. 3.5-б)

Взаимное влияние почти не наблюдается. Заметно отличие только в углублении кольцевой структуры в области между возмущениями - она немного глубже, чем на соответствующем одном возмущении с параметрами A=1, b=0.3 (рис. 3.5-а, 3.6). Приводится данное сравнение, т.к. у одного возмущения с параметрами A=1, b=0.5 (рис. 3.5-г) кольцевая структура на данном срезе z ещё незаметна.

 

Frame 001
Created with Tecplot 9.0-2-3

    а         |                б              |              в            |            г         

 

Рис. 3.5. Взаимное влияние двух возмущений

с одинаковой амплитудой A=1 друг на друга.

а – одно возмущение b=0.3;

б – два возмущения b=0.3 и b=0.5 на расстоянии L=1;

в – два возмущения b=0.3 и b=0.5 на расстоянии L=0.6;

г – одно возмущение b=0.5.

 

2.2. L=0.6 (рис. 3.5-в)

Взаимное влияние существенно. Амплитуда более тонкого (b=0.3) возмущения существенно увеличилась по сравнению с амплитудой одного такого возмущения на том же срезе. Этот факт можно объяснить тем, что при таком сближении возмущений поднимается «фон» для более тонкого (b=0.3) за счет возмущения с большей шириной – это можно увидеть при сравнении начальных условий для обоих рассмотренных случаев на рис. 3.5-б,в. Опять, при сближении углубление кольцевой структуры между возмущениями стало ещё более существенным по сравнению с предыдущим (L=1) случаем (рис. 3.5-в, 3.7).

 

 

   Frame 001
Created with Tecplot 9.0-2-3

Рис. 3.6. Вид сверху.

Два возмущения

A=1, b=0.3 и A=1, b=0.5

на расстоянии L=1.

Соответствует рис. 3.5-б

 

Frame 001
Created with Tecplot 9.0-2-3

 

Рис. 3.7. Вид сверху.

Два возмущения

A=1, b=0.3 и A=1, b=0.5

на расстоянии L=0.6.

Соответствует рис. 3.5-в

 

 

Ни в одном из описанных случаев не было выявлено существенных изменений в смещении центров возмущений. В поисках возможных изменений центров возмущений и интуитивно понятных взаимопоглощений было также рассмотрено влияние возмущений друг на друга без наличия фона.


Случай 3 – аналогично рассмотренному выше случаю, возьмем сначала возмущения с одинаковой шириной (b=1) и разными амплитудами (A1=3, A2=4).

3.1. L=1 (рис. 3.8-б)

Взаимное влияние почти не наблюдается.

3.1. L=0.6 (рис. 3.8-в)

Взаимное влияние существенно. Амплитуды обоих возмущений существенно увеличились по сравнению с амплитудой таких же возмущений на том же срезе, рассмотренных отдельно. Ввиду взаимного увеличения фона при наложении возмущений друг на друга наблюдается более быстрый рост Ψ.

 

Frame 001
Created with Tecplot 9.0-2-3

    а         |                б              |              в            |            г         

 

Рис. 3.8. Взаимное влияние двух возмущений

с одинаковой шириной b=1 друг на друга.

а – одно возмущение A=3;

б – два возмущения A=3 и A=4 на расстоянии L=2;

в – два возмущения A=3 и A=4 на расстоянии L=1.7;

г – одно возмущение A=4

 

 

 

Случай 4 – при рассмотрении возмущений с одинаковой амплитудой (A=3) и разными ширинами (b=1, b=0.8) имеет место значительное сближение центров и поглощение более широким возмущением (b=1) более узкого (b=0.8). Равномерные срезы по z центрального профиля таких возмущений представлены на рис. 3.9. При отсутствии фона не наблюдается образование кольцевых структур.

 

 

Frame 001
Created with Tecplot 9.0-2-3

Рис. 3.9.

Развитие двух возмущений

с параметрами A=3 b=1 и A=3 b=0.8.

Графики представляют собой границу

продольного сечения в разные моменты времени,

распределенные равномерно



Заключение

Целью настоящей работы была разработка эффективных алгоритмов для расчета явлений филаментации, способных позволить в будущем произвести математическое моделирование процесса филаментации в активных элементах реальных лазерных систем. Основные требования к алгоритму – это правильная передача инкремента и обеспечение консервативности, т.е. сохранение энергии. Конструирование алгоритма с такими свойствами для процессов, в которых наблюдается неустойчивость, является одной из самых сложных задач вычислительного анализа. При использовании гидродинамической аналогии неявная консервативная схема на достаточно подробной сетке в лагранжевых переменных обеспечивает строгое сохранение энергии. При использовании эйлеровых координат в одномерных моделях, а также в том случае, если при расчете используются комбинированные сетки, необходимо использовать схемы с весами. В результате точность сохранения полной энергии составляет доли процента. При использовании на прямоугольной сетке консервативных схем с весами закон сохранения полной энергии соблюдается также с высокой точностью. Однако оценка величины инкремента (рис. 2.7.a) показывает, что для схемы на прямоугольной сетке эта величина несколько меньше, чем для других методов расчета. Поэтому в тех случаях, когда необходимо исследование развитой стадии неустойчивости в двумерных задачах целесообразно использовать разностные схемы на комбинированных сетках.

Результаты моделирования различных режимов филаментации показали, что начальная стадия неустойчивости с большой степенью точности описывается моделью Беспалова–Таланова. Что касается развитой стадии, то получены режимы с образованием колец (рис. 3.2). Рассмотрены также процессы взаимодействия филамент, описывающие нелинейную стадию процессов неустойчивости. Показано, что сближение филамент на фоне плоской волны может по-разному влиять на развитие каждой из них в зависимости от их параметров. Например, филаменты с разными амплитудами при одинаковой ширине будут расти медленнее при сближении. В то же время сближение пучков на нулевом фоне способствует их более быстрому росту, что может быть объяснено теорией Беспалова–Таланова, если рассматривать взаимодействующие пучки, как новое начальное состояние.

 


Литература

 

[1]. R.Y.Chiao, E.Garmire, С.H.Townes «Self-trapping of optical beams» Phys. Rev. Lett., V.13, №15, p.479, 1964

 

[2]. В.И.Беспалов, В.И.Таланов «О нитевидной структуре пучков света в нелинейных жидкостях» Письма в ЖЭТФ. Т.3, с.471-476, 1966

 

[3]. Н.Б.Баранова, Н.Е.Быковский, Б.Я.Зельдович, Ю.В.Сенатский «Дифракция и самофокусировка излучения в усилителе мощных световых импульсах» Квантовая электроника, Т.1, №11, с.2435-2449, 1974

 

[4]. Л.М.Дегтярев, В.В.Крылов «Метод численного решения задач динамики волновых полей с особенностями» Журнал вычислительной математики и математической физики, Т.17, ноябрь-декабрь 1977

 

[5]. Л.А.Большев, Л.М.Дегтярев, А.М.Дыхне и др. «Численное исследование мелкомасштабной самофокусировки импульсов света в усилителях на неодимовом стекле» Препринт ИПМ РАН им. М.В.Келдыша №109, 1979

 

[6]. Л.М.Дегтярев, В.В.Крылов «Гидродинамическое описание самофокусировки пучков света в кубической среде» Сборник научных трудов под редакцией А.А.Самарского «Изучение гидродинамической неустойчивости численными методами», Москва, 1980

 

[7]. А.Н.Тихонов, В.Я.Арсенин, В.И.Павлов, А.Х.Пергамент «Исследование расходимости излучения в мощных лазерных усилителях на активных элементах прямоугольно сечения» Препринт ИПМ РАН им. М.В.Келдыша №41, Москва, 1981

 

[8]. С.Н.Власов, В.И.Таланов «Самофокусировка волн» ИПФ РАН, Нижний Новгород, 1997

 

[9]. G.Fibich, B.Ilan «Vectorial and random effects in self-focusing and in multiple filamentation» Phisica D, V.157, p.113-147, 2001