«Математический синтез оптических наноструктур Учебное пособие МОСКВА Российский университет дружбы народов 2008 г. Аннотация Данное пособие является завершающим в цикле учебнометодических пособий для магистерской ...»

К.П. Ловецкий, Л.А. Севастьянов,

М. В. Паукшто, О.Н. Бикеев.

Математический синтез

оптических наноструктур

Учебное пособие

МОСКВА

Российский университет дружбы народов

2008 г.

Аннотация

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

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

Курс «Математический синтез оптических наноструктур» является составной частью магистерской программы «Оптика наноструктур».

Магистерская программа «Оптика наноструктур» реализуется в рамках направления «Прикладная математика и информатика» и направления «Прикладная математика и физика», а возможно, и других направлений. В составе магистерской программы «Оптика наноструктур» курс «Математический синтез оптических наноструктур» является обязательным. Для других магистерских программ этот курс предлагается по выбору без привязки к семестру или факультативно на усмотрение методической комиссии программы. Курс носит теоретический и практический характер.

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

Задачей курса «Математический синтез оптических наноструктур»

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

Инновационность курса Курс является инновационным по содержанию и по литературе, он включает в себя последние научные достижения в области решения задач объектов не превышают либо сравнимы с длиной волны оптического излучения. Эта область знаний интенсивно развивалась в последнее время, но лишь недавно были созданы устойчивые алгоритмы и разработаны численные методы решения задач для многослойных решеток. Следует дифракционных решеток с характерными размерами больше длины волны оптического излучения устойчивые методы решения известны с середины геометрией, такие как двумерные решетки с произвольным профилем, математические модели взаимодействия излучения с веществом в наномасштабах, а затем с их помощью проектировать новые эффективные устройства в высокотехнологичных областях медицины, энергетики, инфокоммуникаций и приборостроения.

предполагают использование традиционных методик преподавания, использованием кредитной системы оценки знаний.

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

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

сформировалось лишь в последние 10 – 15 лет. Поэтому наблюдается сильный дефицит учебно-методической литературы не только в России, но и во всем мире. Разрабатываемые в рамках инновационной программы «Оптика наноструктур» учебные пособия восполнят в некоторой степени этот пробел и составят основной список литературы для слушателей курсов. Вместе с ними следует использовать несколько учебников и монографий, вышедших в свет к настоящему времени и перечисленных в списке литературы. Курс базируется на публикациях научных статей мировых лидеров исследований в данной области в научной периодике, диссертационных работах их учеников, включающих работы по непосредственному моделированию, дизайну и последующему изготовлению лабораторных образцов оптических элементов и устройств.

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

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

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

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

Методы оптимизации традиционно подразделяются на три основных класса: методы нулевого порядка, методы первого и второго порядков.

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

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

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

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

1.1. Симплекс-метод Нелдера -— Мида Симплекс-метод Нелдера—Мида, впервые опубликованный в г., приобрел с тех пор необыкновенную популярность при решении задач многомерной безусловной оптимизации. Однако, несмотря на такую популярность, простого и подробного доказательства сходимости метода для многих классов функций все еще нет. Приведем [1] некоторые теоретические результаты сходимости метода Нелдера—Мида для строго выпуклых функций одной и двух переменных. Здесь будут даны теоремы сходимости для функций размерности 1 и некоторые (ограниченные) результаты о сходимости в двумерном пространстве. Для более общего представления о свойствах этого метода рассмотрим также контрпример МакКиннона, в котором при выборе специального начального приближения минимизация строго выпуклой функции не приводит к точке минимума.

Симплекс метод Нелдера—Мида [2] стал в последние годы наиболее распространенным методом для решения задач оптимизации без ограничений. Не следует путать его с (возможно?) гораздо более известным методом — симплекс-алгоритмом Данцига для решения задач последовательном построении и использовании симплексов в многомерных пространствах, но во всем остальном – это совершенно различные алгоритмы, ничем между собой не связанные.

Особенно популярен алгоритм Нелдера—Мида при решении прикладных задач, возникающих в физике, химии, медицине. Это, в первую очередь, объясняется его простотой как в программной реализации, так и в понимании. Скорость сходимости метода в практических приложениях сравнима, а часто и превосходит градиентные алгоритмы первого порядка (требующие вычисления первых производных функции по всем аргументам). К тому же по вычислительным затратам он гораздо экономичнее.

Метод Нелдера—Мида предназначен для минимизации функции n действительных переменных с использованием лишь вычисляемых на каждом шаге значений минимизируемой функции (метод нулевого порядка). Для него не требуется и в нем не используется (явно или неявно) информация о производных минимизируемой функции. Таким образом метод Нелдера—Мида относится к общему классу прямых методов поиска минимума функций. Как и многие другие методы, относящиеся к подмножеству прямых методов, метод Нелдера—Мида на каждом шаге итеративного процесса хранит невырожденный симплекс – геометрический объект в n -мерном пространстве ненулевого объема, являющийся выпуклой оболочкой, натянутой на n 1 вершину.

Каждая итерация прямого симплекс-метода поиска минимума начинается с построения симплекса, который задается своими n 1-ой вершинами и вычисляемыми в этих вершинах значениями функции. Затем многогранник дополняется одной либо несколькими точками вместе со значениями функции в них. Одна или несколько вершин после этого отбраковывается. Итерационный процесс завершается тогда, когда вершины симплекса и вычисленные в них значения функции при сравнении с предыдущей итерацией удовлетворяют некоторым условиям сходимости.

Среди алгоритмов такого типа метод Нелдера—Мида, пожалуй, самый экономичный – число вычислений значений минимизируемой функции минимально. На практике обычно требуется одно-два вычисления функции для построения нового симплекса.

Метод Нелдера—Мида весьма широко используется для решения прикладных задач.

Приведем известные [1] результаты о сходимости оригинального метода (сформулированного авторами).

Для функций одной переменной алгоритм сходится к точке минимума. Скорость сходимости M-линейна в том случае, скоростью сходимости подразумевается, что существует такое число M, не зависящее от минимизируемой функции, что диаметр симплекса уменьшается не менее чем в два раза В двумерном пространстве значения минимизируемой функции во всех вершинах симплекса сходятся к одному и тому же значению (в случае стандартного симплексметода).

В двумерном пространстве в случае стандартного метода Нелдера—Мида диаметр симплексов стремится к нулю при Интересно отметить, что последнее утверждение не предполагает, что симплексы сходятся в одну точку. До сих пор не известны случаи, когда бы алгоритм не приводил к окончательному решению в виде одной точки, однако доказательства этого факта пока нет.

1.1.1. Алгоритм Нелдера—Мида Алгоритм Нелдера—Мида был предложен [2] в качестве метода минимизации действительной функции f (x) для x R n. Для полного задания симплекс-метода необходимо определить четыре константы:

(expansion) коэффициент сжатия (shrinkage). Согласно оригинальному алгоритму, эти константы должны удовлетворять условиям:

формулировке, однако в описании реализации алгоритма присутствует.

Практически во всех стандартных реализациях метода коэффициенты выбираются равными:

При рассмотрении одномерного случая воспользуемся общими условиями (1.1), а при доказательстве сходимости для двумерного случая ограничимся условиями (1.2).

1.1.2. Описание алгоритма невырожденный симплекс, состоящий из (n 1) -ой вершины, каждая из которых является точкой в R n. Будем далее предполагать, что каждая итерация начинается с упорядочения вершин x1( k ),..., xnk 1 таким образом, что где через fi ( k ) обозначено значение целевой функции f ( xi( k ) ) в точке xi( k ).

На k -ой итерации строится новый симплекс для k 1 -ой итерации, такой, что. Поскольку мы ищем минимум функции f, то считаем x1( k ) лучшей вершиной, а xnk 1 - наихудшей вершиной. Аналогично, значение f ( xnk 1 ) будем считать наихудшим значением функции.

Опишем одну итерацию симплекс-метода Нелдера—Мида. Не выкладках. Результатом каждой итерации будет:

— единственная новая вершина симплекса, замещающая xnk 1 во множестве вершин при переходе к следующей итерации;

— при проведении сжатия симплекса вычисляется n новых точек, которые вместе с x1 образуют новый симплекс на следующей Опишем одну итерацию алгоритма Нелдера—Мида.

1. Сортировка. Отсортировать вершины, используя правило отсечения элементов (смотри ниже) так, чтобы они удовлетворяли условию 2. Отражение. Вычислить точку отражения x r из условия где x x / n - центр тяжести n лучших точек (т. е. всех вершин то направление отражения признается удачным и делается попытка растянуть многогранник в этом направлении. Для этого вычисляется точка растяжения xe (expansion point) успешным, точку xe добавляется к множеству точек симплекса; иначе (если f e f r ), точка xr. добавляется к множеству точек симплекса. На этом итерация заканчивается и его надо сжать. Сжатие (стягивание, усадка) проводится между x и лучшей из точек xn 1 и xr.

чем xn 1 ) надо провести внешнее сжатие: вычислить симплекса и заканчиваем итерацию, иначе переходим к пункту 5 и проводим операцию усечения.

точек симплекса и заканчиваем итерацию, иначе переходим к пункту 5 и проводим операцию усечения.

5. Усечение. Вычисляем значения функции f в n точках на следующей итерации составляются из точек x1, v2,..., vn 1.

На рис. 1 и 2 продемонстрированы эффекты отражения, растяжения, сжатия и усечения для симплекса в двумерном случае. Исходные симплексы обозначены пунктирной линией.

Рис. 1. Отраженная и продвинутая точки симплекс-метода в Рис. 2. Симплексы метода Нелдера—Мида после внешнего Отметим, что во всех случаях, кроме усечения, новая вершина всегда лежит на прямой, проходящей через точки x и xn 1.

В случае, когда в вершинах симплекса целевая функция принимает одинаковые значения, мы будем упорядочивать вершины (при добавлении Упорядочивание точек без усечения. В том случае, когда итерация заканчивается пунктами 2 или 3 (без поведения усечения), точка xnk исключается из множества вершин на следующей итерации. Введенная же в это множество на k -ой итерации точка, обозначаемая v ( k ), становится новой вершиной и занимает позицию с номером j 1 во множестве все же остальные вершины сохраняют свой относительный порядок и на новой итерации.

Упорядочивание точек после операции усечения. В случае, когда итерация заканчивается усечением многогранника, всего лишь одна точка x1( k ). Здесь задействовано лишь одно правило: если во вновь вычисленном наборе вершин имеется несколько точек, значения целевой функции в которых совпадают, то в новый набор переносится точка x1( k ). То есть если тогда x1( k x1( k ). Кроме того, после усечения обязательно проводится упорядочивание вершин в соответствии с применяемым правилом сортировки.

1.1.3. Свойства метода Нелдера—Мида Опишем основные свойства алгоритма оптимизации Нелдера— Мида, следующие непосредственно из определения итерационного процесса.

1. На каждой итерации метода Нелдера—Мида необходимо одно вычисление функции в новой вершине, если итерация заканчивается шагом 2, два вычисления целевой функции, если итерация заканчивается шагами 3 или 4, и n 2 вычисления целевой функции, если итерация заканчивается шагом 5.

2. Шаг «отражение» назван так по той причине, что отраженная точка относительно точки x (центра тяжести всех вершин, кроме лучшей) приводить к увеличению значений в каждой точке по сравнению с предыдущим шагом, за исключением значения f1 функции в точке x1, т. е. возможно, что fi ( k fi ( k ) для 2 i n 1. Дополнительно можно отметить, что в случае 4a (внешнего сжатия) по алгоритму требуется вычисленное значение функции в точке x r значительно улучшает случае f ( xc ) f ( x1 ) выбирается точка x c, в противном случае – точка xr. Ныне обычно выбирается лучшая из точек xc и xr, если обе дают улучшение ситуации по сравнению с точкой x1.

Удивительно, но строгого доказательства сходимости в R n для столь простого с первого взгляда метода до сих пор нет. Имеются довольно сильные результаты о сходимости в одномерном пространстве для строго выпуклых функций и гораздо более слабые — о сходимости процесса для функций двух переменных [1].

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

1.2. Оптимизация методом наименьших квадратов процедуры линейного поиска получили широкое применение. Эти методы также используются и в подпрограммах нелинейной оптимизации методом наименьших квадратов. В задачах на метод наименьших квадратов квадратов:

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

Эти задачи также широко распространены в теории управления, где в конечном итоге необходимо получить некую y( x, t ), соответствующую Данная задача может быть сформулирована как:

где y( x, t ) и (t ) есть некие скалярные функции.

квадратурных формул уравнение (1.9) может быть сформулировано как задача на метод наименьших квадратов:

где - y и включают в себя веса квадратичной схемы. Отметим, что в данной задаче под вектором F ( x) понимается:

В задачах данного типа невязка F ( x), по-видимому, должна быть наименьшей в точке оптимума, поскольку согласно общепринятой практике необходимо провести как можно ближе к реальной траектории.

Хотя приведенная функция для метода наименьших квадратов (уравнение (1.9)) оптимизации без наличия ограничений, определенные характеристики данной задачи часто могут быть использованы для улучшения итеративной эффективности данной методики решения. Градиент и матрица Гессе для задачи метода наименьших квадратов (1.9) имеют особую структуру.

После обозначения матрицы Якобиана для F ( x) размерностью m n через J ( x), вектора градиента функции f ( x) через G ( x), матрицы Гессе для f(x) через H ( x) и матрицы Гессе для каждой Fi ( x) через H i ( x) получим где F ( x) стремится к нулю при стремлении xk к точке решения, то и сама матрица стремится к нулю. Таким образом, при небольших значениях F ( x) в точке решения одним из наиболее эффективных методов является использование направления Ньютона—Гаусса в качестве основы для процедуры и оптимизации.

1.2.1. Метод Ньютона—Гаусса находится на каждой итерации k решением линейной задачи минимизации по методу наименьших квадратов:

эквивалентом направления Ньютона при пренебрежении члена Q(x).

Направление поиска d k может быть использовано в качестве составляющей стратегии линейного поиска, что обеспечивает условие: на каждой итерации идет уменьшение функции f(x).

Рассмотрим возможные преимущества метода Ньютона—Гаусса. На рисунке 3 представлена траектория поиска минимума для функции Розенброка при использовании постановки задачи как метода наименьших квадратов. Метод Ньютона—Гаусса дает сходимость после 48 обращений к расчету функции при конечно-разностном расчете градиента по сравнению с 140 итерациями для BFGS метода без наличия ограничений.

Рис. 3. Применение метода Ньютона—Гаусса для функции Розенброка В методе Ньютона—Гаусса часто встречается ряд проблем в ситуации, когда член второго порядка Q(x) в уравнении (1.11) достаточно значителен по своей величине. Методом, который преодолевает такие трудности, является метод Левенберга—Марквардта.

1.2.2. Метод Левенберга—Марквардта В основу метода Левенбрга-Марквардта [3-6] положено направление поиска, которое находится при решении системы линейных уравнений:

где скаляр задает как величину, так и направление параметра d k. Когда равен нулю, то направление d k будет идентично этому же параметру из то d k стремится к вектору с нулевыми компонентами и направлению наискорейшего спуска. В данном случае предполагается, что для справедливым. Следовательно, член может быть контролируемым с целью обеспечения спуска в случае необходимости учета членов второго порядка, которые, в свою очередь, заметно ограничивают эффективность метода Ньютона—Гаусса.

Отсюда следует, что метод Левенберга—Марквардта основан на направлении поиска, являющегося сочетанием направления Ньютона— Гаусса и наискорейшего спуска (рис. 4). Решение для функции Розенброка сходится после 90 обращений к расчету функции по сравнению с 48 для метода Ньютона—Гаусса. Такая низкая эффективность отчасти объясняется тем, что метод Ньютона—Гаусса обычно более эффективен в случае, когда в решении невязка равна нулю. Однако такая информация не всегда является заранее доступной, и повышенная устойчивость метода Левенберга—Марквардта компенсирует его иногда имеющую место слабую эффективность.

Рис. 4. Метод Левенберга—Марквардта для функции Розенброка 1.2.3. Реализация метода Ньютона—Гаусса Метод Ньютона—Гаусса реализуется с помощью стратегии полиномиального линейного поиска, аналогичного тому, что было рассмотрено в разделе применительно к оптимизации без ограничений.

При решении линейной задачи методом наименьших квадратов (уравнение (1.8)) путем использования разложения QR для J ( xk ) и подобного разложения для F ( xk ) можно избежать дополнительных возмущений.

Такой подход отличается от применения инверсии явной матрицы J ( xk ) Е J ( xk ), где возможно проявление неожиданных ошибок.

дополнительные меры. Эти меры включают в себя корректировку алгоритма метода Левенберга—Марквардта в случае, когда длина шага обусловленности J ( xk ) будет меньше 10. Под числом обусловленности в данном случае понимается отношение наибольшего сингулярного значения к наименьшему.

1.3. Алгоритм Левенберга—Марквардта Алгоритм Левенберга—Марквардта (Levenberg-Marquardt Algorithm, LMA) является наиболее распространенным алгоритмом оптимизации. Он превосходит по производительности метод наискорейшего спуска и другие методы сопряженных градиентов в различных задачах. Изначально считалось, что LMA – это комбинация простейшего градиентного метода и метода Гаусса-Ньютона, однако впоследствии выяснилось, что данный интервалов. В работе [3] дается интуитивно понятное объяснение LMA.

1.3.1. Постановка задачи LMA решает задачу нелинейной минимизации методом наименьших минимизировать, выглядит следующим образом:

где x ( x1, x2,..., xn ) – вектор, а rj – функция отображения из R n в R.

Функцию rj называют невязкой в предположении, что m n.

Для простоты функция f представляется вектором невязки вида:

Рассмотрим линейный случай, когда каждая функция rj линейна.

Здесь якобиан равен константе, r можно представить как гиперплоскость в получим xmin ( J T J ) 1 J T r, т.е. решение системы нормальных уравнений Возвращаясь к общему (нелинейному) случаю, получим:

если функции rj можно аппроксимировать линейными приближениями (т.е. rj ( x) малы) или если rj ( x) малы сами по себе. Тогда гессиан, как и в линейном случае, будет равен:

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

1.3.2. LMA как комбинация простейшего градиентного метода и метода Ньютона— Гаусса Простейший градиентный метод – это наиболее интуитивно понятный способ нахождения минимума функции. Вычисление параметра на очередном шаге выполняется путем вычитания градиента функции, умноженного на заданный положительный коэффициент:

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

Ситуацию можно улучшить, если учитывать информацию о кривизне и градиенте, т. е. вторые производные. Один из способов сделать это – Раскладывая градиент f в ряд Тейлора вокруг текущего состояния x0, получим квадратичной вблизи x0 ) и решая задачу минимума, приравняв левую часть уравнения (1.18) к нулю, получим правило вычисления параметра на очередном шаге по методу Ньютона:

Поскольку метод Ньютона напрямую использует предположение о квадратичности (пренебрегая членами более высоких порядков при разложении в ряд Тейлора), нет необходимости точно вычислять гессиан, а достаточно использовать его аппроксимацию (1.16). Главное достоинство такого подхода – быстрая сходимость. Однако скорость сходимости зависит от начального положения (если быть более точным - от линейности вокруг начального положения).

Легко заметить, что простейший градиентный метод и метод Ньютона— Гаусса дополняют друг друга с точки зрения предоставляемых преимуществ. Основываясь на этом наблюдении, Левенберг предложил алгоритм, в котором правило вычисления параметра есть комбинация правил (1.17) и (1.19):

где H – матрица Гессе, вычисленная в точке xi. Данное правило используется следующим образом: если на очередной итерации невязка f ( x) работает, и мы уменьшаем (обычно в 10 раз), чтобы понизить увеличивается, необходимо следовать направлению градиента, и мы увеличиваем (во столько же раз).

Таким образом, алгоритм Левенберга представляется в виде последовательности действий:

1. Вычислить параметр на очередной итерации по правилу(1.20).

2. Оценить невязку в новом векторе параметров.

3. Если в результате вычисления параметра невязка увеличилась, вернуться на шаг назад (т.е. восстановить прежние значения весов) и увеличить в 10 раз. Затем повторить выполнение, начиная с шага 1.

4. Если в результате вычисления параметра невязка уменьшилась, принять текущий шаг (т.е. оставить новые значения весов) и уменьшить Недостатком данного алгоритма является то, что если значение велико, вычисленная матрица Гессе никак не используется. Однако можно извлечь некоторую выгоду из второй производной даже в этом случае, масштабируя каждый компонент градиента согласно кривизне. Это должно привести к увеличению шага вдоль направлений, где градиент мал, так что классическая проблема впадины больше не возникнет. Этот ключевой момент был замечен Марквардтом. Он заменил единичную матрицу в формуле (1.20) на диагональ гессиана, получив таким образом следующее правило:

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

Стоит отметить, что хотя LMA является не оптимальным, а лишь эвристическим методом, он очень хорошо работает на практике.

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

1.3.3. LMA как метод доверительных интервалов Изначально LMA был представлен Марквардтом так, как показано в предыдущем разделе, где нахождение минимума производится за счет изменения параметра. Позже было предложено рассматривать LMA как метод доверительных интервалов.

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

В алгоритме, основанном на доверительных интервалах, сначала строится модель m( k ), которая аппроксимирует функцию f на конечном интервале вокруг x( k ). Этот интервал называется доверительным интервалом. На каждой итерации с помощью некоторых эвристик интервал изменяют. Модель m( k ) чаще всего представляет собой квадратичную функцию, полученную разложением функции f в ряд Тейлора вблизи x( k ) :

где H – гессиан или приближение гессиана. Подзадача для нахождения размера очередного шага выглядит следующим образом:

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

Решение задачи (1.23) удовлетворяет следующей теореме Теорема. p* – общее решение задачи Легко заметить, что уравнение (1.24) есть ни что иное, как уравнение наоборот. Следовательно, мы приходим к тому же правилу вычисления параметра, что и в предыдущем разделе.

Эвристика для вычисления размера доверительного интервала обычно зависит от отношения полученного и предсказываемого изменения f, т.е.

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

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

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

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

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

Другой пик интереса к изучению распространения волн в одномерных периодических структурах был стимулирован интенсивным развитием в середине 60-х годов и начале 70-х годов таких областей прикладной физики, как оптическая голография, акустооптика, интегральная и оптическая электроника и др. Обращение вновь к одномерным задачам было в этом случае продиктовано стремлением на новом уровне (определяемом возникшими задачами) и с новыми возможностями (в распространения и рассеяния волн в одномерных структурах с целью практического их использования. И хотя к реальным объектам модель одномерной периодической среды чаще всего неприменима, результаты, полученные на ее основе, полезны для качественного анализа и приближенного описания физики протекающих процессов. Таким образом, одномерная задача приобрела характер ключевой или базовой модели, предназначенной, с одной стороны, для обогащения нашего понимания происходящих в периодической структуре процессов и, с другой стороны, для разработки и тестирования количественных методов их анализа.

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

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

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

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

среде с неоднородной диэлектрической проницаемостью описывается следующим однородным уравнением:

Поскольку среда, как мы договорились, неоднородна только в одном направлении, например, в направлении оси Z, то общее решение векторного уравнения (2.1) всегда может быть представлено спектром частных решений, каждое из которых будет однородным в каком-либо направлении в плоскости, нормальной оси Z. Если для каждого такого решения выбрать систему координат таким образом, что направление, в котором решение однородно, будет совпадать, например, с осью X, то векторное уравнение (2.1) сведется к двум линейно-независимым магнитного полей волны на ось X:

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

Как видно из уравнений (2.2) и (2.3), задача о собственных решениях для электромагнитной волны (задача о собственных волнах) в среде с одномерной модуляцией сводится к двум скалярным уравнениям.

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

Ограничим дальнейшее рассмотрение анализом уравнения:

при гармонической модуляции диэлектрической проницаемости среды Здесь — среднее значение диэлектрической проницаемости; m – глубина модуляции; K — пространственная частота модуляции (K = );

— пространственный период модуляции (рис. 5).

Рис. 5. Среда с одномерной гармонической модуляцией Общее решение уравнения (2.4) будем разыскивать в виде спектра решений с разделяющимися переменными:

При таком выборе решения зависимость от y в любом сечении z = const. представляется интегралом по спектру Фурье.

Подставляя (2.5) в (2.4) и учитывая ортогональность функций exp(i y), получим привычной форме уравнения Матье:

Поскольку в (2.5) может принимать любые действительные пропорциональны постоянным распространения плоских волн — прямой (+b) и обратной (-b) — в однородной среде с диэлектрической Мнимые значения при этом соответствуют неоднородным плоским волнам, т.е. волнам прямым, с экспоненциально уменьшающейся амплитудой и волнам обратным, с экспоненциально нарастающей амплитудой. Следовательно, ветвь b в четвертом квадранте на рис. соответствует прямым плоским волнам (однородным и неоднородным), а ветвь b во втором квадранте – обратным волнам.

В соответствии с теоремой Флоке [8] решение уравнения (2.7) будем разыскивать в виде произведения экспоненты на периодическую функцию с периодом. Последняя может быть представлена рядом Фурье:

Решение (2.8) можно трактовать как волну, состоящую из суммы амплитудами an. Подставляя решение (2.8) в уравнение (2.7), получим однородную бесконечную систему уравнений относительно амплитуд an :

где Условием совместности уравнений системы (2.9) является равенство нулю ее определителя:

Уравнение (2.10) является искомым дисперсионным уравнением, следующему виду [8]:

свое, линейно-независимое от других, k е решение уравнения (2.7).

Общее решение уравнения (2.7) при этом можно записать как сумму линейно-независимых решений с произвольными коэффициентами частного решения.

также будут решениями уравнения (2.11). Заметим, что изменение решению для V ( ) необходимо добавить еще такое же решение со знаком уравнением (2.11) с точностью до слагаемого, равного целому четному соответствует решению уравнения (2.4) в виде плоской волны с противоположном направлении по оси Z.

2.2. Анализ дисперсионного уравнения. Обсуждение физического смысла решений степеням q.

Ограничим бесконечный определитель (2N+1) столбцами и строками Раскрывая определитель, начиная с верхней строки, можно получить следующий результат:

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

следующее значение:

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

уравнение (2.11), получим Из уравнения (2.14), учитывая малость q, для всех b, исключая Знаки перед квадратными корнями в выражениях (2.16) – (2.18) минимальной. В окрестности b = 1 из этого условия следует, что при q знак в (2.16) должен быть «минус», если b > 1, и «плюс», если b < 1. Можно показать, что противоположные знаки соответствуют экспоненциально убывающими решениями. При q 0, как это видно из (2.17), границей между гармоническими и экспоненциально убывающими рассчитанные по формулам (2.15) – (2.18), представлены на рис. 7.

Остановимся на физическом смысле полученных решений. Для определенности допустим, что поле внутри модулированной среды возбуждается на плоской границе модулированного и немодулированного непрерывности тангенциальных компонент поля следует, что решения среды вдоль оси Y, можно утверждать, что такого вида решение вдоль оси Y будет сохраняться в любом другом сечении Z.

Качественно картину распространения волн в модулированной среде при Z > 0 среда однородна. Решение для поля в этом случае представляет величина b является мнимой, что соответствует полному внутреннему является экспоненциально убывающим.

В первом приближении будем считать, что рассеяние плоской волны на гармонической неоднородности является очень слабым, и поэтому амплитуда которой определяется нулевым приближением. Гармоническую неоднородность можно рассматривать как две встречных плоских нулевыми фазовыми скоростями. В результате рассеяния электромагнитной волны на этих двух волнах неоднородности возникают две вынужденные волны с компонентами волновых векторов очевидно, должны быть пропорциональны q. Таким образом, полное поле в первом приближении представляет собой плоскую рассеиваемую волну, которую в дальнейшем будем называть нулевой пространственной гармоникой, и две рассеянные плоские волны, которые мы будем называть плюс и минус первыми пространственными гармониками полного поля в модулированной среде (рис. 8).

На рис. 8а приведено схематическое изображение возбуждения модулированной среды плоской волной, падающей из немодулированного диэлектрика. На рис. 8б изображены направления распространения нулевой и плюс первой гармоник при b 0 и b 1. На рис. 8в представлены направления распространения нулевой и плюс первой гармоник при b = 1, а на рис 8г - направления распространения нулевой и плюс первой гармоник при b = 0.

Рис. 8. Геометрия волновых векторов решетки.

В следующем, втором приближении необходимо учесть вторичное рассеяние, т.е. рассеяние +1 и -1 гармоник на неоднородности. При вторичном рассеянии возникают плоские волны с постоянными пространственные гармоники с амплитудами, пропорциональными q 2, и, соответствует учету обратного влияния рассеянного поля на рассеиваемое поле. Это влияние приводит к сдвигу постоянной распространения нулевой гармоники относительно значения, определяемого в нулевом приближении, на величину, пропорциональную q 2. В нашем решении это имеет место для всех b, за исключением точек b 0, b 1 (см. выражение В окрестности точки b 1 значение волнового вектора -1 гармоники по абсолютной величине оказывается равным постоянной распространения нулевой гармоники, а направление его проекции на ось Z обратно соответствует условию Брэгга, при котором поля, отраженные от каждого периода неоднородностей, складываются синфазно. В этом случае амплитуда -1 гармоники (отраженной волны) оказывается максимально возможной и равной амплитуде рассеиваемой (отражаемой) нулевой гармоники.

Решение для в окрестности этой точки оказывается комплексным (см. (2.16)). Мнимая часть соответствует экспоненциальному убыванию амплитуды поля нулевой пространственной гармоники, которое можно объяснить качественно тем, что при отражении от каждого периода неоднородности часть мощности падающей волны (нулевой гармоники) рассеянное поле на том же периоде оказывается точно в противофазе с полем падающей волны (при b = 1). Следовательно, при обратном рассеянии нет возврата энергии из отраженной волны в падающую, как это имеет место при b 1, и амплитуда падающей волны на каждом периоде неоднородности уменьшается на величину, пропорциональную q. При небольшой расстройке от b = 1 экспоненциальная зависимость сохраняется в области, ширина которой равна q (см. рис. 7). Полоса b в окрестности b = 1, в которой Im 0, является полосой брэгговского отражения частоты и угла падения, то, следовательно, существуют частотная и угловая полосы брэгговского отражения, пропорциональные глубине модуляции диэлектрической проницаемости. Заметим, что действительная часть постоянной распространения нулевой гармоники внутри полосы брэгговского отражения остается постоянной и равной.

проницаемости кроме полосы брэгговского отражения первого порядка при b = 1 существуют полосы отражения второго (в окрестности b = 2), третьего (при b = 3) и т.д. порядков. Мнимые части постоянных распространения и ширина полосы отражения в окрестностях этих точек Дисперсионное уравнение (2.14), учитывающее только первый порядок малой величины 2q, описывает брэгговское отражение только первого порядка. Для описания более высокого порядка в дисперсионном вычислении ( 0).

Приведенная выше физическая картина объясняет возбуждение в модулированной среде одного из решений волнового уравнения, которое описывается набором пространственных гармоник. Это решение при пространственных гармоник, кроме амплитуды нулевой, стремятся к нулю, а постоянная распространения нулевой гармоники стремится к постоянной распространения плоской волны в однородном диэлектрике. Это решение называется прямым нулевым. Возбуждение нулевого обратного и решений с другими номерами (прямых и обратных) можно пояснить на примере модулированной среды, ограниченной по оси Z длиной L (рис. 9).

Пусть слева на границу модулированной среды при Z = 0 под углом падает плоская волна. Эта волна возбуждает в модулированной среде спектр пространственных гармоник (нулевое решение), каждая из которых распространяется со своей постоянной распространения в положительном направлении оси Z. На другой границе при Z = L происходит частичное отражение, при котором компоненты спектра нулевого решения являются источниками обратных решений – нулевого, ±1 и т.д. При частичном отражении обратных решений от границы при Z = 0 кроме нулевого решения возбуждаются ±1, ±2 прямые решения и т.д. Амплитуды решений определяются из граничных условий при Z = 0 и Z = L, а соотношения между амплитудами гармоник в каждом решении — значениями b и q, т.е.

параметрами модулированной среды.

Рис.9. Ограниченная модулированная среда и направление волновых векторов пространственных гармоник в решениях с различными 2.3. Амплитуды пространственных гармоник Решение задачи о распространении электромагнитных волн в гармонически модулированной среде, как это показано выше, сводится к уравнению Матье, частными решениями которого являются функции:

представлены рядом Фурье:

где akn - амплитуды спектра Фурье, в рассматриваемой задаче имеющие смысл амплитуд пространственных гармоник k-го решения.

Значения вычисляются из однородной системы (2.9) при 2k и b, определяемом дисперсионным уравнением. Из системы (2.9) akn определяются с точностью до постоянного множителя. Если разыскивая в виде разложения в ряд по степеням q, для akn, учитывая слагаемые с наименьшей степенью q, имеем:

где S = 1, 2, 3, …. При S = 0 akk = 1.

Вычислим akn для частного случая, когда b считается заданным. Для уравнения.

определяется выражением (2.15). Подставляя (2.15) в (2.19), имеем:

В окрестности b = 1, учитывая (2.16) Как видно из выражений (2.20)-(2.22), амплитуды пространственных гармоник k-го решения убывают с ростом S – номера гармоники является лишь окрестность b = 1, в которой амплитуда (k-1)-ой гармоники k-го решения равна амплитуде k-ой гармоники того же решения. На самом Здесь знак плюс перед квадратным корнем соответствует b 1, а знак q является полосой брэгговского отражения. В пределах этой абсолютной величине.

Рис. 10. Зависимость отношения амплитуд (k-1)-ой и k-ой гармоник от расстройки относительно брэгговского условия Зависимость ak,k от расстройки b 1 в окрестности брэгговского условия b 1 получена в предположении, что величина b является заданной. Как уже отмечалось, дисперсионное уравнение определяет лишь связь между b и, и вопрос о том, какую из этих величин можно считать заданной, может быть решен с учетом условий возбуждения поля в периодически модулированной среде. Одним из возможных вариантов возбуждения поля, при котором заданной можно считать величину b, является возбуждение периодически модулированной среды на плоской границе, нормальной вектору решетки (см. рис. 8).

Из анализа выражений для амплитуд спектра Фурье любого k-го диэлектрической проницаемости ( q 1) и при достаточно большой электромагнитного поля в модулированной среде являются почти гармоническими. В окрестности брэгговского условия, где амплитуды k-ой и (k -1)-ой гармоник k-го решения сравнимы по величине, а амплитуды остальных L – в виде плоской прошедшей волны. В модулированной области 2 (0< Z < L) поле записывается в виде прямого и обратного решений, каждое из которых Неизвестные амплитуды отраженной и прошедшей волн, а также амплитуды прямого и обратного решений в модулированной области находятся из условий непрерывности компонент поля, тангенциальных границам раздела областей.

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

Итак, запишем поле в области 1 ( при Z L) где Aпрош - амплитуда прошедшей волны в точке y = 0, z = L.

Удовлетворяя граничным условиям:

получим следующую систему алгебраических уравнений относительно неизвестных амплитуд A, B, Aотр и Aпрош :

KL KL KL

KL KL KL

KL KL KL

KL KL KL

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

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

Для амплитуд прямого и обратного решений в модулированной среде из (2.28) получаем следующие выражения:

пропускания – минимальные.

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

Колебания коэффициента отражения за пределами брэгговского отражения ( 1) можно объяснить следующим образом. Поле в пространственной частотой, равной разности постоянных распространения указанных гармоник - qK 1. При изменении частоты или угла падения период биений изменяется. Если на длине зеркала укладывается нечетное число четвертей биений, то коэффициент отражения принимает максимальное значение. Ситуация здесь аналогична той, которая имеет место при отражении от плоскопараллельной пластинки с толщиной, кратной нечетному числу четвертей длин волн, когда отражения от краев пластинки складываются в фазе. Поскольку амплитуда биений тем больше, чем ближе к единице, то первый максимум коэффициента отражения за значения Г m убывают пропорционально Значения соответствуют таким значениям частоты или угла падения падающей волны, при которых на длине зеркала укладывается целое число полупериодов биений. При этом коэффициент отражения оказывается равным нулю.

Рис. 12. Зависимость коэффициента отражения брэгговского зеркала от расстройки при различных значениях параметра Следует обратить внимание на распределение поля при расстройках в модулированной области (2.31) и (2.32), имеем Подставляя эти выражения в выражение (2.26) для поля в Первое слагаемое в фигурных скобках (2.33) представляет собой стоячую волну, амплитуда которой изменяется вдоль z по гармоническому закону. При этом амплитуда стоячей волны равна нулю на краях зеркала.

Второе слагаемое описывает бегущую волну с амплитудой, медленно меняющейся вдоль z. Из приведенного выражения видно, что в пучностях стоячей волны амплитуда поля в модулированной области всегда превышает амплитуду падающей волны. Особенно четко это заметно при приблизительно равен амплитуде падающей волны, а пучности стоячей волны в q раз больше. Это соответствует значительной запасенной энергии внутри модулированной среды. При этом запасенная энергия при n = 1 оказывается наибольшей. Значительное накопление энергии в модулированной области при позволяет рассматривать значения как резонансные, а участок модулированной среды — как проходной резонатор. В точках резонанса коэффициент пропускания такого резонатора, как и любого другого проходного резонатора без потерь, равен единице.

Тема 3. Метод Тихонравова синтеза многослойного зеркала В работе [9] на основе идей регуляризации дается корректная постановка задачи синтеза многослойных оптических систем с заданными коэффициентами пропускания. Приводятся примеры синтеза некоторых практически важных систем.

последнее время все большее внимание исследователей. Новые возможности для изучения и решения этих задач дает принцип регуляризации [10-12].

К этому же классу относится задача о синтезе оптических систем с заданными характеристиками (энергетическими коэффициентами пропускания и отражения). Такие системы могут быть получены путем создания многослойных покрытий, наносимых в виде тонких пленок на поверхность «подложки» (см., например, [13]). Достаточно полно были исследованы системы с небольшим числом слоев равной оптической толщины (произведения показателя преломления на толщину слоя).

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

Работа [9] посвящена именно этим вопросам, рассматриваемым на диэлектрических пленок при нормальном падении световой волны.

3.1. Физическое содержание рассматриваемой задачи синтеза Рассмотрим пластинку 0 z H с показателем преломления n n( z ), ), на которую падает нормально монохроматическая световая волна длины 2 c /. Амплитуда электрического поля волны E ( z) определяется условиями следующей краевой задачи:

где E0 - амплитуда падающей волны. Если n( z ) – кусочно-непрерывная функция, то задача дополняется однородными условиями сопряжения, соответствующими третьей строке (3.1), на внутренних границах.

Пропускательная способность системы в рассматриваемом случае дается «коэффициентом пропускания»

однозначно определяемым для заданной системы как функция длины рассматривать коэффициент отражения R 1 T.

Рассматриваемая нами задача является обратной по отношению к задаче (3.1) — (3.2) и состоит в том, чтобы по заданной априори функции T ( ) определить n( z ).

Однако такая постановка задачи не является корректной, что характерно для задач синтеза, как это отмечалось в [11], по следующим имеющей требуемый коэффициент пропускания. Во-вторых, данному коэффициенту пропускания (например, близкому в каком-то смысле к T ( ) ) могут отвечать различные системы [12]. Наконец, не всякая система практически реализуема, что также должно быть принято во внимание при постановке задачи.

Мы будем рассматривать задачу синтеза в классе кусочно-постоянных функций n( z ), т. е. для слоистых систем, при условиях типов: показатели преломления слоев, во-первых, имеют определенные для каждого слоя значения ( n j для j го слоя, j 1,2,..., N ), во-вторых, варьируются в заданных пределах: nmin n nmax ; при этом толщины слоев произвольны;

N не является априори заданной величиной.

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

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

Ниже мы рассматриваем случай нормального падения световой волны. Однако это не ограничивает общности постановки задачи. В самом деле, в случае наклонного падения волны достаточно заменить n j на Аналогично этому могут быть рассмотрены системы с поглощением.

3.2. Математическая постановка задачи синтеза Приводимая ниже постановка задачи синтеза оптических систем базируется на принципах, изложенных в [11].

толщины и показатели преломления слоев. Пусть T ( ) - заданная на некоторую функцию T ( ) (см. (3.2)). Тем самым определен нелинейный оператор A :

пропускательной способностью системы слоев, описываемой вектором x.

Пусть E2 N есть 2N -мерное векторное пространство, D2 N - замкнутое выпуклое множество в нем, определяемое так:

Примем за меру близости T ( ) и T ( ) величину системах.

Ясно, что максимально достижимые точности удовлетворяют оценкам, которую назовем предельно достижимой точностью.

Основная задача состоит в том, чтобы приблизить в смысле (3.3) требуемую характеристику с некоторой априори заданной точностью на синтезируемую структуру. Такими дополнительными требованиями являются, в первую очередь, требование минимальности слоев системы N, требованиями устойчивости к внешним воздействиям (при большом числе слоев системы легко разрушаются), возможностями напылительной аппаратуры.

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

множестве D2 N. Причем параметр выбирается из условия существует, так как D2 N - замкнутое выпуклое множество конечномерного пространства. Найденное таким образом решение будет минимальным по 3.3. Дискретная формулировка задачи синтеза В случае системы с произвольным n( z ) для вычисления A( x, ) и явного (через A( x, ) ) выражения его производной можно использовать алгоритм, аналогичный изложенному в [15].

характерных для практики систем с кусочно-постоянной n( z ). В этом представлений оператора A( x, ), принятая в работах по изучению оптических систем [13]. Соответственно где n0, nN n( H ) - заданные показатели преломления обрамляющих покрытие сред (рис. 13), mik - элементы полной характеристической матрицы M ( x), определенные формулой где - угол преломления на границе j -го слоя, легко рассчитываемый Операторная производная A ( x, ), которая нам потребуется для минимизации функционалов, может быть представлена аналитически.

Положим, не меняя обозначений, A( x, ) T (M ( x)), где T - оператор, отображающий матрицу M ( x) в функцию T ( ). Тогда [16] Частная производная матрицы M ( x) по d j будет иметь вид:

а частная производная по n j :

Из (3.6) нетрудно получить, что производная оператора A( x, ), взятая при значении M ( x) и примененная к матрице будет иметь вид:

Если L - одна из матриц или M унимодулярность [13] матрицы M ( x), легко получить Из (3.7) получаем, что градиент первого члена функционала находимыми по формуле где mik,1 и mik,2 - соответственно элементы матриц M Таким образом, F ' ( x) будет представлять собой 2N -мерный вектор с координатами 3.4. Алгоритм минимизации функционала и некоторые результаты расчетов Как показал предварительный анализ, функционал (3.4) обладает большим количеством локальных минимумов, к которым ведут «глубокие наклонные овраги», в результате чего движение к минимуму по методу градиентного спуска происходит как бы по «ломаной линии», направления и величины отрезков которой через один близки друг к другу; чем ближе к минимуму, тем более заметным становится это свойство. В соответствии с такой структурой функционала был построен алгоритм его минимизации, включающий поиск локальных минимумов по случайно заданным начальным приближениям, отбор «подозрительных» на глобальный минимум и уточнение их значений. Полученный таким образом глубокий минимум считается глобальным.

производился по обычному методу проекций градиента [17], при котором n -й итерационный шаг осуществляется по формуле где PD2 N -проектор на множество D2 N ; B( x) F ' ( x), т. е. равен градиенту определяемый из условия минимальности функционала (3.4) по данному направлению градиента.

Процесс уточнения, в соответствии с описанной выше структурой функционала, отличался от (3.9) тем, что для итераций с номерами, кратными трем, Рис. 14. Коэффициенты пропускания и отражения для 11-слойной модели при начальных толщинах Как нетрудно заметить, эта поправка к формуле (3.9) отвечает направлению движения вдоль оврага, а не от «стенки к стенке», как обычно, что позволяет в нашей задаче значительно сократить время счета и существенно ближе подойти к минимуму.

При практической реализации вначале рассматривалась постановка задачи, при которой показатели преломления слоев фиксированы и принимают чередующиеся значения: n1 - для нечетных слоев, n2 = для четных.

На рис. 14 и 15 изображены коэффициенты пропускания и отражения в широком диапазоне длин волн (расчеты проводились с помощью программы MorphoVision [18]) для так называемой отрезающей системы (заданная характеристика, которую необходимо приблизить, такова: T 0, при 6мк, T 1 при 6 мк) в инфракрасной области. Сплошная линия – коэффициент пропускания для оптической системы, штриховая – для коэффициента отражения. Требуемой точности (для синтезированной системы 10 2 ) удалось достичь на 11-слойной системе.

Рис. 15. Коэффициенты пропускания и отражения для 11-слойной Минимизация по суммарной толщине при этом не производилась.

начального приближения ( d H ) и показатели преломления слоев приведены в табл. 1. Основные трудности получения таких систем заключаются в том, что интересующая нас область длин волн весьма широка (от 2-4 до 10 мк), и существующие полуэмпирические методы, достаточно полный обзор которых приведен в [19], не позволяют задать начальное приближение, по которому можно было бы с требуемой точностью получить нужный коэффициент пропускания.

эффективным (однако в нем не учтены дисперсия и поглощение среды).

Для синтезированной системы значение уменьшилось при минимизации в 165 раз – от 1.98 для начальной системы до 0.012 для синтезированной.

составляющих функционала (3.7) в зависимости от величины. Заданная характеристика, которую необходимо приблизить, такова: в области от до 6 мк должно быть полное отражение, в области от 6 до 10 мк должно Показатели преломления для четных слоев 4.00, для нечетных — 1.34.

небольшом увеличении (т. е. можно существенно уменьшить приближения). Заметим, что при дальнейшем увеличении погрешность, наоборот, быстро возрастает, в то время как суммарная толщина уменьшается незначительно. Таким образом, рассматриваемая методика оказывается эффективной и при решении задачи «оптимального» (в указанном выше смысле) синтеза.

Было исследовано влияние на точность реализации требуемой характеристики пропускания варьирования показателей преломления наряду с варьированием толщин. Для этого был проведен синтез 9-слойной отрезающей системы с областью отражения 3-6 мк и областью пропускания 6-10 мк. При этом сначала показатели преломления были фиксированы и принимали для нечетных слоев максимально возможное в этой области длин волн значение n1 4.00, а для четных – минимально возможное n2 1.34 (первое вещество – германий, второе – криолит являются основными веществами, используемыми при практическом конструировании систем, работающих в микронном диапазоне).

Варьировались только толщины. Затем был проведен поиск при варьировании также и показателей преломления в пределах, определяемых величинами n1 и n2. Хотя при этом значения показателей преломления для некоторых слоев полученной системы и приняли значения, отличные от n или n2 (отметим, что для четырех из девяти слоев n приняло одно из крайних значений: n1 или n2 ), точность приближения повысилась лишь очень незначительно (величина уменьшилась лишь в 1.07 раза). Видимо, это связано с тем, что требуемая характеристика имеет широкую область отражения, а для получения широкой области отражения, как это замечено на практике (в случае систем со слоями равной оптической толщины это показано аналитически), необходимо, чтобы показатели преломления соседних слоев отличались друг от друга как можно больше.

Рис. 16. Коэффициенты пропускания и отражения просветляющей 3-слойной системы. Начальная система слоев Рис. 17. Коэффициенты пропускания и отражения просветляющей системы. Оптимизированные толщины слоев Но при синтезе систем с характеристиками другого вида варьирование показателей преломления дает значительный эффект. На рис. 16, изображены коэффициенты пропускания T и отражения R начальной (рис. 16) просветляющей системы и коэффициенты пропускания T и отражения R соответствующей ей полученной (рис. 17) в области длин волн от 2 до 8 мк. Просветление в столь широкой области обеспечивает 2слойное покрытие с показателями преломления и толщинами 1-го и 2-го среды n0 1.00, подложки n 1.56 ). Данная система была получена при синтезе 3-слойных систем. Толщины слоев начальной системы (см. рис.

16) задавались случайным образом, показатели преломления 1-го и 3-го слоев – равными 4.00, а второго 1.34. Интересно отметить, что при минимизации толщина первого слоя обратилась в нуль (система стала двухслойной), а показатель преломления 3-го слоя (для полученной системы – второго) изменился с 4.00 до 1.51.

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

Тема 4. Теоретические аспекты задачи синтеза оптических покрытий 4.1. Качественные закономерности структуры оптических покрытий В течение многих десятилетий актуальной является проблема создания композиционных материалов с новыми свойствами. Наиболее распространенный тип неоднородных структур – слоистые среды, свойства которых меняются вдоль одной координаты. При волновом воздействии на слоистую структуру возникает система отраженных и преломленных волн, которые, взаимодействуя с падающей волной, образуют сложную интерференционную картину. Меняя структуру слоистой среды, можно в широких пределах управлять волновой картиной процесса. На этом основана работа многих приборов и устройств в различных областях физики и техники. Велика роль слоистых структур, применяемых для управления энергетикой волнового процесса.

Слоистые структуры широко применяются в оптике (просветление оптических поверхностей, создание интерференционных зеркал, светоделителей, отрезающих фильтров и т.п.), радиофизике (для преобразования волн в СВЧ- и радиодиапазонах), квантовой электронике, акустике (в ультразвуковой технике, для решения задач звукоизоляции, в качестве согласующих устройств), в теплофизике (для решения задач температурной стабилизации). Упругие слоистые среды применяются в качестве фильтров продольных и поперечных волн, для решения задач вибрационной защиты, создания эффективных демпфирующих устройств.

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

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

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

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

Задачи синтеза как обратные задачи математической физики обладают свойствами некорректности по А.Н. Тихонову. Анализ существующих методов оптимизации и синтеза неоднородных структур при волновых воздействиях показывает, что большинство их них малоэффективно или Многоэкстремальность волновых задач синтеза существенно затрудняет исследование предельных возможностей и характеристик конструкций.

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

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

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

взаимодействующего со слоистой средой, структура которой задается элементом u U. Связь между структурой слоистой среды и ее волновыми свойствами определяется дифференциальным оператором L :

Пусть g - заданная функциональная характеристика волнового поля.

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

на решениях операторного уравнения (4.1) при ограничении:

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

энергетического коэффициента пропускания от частоты:

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

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

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

В соответствии с этим в рамках существующих подходов:

1. Отсутствует возможность объективно оценить, насколько возможности созданных структур, функционирующих в различных областях физики и техники, отличаются от предельно достижимых.

2. Отсутствует возможность эффективного конструирования структур, реализующих предельные возможности.

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

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

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

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

оптимальности в волновых задачах синтеза установлены уникальные свойства оптимальных структур [20, 21].

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

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

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

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

В общей постановке задача оптимального синтеза слоистой структуры с требуемыми свойствами может быть сформулирована следующим образом:

где u - вектор-функция, описывающая структуру оптимизируемой слоистой среды; U - функциональное множество вектор-функций u, соответствующих возможным допустимым параметрам в конструкции; g вектор-функция, компонентами которой является совокупность функциональных характеристик, интересующих проектировщика.

Связи между u и g в общем виде описываются операторным уравнением вида (4.7). Критерий качества характеризует меру близости функциональных характеристик проектируемой структуры к требуемым g.

Задача оптимального проектирования в рассматриваемой постановке заключается в поиске множества решений U*, доставляющих глобальный минимум функционалу (4.2):

т.е. в нахождении всех квазирешений уравнения: Lu g.

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

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

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

4.2. Наклонное падение волны на слоистую структуру Будем считать, что пространственная неоднородность достигается за счет набора слоев из различных однородных и изотропных материалов.

Основой для описания распространения электромагнитных волн в неоднородных средах служат уравнения Максвелла.

Будем рассматривать случай наклонного падения электромагнитной волны под углом на многослойную конструкцию общей толщины l, состоящую из N плоскопараллельных слоев с различными физическими свойствами. Плоскость падения совпадает с плоскостью xz. Ось z направлена в сторону падения, наружная поверхность многослойника, на которую падает волна, совпадает с плоскостью z 0. При таком выборе системы координат электрический и магнитный векторы не будут зависеть от координаты y.

Уравнения Максвелла допускают два независимых решения. В одном из них электрический вектор E перпендикулярен плоскости падения (электрическая волна с горизонтальной поляризацией, ТЕ-волна), в другом электрический вектор E параллелен плоскости падения (электромагнитная волна с вертикальной поляризацией, ТМ-волна). Любое решение системы уравнений Максвелла можно представить в виде линейной комбинации (с комплексными коэффициентами) ТЕ- и ТМ-волн. Отношение коэффициентов этой комбинации определяет поляризацию волны. В соответствии с этим достаточно рассмотреть случай ТЕ- и ТМ-волн.

4.2.1. Наклонное падение ТЕ-волны на систему непоглощающих магнитоэлектрических слоев ТЕ-волна общего вида может быть представлена суперпозицией плоских монохроматических волн, т. е. в виде интеграла Фурье:

полупространства, их которого приходит волна, c - скорость света. Можно магнитоэлектрических слоев сводится к нахождению решения следующей краевой задачи [8,9]:

Здесь - угол, под которым электромагнитная волна выходит из слоистой структуры:

диэлектрической и магнитной проницаемости по толщине слоистой связанными между собой функциональной зависимостью, позволяющей рассматриваемом случае независимым физическим параметром будет являться только. Пусть задан дискретный набор материалов, множество каждого z [0,1] выполнено включение Энергетический коэффициент пропускания ТЕ-волны определяется через решение краевой задачи (4.9):

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

Математически данная задача формулируется следующим образом:

необходимо определить распределение диэлектрической проницаемости значение критерий энергетического коэффициента пропускания от частоты в заданном синтеза (4.9)—(4.11):

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

Оказывается, что конструктивный анализ необходимых условий оптимальности в волновых задачах оптимального синтеза позволяет установить ряд важных качественных закономерностей структуры оптимальных конструкций. А именно:

1) получить верхние оценки для числа различных материалов дискретного набора, которые могут входить в оптимальную конструкцию;

2) получить оценки оптимального числа слоев конструкции;

3) построить систему рекуррентных соотношений, позволяющих априори, до проведения численных расчетов, выделить именно те материалы допустимого набора, которые могут входить в оптимальную конструкцию;

4) установить характер сочленения слоев с различными физическими свойствами в оптимальной конструкции;

5) получить условия, при выполнении которых в оптимальную конструкцию может входить только два материала допустимого набора независимо от количества материалов, составляющих исходный дискретный набор;

6) установить свойство внутренней симметрии во взаимосвязи параметров в оптимальных структурах.

Поэтому справедливо следующее утверждение.

Утверждение 1. Число различных материалов, составляющих оптимальную многослойную систему, не может превышать p q, где p и q - число точек разрыва функций ( )и ( ).

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

Утверждение 2. Физические свойства материалов слоев оптимальной соотношений:

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

Рассмотрим последовательность элементов Утверждение 3. Физические параметры материалов, входящих в состав оптимальной многослойной структуры, являются соседними в последовательности (4.15).

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

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

эффективное сужение множества допустимых вариантов. Поэтому их применение открывает новые возможности для повышения эффективности различных методов поиска оптимального решения.

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

4.2.2. Наклонное падение ТМ-волны на систему непоглощающих магнитодиэлектрических слоев позволяет для случая наклонного падения ТM-волны получить результаты, аналогичные результатам, полученным для случая наклонного падения TEволны. ТМ-волна общего вида также может быть представлена в виде суперпозиции плоских монохроматических волн, т. е. в виде ряда Фурье:

где q( z, ) - спектральная плоскость ТМ-волны. Также можно показать, следующей краевой задачи [22]:

пропускания ТМ-волны определяется через решение краевой задачи (4.16) следующим образом:

Задача оптимального синтеза в вариационной постановке заключается в минимизации критерия на решениях системы (4.16).

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

Утверждение 6. Физические параметры соседних слоев оптимального последовательности (4.19).

4.2.3. Наклонное падение монохроматической электромагнитной волны с вертикальной поляризацией (ТМ-волны) на систему непоглощающих магнитодиэлектрических слоев Для этого случая может быть доказано утверждение, аналогичное утверждению 4.

Утверждение 7. В случае наклонного падения монохроматической электромагнитной волны с вертикальной поляризацией на систему непоглощающих магнитодиэлектрических слоев оптимальная конструкция состоит не более чем из двух материалов независимо от количества проницаемости которых являются соседними в последовательности (4.19).

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

Качественные закономерности оптимальных покрытий установлены для случая поляризованных магнитных электромагнитных волн ТЕ- и ТМтипа, а также для случая монохроматических волн.

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

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

Фотонные кристаллы по характеру изменения диэлектрической проницаемости можно разделить на три основных класса.

Рис. 17. Схематическое представление одномерного фотонного периодически изменяется в одном пространственном направлении (рис.

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

Рис. 18. Схематическое представление двумерного фотонного 2. Двумерные, в которых диэлектрическая проницаемость периодически изменяется в двух пространственных направлениях (рис. 18). На рис. диэлектрическими проницаемостями 1, которые находятся в среде с кубической решетке. Такие фотонные кристаллы могут проявлять свои свойства в двух пространственных направлениях, и форма областей с прямоугольниками, как на рис. 18, а может быть любой (окружности, эллипсы, произвольная и т.д.). Кристаллическая рештка, в которой упорядочены эти области, также может быть другой, а не только кубической, как на приведнном рис. 18.

периодически изменяется в трех пространственных направлениях.

Такие фотонные кристаллы могут проявлять свои свойства в трех пространственных направлениях, и их можно представить как массив объемных областей (сфер, кубов и т.д.), упорядоченных в трехмерной кристаллической рештке.

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

Методы для временной области удобны в отношении динамических задач, которые предусматривают временную зависимость электромагнитного поля от времени. Они также могут быть использованы для расчета зонных структур фотонных кристаллов, однако практически сложно выявить положение разрешенных и запрещенных зон в выходных данных таких методов. Кроме того, при расчете зонных диаграмм фотонных кристаллов используется преобразование Фурье, частотное разрешение которого зависит от общего времени расчета метода. То есть для получения большего разрешения в зонной диаграмме нужно потратить больше времени на выполнение расчетов. Есть еще и другая проблема временной шаг таких методов должен быть пропорционален размеру пространственной сетки метода. Требование увеличения частотного разрешения зонных диаграмм требует уменьшения временного шага, а, следовательно, и размера пространственной сетки, увеличения числа итераций, необходимой оперативной памяти компьютера и времени расчета.

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

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

5.2. Распространение света в фотонных кристаллах 5.2.1. Спектр отражения и пропускания одномерных фотонных кристаллов Метод матриц распространения для расчета отражения и пропускания света в одномерных фотонных кристаллах E0 exp(i(k0 r характеризующихся толщиной d j и, в общем случае, комплексной диэлектрической проницаемостью " или коэффициентом преломления является плоскостью падения волны. Прошедшее и отраженное от многослойной структуры электромагнитное поле является результатом многолучевой интерференции в каждом из ее слоев. Однако в силу линейности волнового уравнения в каждом из слоев поле является положительном и отрицательном направлениях оси z (далее называемыми прямой и обратной волнами):

сохраняется постоянной внутри многослойной структуры, что является следствием трансляционной симметрии вдоль слоев, а нормальная Амплитуды прямой и обратной волн E j и E j в уравнении (5.1) являются комплексными величинами вследствие многолучевой интерференции.

Ограничиваясь рассмотрением стационарного случая, множителями exp(ik x x i t ) далее будем пренебрегать.