Компьютерные методы анализа моделей прикладной электродинамики и нелинейной оптики

17. Денисов A.M., Макеев А. С. Итерационные методы решения обратной задачи для одной модели популяции // ЖВМиМФ. 2004. 44. № 8. С. 1492-1501.

18. Денисов A.M. Существование и единственность решения системы интегральных уравнений первого рода // Диф. ур-ния. 2003. 39. № 9. С. 1201-1208.

19. Krylov A.S., Liakishev А. V. Numerical projection method for inverse Fourier type transforms and its application // Numerical Functional Analysis and Optimization. 2000. 21. N 1. P. 205-216.

20. Blagonravov L.A., Skovorod'ko S.N., Krylov A.S., Orlov L.A., Alekseev V.A., Spil-rain E. E. Phase transition in liquid cesium near 590 К // J. of Non-Crystalline Solids. 2000. 277. N 2/3. P. 182-187.

21. Zhdanov M.S.,Dmitriev V. I., Fang Sheng, Hursan G. Quasi-analytical approximations and series in electromagnetic modeling // Geophysics. 2000. 65. N 6. P. 1746-1757.

22. Дмитриев В. И. О методах решения обратных задач // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2001. № 4. С. 3-7.

23. Дмитриев В. И., Мерщикова Н.А. Синтез магнитотеллурического поля // Физика Земли. 2002. № 11. С. 69-75.

24. БердичевскийМ. Н., Дмитриев В. И., Мерщикова Н.А. Об обратной задаче зондирования с использованием магнитотеллурических и магнитовариационных данных. М: МАКС Пресс, 2000.

25. Zhdanov M.S. Geophysical inverse theory and regularization problems. Elsevier, 2002.

26. Ильинский А. С., Галишникова Т.Н. Математическое моделирование процесса отражения плоской электромагнитной волны от волнистой поверхности // Радио и электроника. 1999. 44. № 7. С. 773-786.

27. Щедрин Б.М. О существовании решения обратных задач восстановления строения вещества по дифракционным данным // Прикладная математика и информатика. № 6. М: Издательский отдел ф-та ВМиК МГУ, 2000. С. 63-71.

28. Burova Е., Shchedrin В. A windows application: program for qualitative phase analysis of poly-crystalline mixtures // Cryst. Reports. 2000. 5. N 2. P. 340-342.

Поступила в редакцию 01.09.04

Е. В. Захаров, А. С. Ильинский

КОМПЬЮТЕРНЫЕ МЕТОДЫ АНАЛИЗА МОДЕЛЕЙ ПРИКЛАДНОЙ ЭЛЕКТРОДИНАМИКИ И НЕЛИНЕЙНОЙ ОПТИКИ

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

1. Дифракция электромагнитных волн на идеально проводящих бесконечно тонких незамкнутых поверхностях (экранах). Математическая постановка задачи состоит в следующем: определить в R3 электромагнитное поле {Е, Н}, гармонически зависящее от времени (exp( — iujt)), удовлетворяющее системе уравнений Максвелла

rot Н = —iuieE, rot Е = iufiH и краевому условию на идеально проводящей бесконечно тонкой поверхности S с краем (ребром) L:

Ё X п = —Ё° X п

s

где Е° — вектор электрической напряженности первичного поля, возбуждаемого в Д3 плоской электромагнитной волной или системой локальных дипольных источников. Кроме того, выполнены условия излучения Зоммерфельда на бесконечности и условия Майкснера [1] на ребре Ь.

Эта задача эквивалентна интегральному уравнению 1-го рода по поверхности 5, численное решение которого для плоского случая (¿'-поляризованное поле) было получено В. И. Дмитриевым и Е.В. Захаровым в конце 60-х гг. [2], а примененная вычислительная схема получила название "метода саморегуляризации". В последующие годы на основе этой технологии были решены многочисленные задачи дифракции и теории антенн [3, 4], которые сводились к одномерным интегральным уравнениям Фредгольма 1-го рода, ядра которых имели логарифмическую особенность.

В начале 80-х гг. в цикле работ А. Г. Давыдова, Е. В. Захарова, Ю. В. Пименова была разработана новая вычислительная технология решения одномерных и многомерных интегральных уравнений, названных гиперсингулярными, так как они содержали интегральные операторы, понимаемые в смысле конечной части по Адамару [5-8].

В трехмерном векторном случае поставленная электродинамическая задача эквивалентна следующему гиперсингулярному интегральному уравнению относительно вектора плотности электрических

токов ]{М) на поверхности 5 (Мо £ 5):

*

где значок (*) указывает, что интеграл следует понимать в смысле конечной части по Адамару, а Сг(М, Р) = ехр(гк11мр) — фундаментальное решение уравнения Гельмгольца в пространстве.

Вопросам разрешимости уравнения (1) посвящена монография [9], а вопросы сходимости соответствующих вычислительных схем (схема кусочно-постоянной аппроксимации и коллокации) рассматривалась в работах [9-12]. На основе численного решения уравнения (1) в последнее десятилетие было проведено моделирование широкого спектра моделей теории антенн [13-16]. Из последующих результатов отметим исследование по модификации методов и алгоритмов для экранов малых электрических размеров [17].

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

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

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

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

Для сведения исследуемых краевых задач к интегральным уравнениям использованы формулы Грина, двумерные периодические фундаментальные решения для неоднородных уравнений Гельмголь-ца в областях над и под периодической прозрачной границей д\)2 (М, Р), а также учитывались свойства поверхностных потенциалов. В результате для решения задачи отражения поля плоской трехмерной волны от прозрачной границы раздела сред получена система четырех интегральных уравнений относительно неизвестных компонент электрического и магнитного полей на границе раздела сред и их нормальных производных [23-25]. Если на периодической границе раздела сред выполняются условия Леонтовича, то вместо системы четырех интегральных уравнений получаем две независимые системы относительно неизвестных компонент искомого поля [25].

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

017

дп

вида

й(М) =

2тг

й(Р)—(д2(М1Р)-д1(М1Р)) +

1 + а дй(М) ~2 дп

+ Сдг

1

2п

дпг

й

ЖзР + йо(М), Мб 5о,

д2

дпмдпI

(</2(М, Р) - (М, Р))-

дй(Р) д

дпр дп

м

{ад2

й8р +

Мм)

дп

м е 50,

м

где а = 1, й(Р) = и(Р) = Е^\р) для Р-поляризации и а = е2/е\, й(Р) = = Н^'(Р) для ^-поляризации. Здесь 5о — период отражательной поверхности, а е\)2 — диэлектрические проницаемости сред над и под границей раздела двух сред.

Для численного решения система интегральных уравнений (2), (3) сводилась к решению систем линейных алгебраических уравнений. Поскольку ядра этих уравнений либо имеют логарифмическую особенность, либо требуют доопределения в точках совпадения аргументов и ряды плохо сходятся, были разработаны эффективные алгоритмы их вычисления. Разработанные алгоритмы были реализованы, и это позволило провести серию вычислительных экспериментов. В работах [26, 27] для этого случая исследуются характеристики отражения и проведен анализ распределения электромагнитного поля на границе раздела двух прозрачных сред. Рассмотрены случаи, когда отраженное поле имеет как одно-, так и многоволновой характер.

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

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

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

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

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

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

В настоящем обзоре рассматриваются методы построения вычислительных алгоритмов для задач дифракции на локально-неоднородном теле [29-31]. Будем рассматривать дифракцию заданного поля {Eq, созданного произвольной системой сторонних источников колебаний частоты / = и/2тг в

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

Для применения вычислительных методов важно свести краевую задачу в неограниченном пространстве к исследованию задачи в ограниченной области. Для такой формулировки задачи нам потребуется понятие нормальных сферических волн. Хорошо известно [31], что система однородных уравнений Максвелла в неограниченном пространстве с характеристиками £q, имеет счетную систему решений удовлетворяющих условиям излучения. Эта система волн может быть получена

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

Фп(М) = (I (k0r)Ykm (в, у), М=М(г, в, <р),

где СКХ) = л/ШН-1+1/2(х) — сферические бесселевы функции, ко = u^/So/Jo, Y™(6, tp) — сферические функции. Под индексом п здесь и далее будем понимать произвольно упорядоченную последовательность пар целых чисел (к,тп). Поля нормальных волн разделены на два класса решений — поля "магнитного" типа и поля "электрического" типа, которые выражаются через скалярную функцию ^„(М) следующим образом:

для полей "магнитного" типа

Ёп = Ё^ = rot(гФп1г), Нп = Н^ = - т-i— rot rot(гФп1г);

lOJ/I о

для полей "электрического" типа

Ёп = ЁW = rot rot{rWjr), Нп = НW = rot(r^„Tr),

к q ioj/I о

где ir — единичный вектор сферической системы координат.

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

Ё{г,в,<р) = ^ а

n,j

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

г 1) (Л 1

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

а"= Т / ^ ^ Лз = Д* / ^ ^ (1'?1

а условия излучения (4) сформулировать в виде импедансных условий на сфере 5д, объемлющей тело £):

= Е ¿г /[НК]!ГС1.3Р-ЁП{М) = [г(м,р)[н(р)Ъ(р)] ¿8р,

„ л Рп 3 Бп >

*

= — К{Р)Ёп{М).

Рп

п,3

п,3

Матрица Z(M, Р) является суммой диадных произведений векторов Е* <8> Еп, а /Зп — нормировочные коэффициенты. Таким образом, на любой сфере 5д, объемлющей тело £), соотношения (5) дают систему нелокальных линейных граничных условий, связывающих касательные составляющие векторов Ё и Н.

Внешняя задача дифракции эквивалентна краевой задаче для системы уравнений Максвелла внутри сферы 5д с условиями (5) на поверхности сферы.

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

Произвольный вектор, касательный к единичной сфере с центром в точке О, представим в виде линейной комбинации вектор-функций электрического и магнитного типа. Для удобства будем использовать две системы вектор-функций е„(#,</?) и ]гп(в,ср). Они совершенно равноправны на единичной сфере, а их обозначение связано с тем, что они пропорциональны касательным составляющим сферических волн на единичной сфере.

Приближенное решение краевой задачи дифракции в области г ^ Я определено в виде конечной суммы

Е? = Е <£(г)ёп(е}(р), н» = £

п п

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

У [хо1Ём - к*п%твйвй<р = О,

5г п= 1,...,ЛГ.

/ (кй Нм - ше(М)Ём^ е*п вт в <1в = / (/е„) (Ш,

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

Легко установить [31], что краевая задача для определения коэффициентов с^ и Ь^ однозначно

разрешима при любом И, а при А^ —> оо функции , Й14^ стремятся к решению исходной краевой

задачи в обобщенном смысле. Матрицы системы (6) полные, если е(М) и д(М) отличны от £q и Краевая задача для системы (6) при небольшом N может быть решена численно. Для решения таких краевых задач разработаны специальные экономичные методы решения [32], пригодные для систем дифференциальных уравнений высокой жесткости.

Однако порядок систем уравнений N определяется величиной безразмерного параметра k^R = = 2itR/\q. Как правило, N > k^R. Это связано со скоростью сходимости разложений по системе нормальных волн, которая подробно исследована для сферических неоднородностей, начиная с работы [1]. Разработка вычислительных методов для большого числа сферических гармоник приводит не только к увеличению объема вычислений, но и к накоплению ошибок.

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

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

ди ^ ( д2 и д2и\ т . . . ^ о

~dt+u~ \Щ + Щ)= а{и)} x =

Уравнение (7) описывает процесс релаксации и диффузии фазовой модуляции и = u(x,t]G) световой волны в тонком слое нелинейной среды, ограниченной апертурой Q. Наведенные неоднородности показателя преломления слоя считаются пропорциональными функции Iq(u), задаваемой суперпозицией отображения и —> 1(и) фазы в интенсивность, которая учитывает конкретный способ сложения входной волны и волны обратной связи, с преобразованием G = G(y) : R2 —> R2 пространственных аргументов. Управление преобразованием осуществляется с помощью призм, зеркал, а также современных управляющих устройств оптоэлектроники. Варьированием преобразования аргументов в эксперименте можно получить целый спектр явлений самоорганизации светового поля: многолепестковые и ротационные волны, оптические спирали, стационарные структуры. Богатая пространственно-временная динамика системы и ее чувствительность к выбору параметров преобразования аргументов обусловливают интерес к изучению соответствующих моделей, описываемых ФДУ с преобразованием аргументов, и постановке новых типов задач управления такими преобразованиями.

В случае гладкого обратимого преобразования G области Q на себя с единичным якобианом (например, отражения от координатных осей, поворота в круге с центром в начале координат) суперпозиция описывается формулой

IG

где G~l — обратное преобразование. В моделях с интерференцией волн (/(и) = К(1 + 7 cos и)) и поворотом аргументов в круге или кольце в [35, 36] дано описание вращающихся многолепестковых волн на основе теории бифуркации Андронова-Хопфа. В [45] установлено существование устойчивых стационарных структур. В [37] доказано существование глобального аттрактора и получены оценки его размерности. Аналогичные результаты для моделей с дифракцией, учитываемой дополнительным уравнением типа Шредингера, получены в [38, 41]. Модели оптических систем с запаздывающей обратной связью рассмотрены в [43].

С точки зрения задачи управления преобразованием аргументов важно максимально расширить класс допустимых преобразований и отказаться прежде всего от обременительного условия обратимости G(y). При этом формула (8), вообще говоря, теряет смысл даже для случая бесконечно дифференцируемых отображений (например, для G(y) = у0, где у0 — заданная точка из и необходимо обобщенное толкование (см. также [39]):

g~1

где С-1 (О) = {у : у Е Сг(у) Е — полный прообраз множества при отображении С (у). В [42] доказаны теоремы существования и единственности начально-краевой задачи Дирихле для (7), (9) в терминах пространств Соболева с дробными порядками гладкости, установлено существование аттрактора соответствующей динамической системы. Это позволило предложить постановки задач оптимального управления на классе произвольных измеримых по Лебегу преобразований аргументов. Адаптированная к (9) проекционно-разностная схема (ПРС) аппроксимации задачи (7) предложена в [44]; там же получена оценка скорости сходимости ПРС и соответствующих аппроксимаций задачи управления. Проведенное в [40, 46] численное моделирование задачи управления пространственным преобразованием показало эффективность данного метода управления в контуре обратной связи. Полученные оптимизированные преобразования, рассчитанные для актуальных с точки зрения приложений задач формирования локализованных структур, являлись именно необратимыми. Эти результаты подтверждают необходимость дальнейших исследований в области управления преобразованием аргументов.

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

Изучение процесса распространения светового импульса в одномерном фотонном кристалле проводится на основе следующего волнового уравнения:

d2E(z,t) n2(z) d2E(z,t) 4тг32Р„;

У ¿¿2 Рш = е

Здесь E(z,t) — напряженность электрического поля; z — координата, вдоль которой распространяется световой импульс; Lz — длина среды; t — время; n(z) — показатель преломления среды; с — скорость света; е(Е) — нелинейная добавка к диэлектрической проницаемости, a P„i описывает нелинейную поляризацию среды.

Волновое уравнение (10) с помощью традиционного подхода приводится к одному уравнению Шре-дингера для одноволнового случая [49-52] или к системе уравнений Шредингера для двухволнового взаимодействия фемтосекундных импульсов в фотонном кристалле [47, 48].

Для рассматриваемого уравнения построены инварианты и предложены консервативные разностные схемы, построенные на основе подхода, предложенного В. А. Трофимовым и основанного на отказе от выделения направления распространения световой волны. Компьютерное моделирование продемонстрировало преимущества этого подхода для данного класса задач [48, 49], также показана возможность локализации части световой энергии внутри одномерного нелинейного фотонного кристалла. В зависимости от длительности фемтосекундного импульса и его интенсивности в кристалле может быть локализовано до 50% прошедшей в него энергии. Данный эффект может быть использован для хранения информации в трехмерных оптических устройствах памяти [50, 51]. При анализе динамики взаимодействия субимпульсов, локализованных внутри некоторых слоев одномерного нелинейного фотонного кристалла, показано, что они движутся и взаимодействуют подобно солитонам. Описаны возможные сценарии их взаимодействия: пространственно-периодический режим; слияние нескольких близких по характеристикам субимпульсов в один более высокоинтенсивный и медленно распространяющийся субимпульс; формирование неподвижного импульса в слое фотонного кристалла [52]. Отметим, что исследования, результаты которых приведены в этом пункте, выполнены совместно с лабораторией математического моделирования в физике факультета ВМиК МГУ.

СПИСОК ЛИТЕРАТУРЫ

1. Хенл X., Мауэ А., Вестпфаль К. Теория дифракции. М.: Мир, 1964.

2. Дмитриев В.И., Захаров Е.В. О численном решении некоторых интегральных уравнений Фредгольма 1-го рода // Вычисл. методы и программирование. Вып. 10. М.: Изд-во МГУ, 1968. С. 49-54.

3. Захаров Е.В., Пименов Ю.В. Численный анализ дифракции радиоволн. М.: Радио и связь, 1982.

4. Захаров Е. В., Пименов Ю. В. Численный анализ дифракции электромагнитных волн на идеально проводящих экранах. Препринт ИРЭ АН СССР. № 23(395). М., 1984.

5. Давыдов А.Г., Захаров Е.В., Пименов Ю.В. Метод решения задач дифракции электромагнитных волн на бесконечно тонких цилиндрических экранах // ДАН СССР. 1981. 261. № 2. С. 338-341.

6. Давыдов А.Г., Захаров Е.В., Пименов Ю.В. Метод численного решения задач дифракции электромагнитных волн на незамкнутых поверхностях вращения // ДАН СССР. 1983. 269. № 2. С. 329-333.

7. Давыдов А.Г., Захаров Е.В., Пименов Ю.В. Метод численного решения задач дифракции электромагнитных волн на незамкнутых поверхностях произвольной формы // ДАН СССР. 1984. 276. № 1. С. 96-100.

8. Давыдов А.Г., Захаров Е.В., Пименов Ю.В. Метод численного решения задач дифракции электромагнитных волн на поверхностях с анизотропной проводимостью // ДАН СССР. 1985. 280. № 2. С. 333-337.

9. Ильинский А. С.,Смирнов Ю.Г. Дифракция электромагнитых волн на производящих тонких экранах: (Псевдодифференциальные операторы в задачах дифракции). М.: Радиотехника, 1996.

10. Захаров Е.В., Собянина И.В. Об одномерных интегродифференциальных уравнениях задач дифракции на экранах // ЖВМиМФ. 1986. 27. № 4. С. 632-636.

11. Давыдов А.Г., Захаров Е.В., Халеева И.В. Интегральные уравнения с ядрами типа Ада-мара в задачах дифракции // Актуальные вопросы прикладной математики. М.: Изд-во МГУ, 1989. С. 118-127.

12. Захаров Е.В., Халеева И. В. Гиперсингулярные интегральные уравнения 1-го рода задач дифракции электромагнитных волн на незамкнутых поверхностях // ЖВМиМФ. 1993. 33. № 2. С. 313-318.

13. Davidov A.G., Zakharov E.V., Pimenov Yu.V. Numerical analysis of fîelds in the case of electromagnetic excitation of unclosed surfaces // J. of Comm. Tech, and Elect. 2000. 45. Sup. 2. P. 247-259.

14. Давыдов A.Г., Захаров E.B., Пименов Ю.В. Гиперсингулярные интегральные уравнения в вычислительной электродинамике // Прикладная математика и информатика. № 9. М.: МАКС Пресс, 2002. С. 5-22.

15. Ilinskii A.S., Berezhnaya I. V., Alkhovski Е.А. Investigations of HF and VHF ground-based antennas for communication systems // Electromagnetics. 1999. 19. N 2. P. 171-185.

16. Ильинский A. С. Интегральные уравнения вибраторных антенн // Труды Математического центра им. Н. И. Лобачевского. Казань, 2002. 13. С. 5-13.

17. Давыдов А.Г., Захаров Е.В., Пименов Ю. В. О некоторых вычислительных аспектах применения метода гиперсингулярных интегральных уравнений в электродинамике // Прикладная математика и информатика. № 14. М.: МАКС Пресс, 2003. С. 11-17.

18. Дмитриев В.И., Захаров Е.В. Метод решения задач электродинамики неоднородных сред // ЖВМиМФ. 1970. 10. № 6. С. 1458-1464.

19. Захаров Е.В., Ильин И. В. Интегральные представления электромагнитных полей в неоднородной слоистой среде // Изв. АН СССР. Сер. Физика Земли. 1970. № 8. С. 62-71.

20. Захаров Е.В. Метод решения задач теории МТЗ для сред с горизонтально периодической структурой // Изв. АН СССР. Сер. Физика Земли. 1975. № 10. С. 53-56.

21. Захаров Е.В. О единственности и существовании решений интегральных уравнений электродинамики неоднородных сред // Вычисл. методы и программирование. Вып. 24. М.: Изд-во МГУ, 1975. С. 37-43.

22. Захаров Е.В. Математическое моделирование в электромагнитном каротаже. Л.: Недра, 1979.

23. Галишникова Т. Н., Ильинский А. С. Численные методы в задачах дифракции. М.: Изд-во МГУ, 1987.

24. Галишникова Т. Н., Ильинский А. С. Рассеяние плоской волны на волнистой поверхности / / Математические модели естествознания. М.: Изд-во МГУ, 1995. С. 86-111.

25. Ильинский А. С., Галишникова Т.Н. Математическое моделирование процесса отражения плоской электромагнитной волны от волнистой поверхности // Радиотехника и электроника. 1999. 44. № 7. С. 773-786.

26. Ильинский А.С., Галишникова Т.Н., Бережная И. В. Характеристики отражения от волнистой границы раздела прозрачных сред // Прикладная математика и информатика. № 3. М.: МАКС Пресс, 1999. С. 33-42.

27. Ильинский А.С., Галишникова Т.Н., Бережная И.В. Исследование характеристик отражения электромагнитного поля от волнистой периодической границы раздела прозрачных сред // Прикладная математика и информатика. № 6. М.: МАКС Пресс, 2000. С. 5-11.

28. Ильинский А. С., Галишникова Т.Н., Бережная И. В. Сравнение двух математических моделей в задаче дифракции iï-поляризованной волны на нерегулярной границе раздела сред // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2001. № 1. С. 8-13.

29. Свешников А.Г. Неполный метод Галеркина // ДАН СССР. 1977. 236. № 5. С. 1076-1079.

30. Свешников А.Г. Дифракция на ограниченном теле // ДАН СССР. 1969. 184. № 1. С. 63-65.

31. Свешников А. Г., Ильинский А. С. Прямой метод для задачи дифракции на локальном неоднородном теле // ЖВМиМФ. 1971. 11. № 4. С. 960-968.

32. Ильинский А. С., Быков А. А. Об одном методе решения краевых задач для систем линейных дифференциальных уравнений // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 1980. № 1. С. 11-19.

33. Апельцин В.Ф., Ильинский А. С., Сабитов Б. Р. Обоснование модифицированного неполного проекционного метода для задач рассеяния от гидрометеоров // ЖВМиМФ. 1986. 26. № 10. С. 1535-1551.

34. Ильинский А. С., Некрасов JI.M. Численный метод решения задач дифракции на неоднородном диэлектрическом цилиндре и его обоснование // ЖВМиМФ. 1995. 35. № 1. С. 53-70.

35. Разгулин А. В. Об автоколебаниях в нелинейной параболической задаче с преобразованным аргументом // ЖВМиМФ. 1993. 33. № 1. С. 69-80.

36. Разгулин А. В. Ротационные волны в оптической системе с двумерной обратной связью // Матем. моделир. 1993. 5. № 4. С. 105-119.

37. Vorontsov M.A.,Razgulin A.V. Properties of global attractor in nonlinear optical system having nonlocal interactions // Photonics and Optoelectronics. 1993. 1. N 2. P. 103-111.

38. Разгулин А. В. О конечномерном аттракторе нелинейной оптической системы // Изв. РАН. Сер. физич. 1996. 60. № 3. С. 127-135.

39. Потапов М.М. Уравнение нелинейной оптики с преобразованиями пространственной независимой переменной в роли управляющих воздействий // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 1997. № 3. С. 13-16.

40. Razgulin A.V. Localized and periodic patterns in nonlinear optical system with controlled transforms of spatial arguments // Proc. SPIE. 1998. 3733. P. 211-217.

41. Акимова И.Г., Разгулин А.В. Ротационные волны в оптической системе с дифракцией и поворотом пространственных аргументов // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 1999. № 2. С. 20-25.

42. Разгулин А. В. Об одном классе функционально-дифференциальных параболических уравнений нелинейной оптики // Диф. ур-ния. 2000. 36. № 3. С. 400-407.

43. Razgulin A.V. Finite-dimensional dynamics of distributed optical systems with delayed feedback // Computers and Mathematics with Applications. 2000. 40. N 12. P. 1405-1418.

44. Разгулин А. В. Аппроксимация задачи управления преобразованием аргументов в нелинейном параболическом уравнении // ЖВМиМФ. 2001. 41. № 12. С. 1844-1856.

45. Чушкин В. А., Р аз г у л и н А. В. Стационарные структуры в функционально-дифференциальном уравнении диффузии с отражением пространственного аргумента // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2003. № 2. С. 13-20.

46. Разгулин А. В., Саввина С. С. О численной оптимизации двумерного преобразования аргументов в функционально-дифференциальном уравнении диффузии // Прикладная математика и информатика. № 15. М.: МАКС Пресс, 2003. С. 26-35.

47. Mantsyzov B.I.,Fedotov M. V., Pospelova A. A. Nonlinear solitary waves in multidimensional resonant photonic bandgap structures // Proc. SPIE. 1999. 3736. P. 211-220.

48. Трофимов В.А., Терешии Е.Б., Федотов М.В. Консервативная разностная схема для задачи двухволнового взаимодействия фемтосекундных импульсов в фотонном кристалле // ЖВМиМФ. 2003. 43. № 10. С. 1550-1553.

49. Трофимов В. А., Терешин Е.Б., Федотов М.В. Исследование разностных схем для задач самовоздействия фемтосекундного импульса в фотонном кристалле // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2003. № 1. С. 20-27.

50. Трофимов В. А., Терешин Е. Б., Федотов М.В. О возможности локализации световой энергии в нелинейном фотонном кристалле // Оптика и спектроскопия. 2003. 95. № 1. С. 106-109.

51. Трофимов В. А., Терешин Е.Б., Федотов М.В. Локализация световой энергии последовательности фемтосекундных импульсов в одномерном фотонном кристалле // ЖТФ. 2004. 74. № 5. С. 66-70.

52. Трофимов В.А., Терешин Е.Б., Федотов М.В. Солитоноподобное распространение световых импульсов в нелинейном фотонном кристалле // Оптика и спектроскопия. 2004. 97. № 5. С. 823-832.

Поступила в редакцию 18.06.04

В. Б. Андреев, А. В. Гулин, Е. С. Николаев, Б. Н. Четверушкин ТЕОРИЯ ЧИСЛЕННЫХ МЕТОДОВ1

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

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

Ьи = — е(р(х)и')' — г(х)и' + д(ж)и = /(ж), 0 < х < 1,

и(0)=^о, и{1)=д1,

где р(х) ^ со > 0, г(х) ^ с\ > 0, д(ж) ^ 0, а е 6 (0,1] — малый параметр. Хорошо известно, что при е —т- 0 решение этой задачи сходится на полуинтервале 0 < х ^ 1 к решению вырожденной задачи

-г(х)у'+ д(х)у = ¡(х), у(1)=д1,

а при малых е не использованное в вырожденной задаче граничное условие приводит к образованию в окрестности точки х = 0 так называемого пограничного слоя, где решение и(х) исходной задачи (1) претерпевает сильные изменения. Наличие малого параметра и порождаемого им пограничного слоя приводит к значительным трудностям при численном решении задачи (1).

Обозначим одну из них на простейшем примере уравнения (1), когда р(х) = г(х) = 1, д(ж) = 0. Введем на отрезке [0,1] равномерную сетку с шагом /г = 1/^ и построим на ней классическую аппроксимацию уравнения (1) с указанными коэффициентами:

_= г = 1,..., N — 1.

1 Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (проекты 02-01-01089, 02-01-00555, 01-01-00488, 99-01-01078, 01-01-00572, 04-01-00633, 01-03-00606) и программы "Университеты России".