Предложенный метод решения уравнения Гельмгольца с однородными граничными условиями проверен для эллиптического сечения. Метод обобщен для уравнения Гельмгольца с неоднородными граничными условиями. Обобщение проверено для круглого и эллиптического сечения.
The laser gain under nonhomogeneous boundary conditions.pdf Введение Большинство газоразрядных лазеров имеют активные элементы в виде цилиндрических трубок. Граничные условия на стенках трубок принято считать нулевыми. Мощность излучения лазера определяется взаимодействием активного вещества с полем, а поскольку поле резонатора (каустика) в измерительных лазерах имеет диаметр в два и более раз меньший, чем диаметр трубки, то для усиления в активной среде можно считать граничные условия, задаваемые на каустической поверхности, отличными от нулевых. Дифракционными потерями при этом можно пренебречь. В данной работе мы рассматриваем только основной тип колебаний в резонаторе, что справедливо для измерительных лазеров (в отличие от случаев, когда нужна максимальная мощность лазера и работают в многомодовом режиме). Пока мы также пренебрегаем пространственной неоднородностью активной среды. Изменение геометрии сечения активного элемента (прямоугольник, эллипс и т.д.) изменяет коэффициент усиления излучения. Могут появиться неоднородные граничные условия. Возникает вопрос, может ли коэффициент усиления быть больше. Судя по результатам работы [1], может. Сегодня мы располагаем более совершенными методами и средствами вычислений. Мы предлагаем метод нахождения среднего (по поперечному сечению граничной поверхности) значения коэффициента усиления измерительного лазера при неоднородных граничных условиях. Вначале мы еще раз проверим наш метод на сечении, для которого можно записать аналитическое выражение для коэффициента усиления, с однородными граничными условиями и потом перейдем к неоднородным. 1. Случай эллиптического сечения трубки с однородным граничным условием В работе [2] мы описали метод нахождения коэффициента усиления лазера при однородных граничных условиях, предложив метод решения однородного уравнения Гельмгольца с однородными граничными условиями. Были проведены расчеты для случая круглого, эллиптического и прямоугольного сечений, причем для случаев сечений в виде круга и прямоугольника было показано, что результаты совпадают с аналитическим решением. Однако и для случая эллиптического сечения также можно проверить наш метод, получив аналитическое выражение для коэффициента усиления. При постоянном вдоль оси трубки сечении геометрическая часть среднего значения по поперечному сечению граничной поверхности коэффициента усиления лазера имеет вид [1] , (1) где функция f, описывающая пространственное распределение коэффициента усиления среды, удовлетворяет уравнению Гельмгольца в области V: (2) с граничным условием . (3) Здесь Г - граница области, в которой идет поиск решения; S0 - площадь поперечного сечения трубки; k0 - усиление на оси системы. Рассмотрим уравнение (2) в эллиптических координатах (u, ν, z), которые связаны с декартовыми следующими соотношениями: , , , , , . Воспользовавшись методом разделения переменных , получаем набор трех уравнений: , , , где c и d - постоянные разделения. Обозначая , получаем из третьего уравнения каноническую форму уравнения Матье (здесь и далее при вводе функций Матье и связанных с ними величин мы придерживаемся обозначений из [3, 4], у других авторов они могут отличаться): , (4) а из второго уравнения - модифицированное уравнение Матье: . (5) Выберем за начало отсчета (u = 0) центр эллиптического сечения. Хотя диаметр каустики («пятно» генерации) зависит от z, но эта зависимость слабая (диаметры пятен на зеркалах измерительных лазеров обычно отличаются в 1.4 раза). И поскольку в данной работе мы ограничиваемся только основным типом колебаний, будем в дальнейшем считать, что зависимости от z у нас нет. Тогда постоянная разделения c = 0 (отсюда ). Как известно [4], любую функцию от (u, ν), непрерывную и однозначную внутри эллипса и обращающуюся в нуль на самом эллипсе, можно в каждой точке, лежащей внутри эллипса, разложить в ряд по решениям уравнений (4) и (5) - в ряд по произведениям функций Матье и . В нашем же случае вследствие уравнения (2) достаточно взять одно слагаемое этого ряда. Поскольку по симметрии нас интересует четное по ν периодическое решение и, как и в работе [2], мы рассматриваем вклад только от минимального значения параметра λ, то в качестве функции f следует взять , где α - нормировочный множитель для соответствия работе [2] (там мы нормировали функцию f на единицу в центре сечения). Пусть эллипс имеет полуоси a и b (b < a) и его уравнение в декартовых координатах имеет вид . Тогда для его уравнения в эллиптических координатах выполняется: , , откуда (т.е. ρ - фокальное расстояние), . В качестве α следует взять . Поскольку , , то для среднего значения коэффициента усиления из (1) получаем формулу (здесь и далее длина трубки равна единице): , (6) где . Из граничного условия (3) следует, что λ должны быть такие (конкретнее, в данном случае мы ищем минимальное такое λ), что . (7) К сожалению, интеграл (6) удается посчитать лишь численно. При вычислении значений функций Матье можно использовать ряд ([3, 4]) . (8) Формула для получается из (8) заменой , т.е. заменой cos на ch. Однако простых аналитических выражений для коэффициентов не существует - для нахождения коэффициентов можно воспользоваться рекуррентными соотношениями из [4, 5], в которые, в свою очередь, также входят параметр q и зависящая от него постоянная разделения d, а также условиями нормировки. Методы вычислений функций Матье до сих пор являются предметом исследований (см., например, [6, 7]). После того, как мы выразили нужное количество коэффициентов и составили функции Матье в виде суммы соответствующего числа слагаемых, нужно численно решить уравнение (7) и найти нужное λ и затем с данным λ вычислить интеграл (6). Это громоздкая вычислительная задача. Тем не менее вычисления можно весьма заметно упростить, если рассматривать эллипс с малым эксцентриситетом - когда b близко к a и соответственно параметр q мал. В этом случае можно воспользоваться степенными рядами для функций Матье при малых q [3]: . (9) По методу, изложенному в работе [2], мы посчитали собственное значение λ и коэффициент усиления для эллипсов с малым эксцентриситетом. Далее мы ограничились всего тремя слагаемыми в (9) для и и вычислили значение коэффициента усиления по формуле (6), численно посчитав интеграл. Результаты сравнения коэффициентов усиления двумя разными методами приведены в таблице (значение полуоси b бралось равным 1). Расчетные значения коэффициента усиления в зависимости от отношения полуосей эллиптического сечения, посчитанные двумя способами a/b λ q k/k0 по методу работы [2] k/k0 по формуле (6) 1.001 2.404 0.0029 0.432 0.430 1.01 2.393 0.029 0.432 0.430 1.02 2.381 0.057 0.432 0.430 1.03 2.370 0.086 0.432 0.431 1.05 2.348 0.141 0.432 0.431 1.07 2.327 0.196 0.432 0.433 1.10 2.298 0.277 0.431 0.435 1.15 2.253 0.409 0.431 0.441 Из таблицы видно хорошее совпадение значений (которое, что было и ожидаемо, ухудшается с ростом q, т.е. когда замена функций Матье на три слагаемых из (9) становится все более грубой), что может служить очередным подтверждением правильности предложенного в работе [2] метода решения уравнений (2), (3). Можно проверить наш метод и для случая эллипса, когда эксцентриситет не является малым. Чтобы уменьшить вычислительные трудности, можно воспользоваться таблицами с уже известными значениями коэффициентов в формуле (8) для каких-то q. В [3, 5] есть значения, например, для q = 5: A00 = 0.540612446; A02 = -0.627115414; A04 = 0.147927090; A06 = -0.017848061; A08 = 0.001282863; A010 = -0.000060723; A012 = 0.000002028; A014 = -0.000000050; A016 = 0.000000001. По формуле (8) мы составили соответствующие суммы в качестве и . С помощью метода из работы [2] нетрудно найти параметры эллипса, дающие нужное значение q = 5: b = 1, a = 2.6917 и соответственно λ = 1.7896. При данных значениях формула (7) удовлетворяется с точностью до 3∙10-6. Мы посчитали коэффициенты усиления по методу работы [2] и по формуле (6) - результат в обоих случаях равен k/k0 = 0.384, что является очередным доказательством правильности нашего метода. 2. Метод решения уравнения Гельмгольца с неоднородными граничными условиями Поскольку в реальных физических задачах возникают неоднородные граничные условия вида , (10) где u - некая известная функция, задающая граничные условия, то нужно обобщить наш метод на решение уравнений (2), (10). Поскольку по сравнению с [2] в нем модифицировались некоторые шаги и формулы, кратко приведем его. Итак, мы ищем решение (2), (10) в виде разложений по точным решениям соответствующего дифференциального уравнения (2) (подобный метод решения граничной задачи в технической литературе часто называют методом Трефтца [8]), т.е. ищем решение в виде , (11) где (r, , z) - цилиндрические координаты (от z зависимости нет), Jk - функция Бесселя неотрицательного целого порядка k (при k = 0 будет одна функция a0J0(λr)), или обозначив для единообразия через ξk соответствующие функции из (11) (т.е. , или ) в виде . Функция (11) точно удовлетворяет уравнению (2), но граничному условию (11) - уже приближенно. Под приближенным выполнением граничного условия мы понимаем следующее - выберем на границе Г некоторые N точек: ξ1, ξ2, ξ3, …, ξN и потребуем минимальности суммы квадратов разностей значений искомой функции и граничной функции в соответствующих точках: . (12) Таким образом, разность на границе Г должна быть максимально возможно близкой к нулю, т.е. на границе Г должна выполняться аппроксимация функции u комбинацией линейно-независимых функций ζk: . (13) Эту аппроксимацию делаем через ортонормированные (на границе Г) функции, поэтому перейдем от набора линейно-независимых функций {ζ1, ζ1, …, ζn} к набору ортонормированных на границе Г функций {Ψ1, Ψ2, …, Ψn}. Cкалярное произведение двух функций ρ и τ на границе Г определяется так: , где rj и j - соответствующие значения полярных координат точки ξj, а норма произвольной функции ρ на границе Г считается следующим образом: . Переход к ортонормированному набору осуществляем модифицированным процессом Грама - Шмидта (поскольку классический метод является численно неустойчивым [9]) - построением ненормированных ортогональных функции {Ỹk} и дальнейшей нормализацией - по алгоритму: , , . И далее последовательно по формулам: , , … , и . Получившийся набор функций Ψi ортонормирован на границе Г, а поэтому функцию u можно разложить в ряд из этих функций: . (14) Причем в силу ортонормированности системы {Ψ1, Ψ2, …, Ψn} коэффициенты разложения βk имеют вид . Сравнивая (13) и (14), имеем равенство , выполняющееся для всех точек j = 1, 2, …, N. И поскольку в силу своего построения каждая из функций Ψк (k = 1, 2, …, n) является линейной комбинацией функций {ζ1, ζ2, …, ζn}, то приравнивая в последнем равенстве коэффициенты при одинаковых функциях , получаем выражения для коэффициентов ak, они весьма громоздки. Отметим следующее отличие - в однородном случае наш метод позволял находить возможные значения собственных чисел λ. Как известно [10, 11], если параметр λ в уравнениях (2), (10) не совпадает ни с одним из собственных чисел уравнений (2), (3), то решение уравнений (2), (10) существует и единственно, а если он совпадает с одним из собственных чисел уравнений (2), (3), то решения может не существовать, а если существует, то оно не единственно. Параметр λ, как и граничное условие (10), определяется физикой задачи. Поэтому в модифицированном алгоритме для поиска решения уравнения Гельмгольца с неоднородными граничными условиями нужно опустить этап нахождения λ и сразу находить значения коэффициентов ak. При необходимости следует увеличивать число точек разбиения N или число коэффициентов разложения n - до тех пор, пока сумма (12) не достигнет нужного нам минимального значения. После того, как мы нашли коэффициенты ak с нужной точностью, мы можем определить среднее значение коэффициента усиления по поперечному сечению граничной поверхности: . Продемонстрируем применение данного метода для сечений разной формы с неоднородным граничным условием. Сразу отметим, что в данной работе авторы пока не затрагивают вопросы, насколько нижеприведенные граничные условия (функция u в примерах ниже) соответствуют реальным в лазере, и не рассматривают вопросы определения конкретной граничной функции на практике - предлагается лишь метод решения однородного уравнения Гельмгольца с неоднородными граничными условиями, позволяющий получать решения с высокой точностью при относительно небольшой вычислительной сложности. 3. Случай круглого сечения трубки с неоднородным граничным условием Для указанного сечения можно получить точное решение. Если выбрать за начало отсчета полярной системы координат центр круга, и взять не зависящее от полярного угла граничное условие, то искомая функция тоже не будет зависеть от полярного угла , уравнение (2) превратится в уравнение для функции Бесселя нулевого порядка, и общее решение уравнения (2) будет иметь вид . Зададим функцию u в граничном условии (10) на окружности радиуса r = a в виде и зададим какое-нибудь λ. Данное λ не будет является собственным значением системы (2), (3), если произведение λ и радиуса a окружности не принимает значений , где λk - k-й корень функции J0(x). Можно взять, к примеру, λ = a = 1, тогда очевидным единственным решением уравнений (2), (10) будет функция (и при этом, как и для однородного граничного условия, у нас есть нормировка функции f на единицу в центре сечения). Тогда по формуле (1) получаем следующее выражение для среднего значения коэффициента усиления: . Для проверки обобщения нашего алгоритма мы провели расчеты для круглого сечения трубки. Поскольку у нас есть симметрия при замене угла углом - и углом π + , то в качестве функции следует взять и функция f n(r, ) примет следующий вид: . (15) Среднее по сечению значение коэффициента усиления запишем как , где a0 введен для соответствия однородному случаю, ибо теперь из-за неоднородного граничного условия в общем случае может быть , а из (15) следует, что значение f n(r, ) в центре круга равно a0. Результаты расчета, как и для однородного граничного условия, дают выполнение равенств a0 = 1, a2 = a4 = a6 = a8 = a10 = a12 = 0 с очень высокой точностью, и при этом значение коэффициента усиления k также отлично согласуется с найденным выше теоретическим. 4. Случай эллиптического сечения трубки с неоднородным граничным условием Как и в п. 1 с однородным граничным условием, рассмотрим уравнение (2) в эллиптических координатах (u, ν, z), выберем за начало отсчета центр эллиптического сечения с полуосями a и b, а граничное условие (10) на границе эллипса зададим в виде , где . Из рассуждений п. 1 тогда следует, что решением такой задачи (2), (10) будет функция . Если параметр q такой, что соответствующее λ (связанное с q по формуле ) не является собственным значением системы (2), (3), то это решение единственно. Возьмем, например, эллипс с полуосями b = 1, a = 2 и параметр . В работе [2] были определены собственные значения системы (2), (3) для эллипса с такими полуосями, и данное значение не является собственным. При таких значениях a, b, λ параметр q = 5 и мы можем воспользоваться значениями из п. 1 выше. Коэффициент усиления будет находиться по формуле (6) (при этом уравнение (7), естественно, не выполняется), и результаты численного интегрирования дают k = 0.185k0. Применим теперь наш метод. Поскольку для эллипса опять сохраняется симметрия при замене угла углом - и углом π + , то в качестве функции f n(r, ) можно взять выражение (14) и среднее по сечению значение коэффициента усиления для эллипса будет . При вычислении коэффициентов разложения βk, требующихся для нахождения коэффициентов ak, а также вычисления левой части (11) следует учесть связь цилиндрических и эллиптических координат точек, лежащих на границе эллипса: . Результаты вычислений по нашему методу также дают k = 0.185k0, что может служить подтверждением его правильности. Заключение В данной работе был проверен наш метод решения уравнения Гельмгольца с однородными граничными условиями в случае сечения, имеющего вид эллипса, также этот метод был обобщен на случай неоднородных граничных условий. Преимуществом данного метода является его численная устойчивость и относительно небольшая вычислительная сложность [2]. Проверка на сечениях в виде круга и эллипса показала непротиворечивость этого обобщения. В дальнейшем мы планируем усложнить рассматриваемую модель с учетом распределения интенсивности поля в резонаторе и применить предложенный метод для расчета при других сечениях трубки и граничных условиях [1] с целью найти оптимальное сечение.
Привалов В.Е. // Изв. вузов. Физика. - 2013. - Т. 56. - № 2/2. - С. 246.
Кожевников В.А., Привалов В.Е. // Научно-технические ведомости СПбГПУ. Физико-математические науки. - 2018. - Т. 11. - № 2. - С. 84-95.
Справочник по специальным функциям с формулами, графиками и математическими таблицами / ред. М. Абрамовиц, И. Стиган. - М.: Наука, 1979. - 832 с.
Мак-Лахлан Н.В. Теория и приложения функций Матьё. - М.: ИЛ, 1953. - 476 с.
Таблицы для вычисления функций Матье; собственные значения, коэффициенты и множители связи. - М.: ВЦ АН СССР, 1967.
Frenkel D. and Portugal R. // J. Phys. A: Math. Gen. - 2001. - V. 34. - P. 3541-3551.
Bibby M. and Peterson A. Accurate Computation of Mathieu Functions. Synthesis Lectures on Computational Electromagnetics #32. - Morgan & Claypool Publishers, 2014. - 124 p.
Михлин С.Г. Вариационные методы в математической физике. - М.: Наука, 1970. - 512 с.
Higham N.J. Accuracy and Stability of Numerical Algorithms. - 2nd ed. - Philadelphia: Society for Industrial and Applied Mathematics, 2002. - 711 p.
Владимиров В.С. Уравнения математической физики. - М.: Наука, 1981. - 512 с.
Свешников А.Г., Боголюбов А.Н., Кравцов В.В. Лекции по математической физике. - М.: Изд-во МГУ, 1993.