Диссертация (Формирование и характеристики плазменных каналов при филаментации фемтосекундного лазерного излучения в воздухе), страница 12

Правда, после первого шага первый полушаг очередной итерации можно объединить с последним полушагом предыдущей, но тогда не удастся узнать значения физическихвеличин после полного шага. Поэтому из приведенных схем использовалась последняя схема (2.79). Она лишь немного уступает в скорости расчета простым схемам (2.76, 2.77), нопозволяет использовать больший шаг интегрирования ∆.В случае стационарной задачи самофокусировки схема метода расщепления остаетсяпрежней, за исключением состава линейных и нелинейных членов (учитываются только дифракция и керровская нелинейность).2.5.3. Учет дифракции и дисперсииНеосесимметричная задачаУчет линейных факторов в неосесимметричной постановке задачи филаментации илистационарной задачи самофокусировки осуществлялся в спектрально-угловом пространстве.— 43 —Уравнение (2.55a) на этапе учета линейных факторов имело вид(, , , )120= ∆⊥ (, , , ) +2∫︁∞(︀)︀Ω (, , , Ω) 2 (0 + Ω) − (0 + 1 Ω)2 Ω ,−∞(2.80)Перейдем от поля (, , , ) и частотного спектра (, , , Ω) (2.9) к его частотноугловому спектру1(, , , ) =(2)3∫︁ ∫︁ ∫︁( , , , Ω)Ω + + Ω.(2.81)Уравнение (2.80) примет вид)︀(︀( , , , Ω)= (−2 − 2 )( , , , Ω) + 2 (0 + Ω) − (0 + 1 Ω)2 ( , , , Ω).(2.82)Нетрудно видеть, что аргументы , и Ω являются параметрами, а само уравнениепревратилось в обыкновенное линейное дифференциальное уравнение.

Его решение при известном начальном условии ( , , , Ω) дается формулой:20−2 − 2 + 2 (0 + Ω) − (0 + 1 Ω)2∆20( , , + ∆ , Ω) = ( , , , Ω) · .(2.83)Решение в пространственно-временном представлении получается применением обратного преобразования Фурье:∫︁ ∫︁ ∫︁(, , + ∆ , ) =( , , + ∆ , Ω) · −Ω − − Ω .(2.84)В случае стационарной задачи, в которой не учитывается дисперсия, в 2.82 следует отбросить дисперсионный член, считая, что (0 + Ω) = 0 + 1 Ω.После дискретизации переменных преобразования (2.81, 2.84) заменяются на дискретныепреобразования Фурье.Осесимметричная задачаУчет дисперсии в осесимметричной задаче происходил аналогично неосесимметричной,но без перехода к угловым гармоникам: 2 (0 + Ω) − (0 + 1 Ω)2∆20.(, + ∆ , Ω) = (, , Ω) · (2.85)Учет дифракции выполнялся отдельно от учета дисперсии, то есть производилось дальнейшее расщепление по физическим факторам.Для радиальной части оператора Лапласа была построена аппроксимация третьего по-— 44 —рядка точности1 (︂)︂⃒ ⃒⃒≈ −1 + + +1 ⃒(2.86)с коэффициентами+1, =(−1 − )(−1 − +1 )−13− =,(+1 − )(+1 − −1 )3− = − − .(2.87)(2.88)(2.89)Используя приведенную аппроксимацию, построим схему Кранка-Николсона для уравнения дифракции(︂)︂(, , )1 20=(2.90) 20)︀1 (︀+1− +1= −1 + + +1 + +1+ +1−1 + +1 ,∆2(2.91)откуда легко получить систему линейных алгебраических уравнений с трехдиагональнойматрицей: +1−1 + ( −40 +140 ) + +1) − +1 ,+1 = − −1 − ( +∆∆ (2.92)где = 1, − 2.

Остается добавить граничные условия.На правом конце использовались нулевые граничные условия второго рода, аппроксимация которых имела вид+1+1(2.93) −2 − −1 = 0.Нетрудно видеть, что при смещении нулевого расчетного узла на полшага 1 = 30 , поэто21.му 0 = 0. Кроме того, можно формально положить −1 = − ∆. Тогда 0 = −0 =2(∆)2В итоге, при = 0 разностное уравнение (2.92) будет иметь вид(︂240−−2(∆)∆)︂+10(︂)︂22402+1+1 = − −+0 − ,22(∆)(∆)∆(∆)2 1(2.94)или после упрощения(︂(︂)︂)︂20 (∆)220 (∆)2+1+10 − 1 = − 1 −0 + 1 .1+∆∆(2.95)Полученная система линейных алгебраических уравнений с трехдиагональной матрицейлегко решается методом прогонки.— 45 —2.5.4.

Учет нелинейных факторов и затуханияНа этом этапе необходимо проинтегрировать систему (2.55) с отброшенными дифракционным и дисперсионным членами:202 2(, , , )= 0 (∆ (, , , ) + ∆ (, , , )) (⃗, ) − 0 (, , , )(, , , ),0(2.96a) (, , , )= ((, , , ))(0 − (, , , )).(2.96b)Как показывают оценки, потери энергии импульсом за один шаг интегрирования крайнемалы, поэтому можно считать, что интенсивность (, , , ) остается «замороженной» длянелинейных процессов. Это позволяет нам проинтегрировать второе уравнение и получитьнаведенную временным слоем в данной точке пространства (, , ) концентрацию электронов (, , , + ∆ ): (, , , + ∆ ) = 0 − −((,,, ))Δ (0 − (, , , )).(2.97)Зная изменение числа электронов, можно определить потери излучения на ионизациюсреды, вычислив коэффициент (, , , ) по формуле (2.45).Керровская добавка к показателю преломления определяется по формуле (2.21).

Длячисленного интегрирования инерционного отклика использовался метод трапеций.Теперь, когда коэффициенты ∆ , ∆ , определены, можно проинтегрировать уравнение (2.96a). Оно допускает точное решение (конечно, в предположении «замороженности»интенсивности):Δ,(, , + ∆ , ) = (, , , ) · (,, , ) · 2−где (, , , ) = −0(∆ + ∆ )∆ .0(2.98)(2.99)нелинейный набег фазы.Подход к интегрирование системы (2.96) в случае осесимметричной задачи остается темже, поскольку во всех выкладках текущего раздела пространственные координаты и являются параметрами и могут быть без последствий заменены единой координатой .2.5.5. Учет турбулентных флуктуаций показателя преломленияАналитическое исследование филаментации лазерного импульса в турбулентной атмосфере затруднено ввиду случайности флуктуаций показателя преломления.

При сравнительнонебольшой мощности излучения можно считать нелинейные факторы малыми добавками ипровести рассмотрение задачи при определенных приближениях, например, в приближении— 46 —заданной формы пучка [173]. В общем случае, однако, такое приближение не оправдывается,и необходимо проводить статистическое численное моделирование.Учет влияния турбулентных флуктуаций показателя преломления атмосферы осуществлялся в рамках модели фазовых экранов [89].

Согласно этой модели, слой атмосферы толщиной (обычно 2–5 метров) заменяется тонким фазовым экраном, расположеннымпосередине слоя, при прохождении которого пучок получает фазовый набег, как от всегослоя:(, , screen + 0, ) = (, , screen − 0, ) · (,) .(2.100)Толщина стягиваемого слоя не должна быть очень большой, в противном случае распространение пучка сквозь слой уже нельзя будет рассматривать как получение фазового набега, поскольку вследствие дифракции исказится поперечное распределение интенсивности.Если расстояние между экранами превосходит максимальный по размеру масштаб флуктуаций (так называемый внешний масштаб турбулентности), то фазовые экраны можно считатьстатистически независимыми.В качестве модели турбулентной атмосферы использовалась модель с модифицированным кармановским спектром пространственных флуктуаций показателя преломления () =0.0332 ( 2+20 )−11/62exp − 2(︂)︂.(2.101)5.922, =, где 0 и 0 — внешний и внутрен00ний масштабы турбулентности.

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

Структурная постоянная 2 определяет интенсивность турбулентных флуктуаций, чем больше ее значение,тем больше амплитуда флуктуаций. Для спокойной атмосферы 2 ≃ 10−17 см−2/3 , величина2 ≃ 10−13 см−2/3 соответствует очень сильной турбулентности [174].Здесь — поперечное волновое число, 0 =2.5.6. Адаптивный шаг интегрированияПравильный выбор шага интегрирования ∆ чрезвычайно важен для получения корректного результата моделирования за приемлемое время.

По отдельности уравнения, содержащие только линейные (2.80) и только нелинейные факторы (2.96), могут быть проинтегрированы с достаточно большим шагом, поскольку известны их точные решения. Однаковзаимное влияние этих факторов не дает возможности так поступить. Нелинейный набегфазы оказывает влияние на дифракцию, приводя, в частности, к самофокусировке, и надисперсию, приводя к фазовой самомодуляции, самосжатию импульса и прочему. Изменение интенсивности при дифракции и дисперсии приводит к нарушению предположения о«замороженности» интенсивности для нелинейных процессов.

Поэтому шаг интегрированиянеобходимо ограничить.— 47 —С другой стороны, шаг интегрирования для области развитого филамента в атмосфере доходил до величин порядка 1 мкм, что гораздо меньше необходимого шага для области, предшествующей филаментации. На этапе самофокусировки, когда нелинейные факторы еще слабы, допустим гораздо больший шаг, позволяющий существенно сократить времявычислений, по крайней мере, для начального этапа распространения импульса. Однако вобщем случае невозможно предсказать, где начнется область интегрирования, требующаямалого шага ∆ . Поэтому шаг интегрирования должен быть не только переменным, но иадаптивным, подстраиваясь под получаемое решение.Для оценки пригодности текущего шага интегрирования использовался простой критерий, основанный на требовании малости нелинейного фазового набега (2.99): max | | 6 0.01.,,Константа 0.01 была подобрана эмпирически.

Если нелинейный набег превосходил пороговоезначение, следующий шаг интегрирования пропорционально уменьшался:∆+1 = ∆max | (, )|,,0.01.(2.102)2.6. Параллельные алгоритмы решения задачи филаментации2.6.1. Оценка объема вычислительной нагрузкиВ случае общей задачи филаментации (2.55) на каждом шаге по эволюционной координате необходимо иметь данные о значениях комплексной амплитуды поля во всех точкахрасчетной сетки ,, = ( , , ), 0 6 6 − 1, 0 6 6 − 1, 0 6 6 − 1. Таким образом, в памяти программы должны храниться и обрабатываться данные о значенияхполя в × × точках сетки.Количество точек по поперечным пространственным координатам определяется, с однойстороны, необходимостью разместить на сетке весь лазерный пучок, так чтобы поле на границе сетки было исчезающе малым, и, с другой стороны, необходимостью иметь высокоепространственное разрешение, чтобы прописать тонкую структуру филаментов и плазменных каналов.