авторефераты диссертаций БЕСПЛАТНАЯ БИБЛИОТЕКА РОССИИ

КОНФЕРЕНЦИИ, КНИГИ, ПОСОБИЯ, НАУЧНЫЕ ИЗДАНИЯ

<< ГЛАВНАЯ
АГРОИНЖЕНЕРИЯ
АСТРОНОМИЯ
БЕЗОПАСНОСТЬ
БИОЛОГИЯ
ЗЕМЛЯ
ИНФОРМАТИКА
ИСКУССТВОВЕДЕНИЕ
ИСТОРИЯ
КУЛЬТУРОЛОГИЯ
МАШИНОСТРОЕНИЕ
МЕДИЦИНА
МЕТАЛЛУРГИЯ
МЕХАНИКА
ПЕДАГОГИКА
ПОЛИТИКА
ПРИБОРОСТРОЕНИЕ
ПРОДОВОЛЬСТВИЕ
ПСИХОЛОГИЯ
РАДИОТЕХНИКА
СЕЛЬСКОЕ ХОЗЯЙСТВО
СОЦИОЛОГИЯ
СТРОИТЕЛЬСТВО
ТЕХНИЧЕСКИЕ НАУКИ
ТРАНСПОРТ
ФАРМАЦЕВТИКА
ФИЗИКА
ФИЗИОЛОГИЯ
ФИЛОЛОГИЯ
ФИЛОСОФИЯ
ХИМИЯ
ЭКОНОМИКА
ЭЛЕКТРОТЕХНИКА
ЭНЕРГЕТИКА
ЮРИСПРУДЕНЦИЯ
ЯЗЫКОЗНАНИЕ
РАЗНОЕ
КОНТАКТЫ


Pages:     | 1 || 3 | 4 |   ...   | 6 |

«МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ИМЕНИ М.В. ЛОМОНОСОВА Факультет вычислительной математики и кибернетики ...»

-- [ Страница 2 ] --

Интересным представляется рассмотреть таблицу 3, содер жащую значения F1 -меры для различных значений параметров программы. По ней можно сказать, что наилучшее качество ре зультатов получается тогда, когда требуемое значение полноты Сборник статей молодых ученых факультета ВМК МГУ, № 10, 58 П. В. Алейников максимально, а требуемое значение точности – около 50%. Этот вывод представляется полезным для работы с программой.

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

Список литературы [1] Дроздов В. В. Метод автоматизированной генерации пра вил синтаксического анализа проектной документации.

М.: МИЭМ, 2010. 121 с.

[2] Тестелец Я. Г. Глава I. Слово и предложение. Структура зависимостей // Введение в общий синтаксис. М.: РГГУ, 2001. С. 61–106. 800 с. ISBN 5–7281–0343-X [3] Ian H. Witten and Eibe Frank Data Mining: Practical Machine Learning Tools and Techniques (Second Edition).

Morgan Faufmann, 2005 ISBN 0120884070.

[4] Синтаксически размеченный корпус русского язы ка: инструкция пользователя [HTML] (http:

//www.ruscorpora.ru/instruction-syntax.html) Сборник статей молодых ученых факультета ВМК МГУ, № 10, Генерация правил для синтаксического анализа [5] Розенталь Д. Э., Теленкова М. А. Словарь-справочник лингвистических терминов. Изд. 2-е. М.: Просвещение, 1976.

[6] Жеребило Т. В. Словарь лингвистических терминов: Изд.

5-е, испр. и доп. Назрань: ООО Пилигрим, 2010.

486с. ISBN 5989931336.

[7] F. Sebastiani Machine Learning in Automated Text Categorization [PDF] (http://nmis.isti.cnr.it/ sebastiani/Publications/ACMCS02.pdf) [8] Арефьев Н. В., Булгаков И. А., Мальковский М. Г. Оцен ка качества и сопровождение синтаксического анализа тора русскоязычных текстов Программные системы и инструменты: Тематический сборник. Под ред. Короле ва Л. Н. №12. М.: Издательский отдел Факультета ВМиК МГУ;

МАКС Пресс, 2011. С. 111–124.

[9] Christopher D. M., Prabhakar R. and Hinrich S. Introduction to Information Retrieval Cambridge University Press, 2008 ISBN 0521865719.

[10] Marshall, R J. Comparison of Misclassication Rates of Search Partition Analysis and other Classication Methods Statistics in Medicine, 2006 №25 p.3787–3797.

[11] Воронцов К. В. Комбинаторный подход к оценке качества обучаемых алгоритмов. Математические вопросы ки бернетики / Под ред. О. Б. Лупанов. М.: Физматлит, 2004. T. 13. С. 5–36.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 60 Сборник статей молодых ученых факультета ВМК МГУ, № 10, УДК 519. Задача управления портфелем пенсионного фонда при гарантированном исполнении обязательств в случае коррелированной динамики цен активов © 2013 г. Н. А. Андреев, А. В. Дружинина nandreev@cs.msu.su, alisa.druzhinina@gmail.com Кафедра системного анализа ВМК МГУ, Лаборатория по финансовой инженерии и риск-менеджменту НИУ ВШЭ 1 Введение В данной статье рассматривается задача нахождения оп тимальной инвестиционной стратегии пенсионного фонда, ра ботающего по схеме с фиксированными выплатами. Особенно стью инвестиционной политики подобного финансового инсти тута являются его обязательства перед участниками, а имен но –– обязательства произвести выплаты заранее установлен ного размера в заранее специфицированные моменты време ни. Следовательно, сумма активов фонда в момент выплаты должна гарантированно превышать сумму его обязательств.

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

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

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

Многопериодные модели появились в ответ на критику среднедисперсионного подхода [5], одним из основных недо статков которого являлась однопериодность модели, что про Управление портфелем пенсионного фонда тиворечит долгосрочной задаче инвестирования пенсионного фонда, который имеет возможность изменять инвестиционный портфель до момента осуществления выплат. В работах [2,7,8] задача инвестирования рассматривается уже как задача по иска оптимального управления, максимизирующего функцию ожидаемой полезности от капитала на конец горизонта пла нирования в предположении о независимости межвременного распределения доходности активов. Взяв в качестве функции полезности функцию с постоянной относительной несклонно стью к риску, авторы получили оптимальную инвестиционную стратегию, называемую в экономике близорукой : агент ведет себя таким же образом, как если бы он последовательно решал задачу инвестирования на один период, т.е. доли вложения в активы остаются неизменными на протяжении всего горизонта планирования. Классическая работа [6] рассматривает задачу оптимального управления в непрерывном времени. Результа том является утверждение, схожее с [8] для дискретного вре мени –– оптимальной является близорукая инвестиционная стратегия в следующих предположениях: используется изоэла стичная функция полезности, цены активов следуют процес су геометрического броуновского движения. В [3] рассматри вается более общий случай данных моделей: автор обобщает предыдущие результаты на класс логарифмических функций ожидаемой полезности с введением межвременных корреляций между доходностями активов, учетом стохастической динами ки доходов, а также неопределенности относительно смертно сти вкладчиков. Задача, поставленная подобным образом, ре шается аналитически при помощи динамического программи рования, однако при большем обобщении (например, вводе в модель торговых или регуляторных ограничений, транзакци онных издержек или при снятии предпосылки о нормальном распределении доходности активов) получение решения в яв ном виде невозможно. Применение численным методов также весьма ограниченно из-за т. н. проклятия размерности.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 62 Н. А. Андреев, А. В. Дружинина Другой распространенный подход, называемый правилом принятия решений, подразумевает априорное задание прави ла, в соответствии с которым стратегия может периодически перестраиваться. Например, в [4] этот подход используется для многопериодной модели управления активами и обязательства ми, где оптимизация проводится по доле активов с фиксиро ванной доходностью. К преимуществам подхода можно отне сти относительную простоту включения транзакционных из держек и дополнительных регуляторных, торговых или корпо ративных ограничений на инвестиционную политику, а также несложную интерпретируемость результатов для лиц, прини мающих решения. Недостатком же подхода является тот факт, что оптимизационная функция чаще всего не унимодальна и для нахождения решения требуется применять глобальные оп тимизационные алгоритмы, что налагает ограничения на ко личество управляющих переменных. Кроме того, правила при нятия решения могут быть по-разному специфицированы (со ответственно, давать различные результаты) и неочевидно, ка ким образом можно сравнивать разные спецификации, а также насколько это корректно.

Среди подходов к управлению активами и обязательствами широко используется стохастическое оптимальное управление, позволяющее включать в анализ коррелированность активов и обязательств, учитывать различные уровни неприятия риска, а также включать в себя торговые, регуляторные и корпоратив ные ограничения на инвестиционную политику. В данной ра боте будет рассмотрена одна из подобных постановок, заклю чающаяся в определении оптимальной стратегии управления портфелем из трех активов при наличии фазового ограниче ния в терминальный момент. Решение ищется в классе непре рывных адаптивных стратегий без учета торговых ограниче ний и транзакционных издержек. Целью работы является по лучение решения в полуаналитическом виде для общего случая совместной динамики цен и с учетом стохастической природы Сборник статей молодых ученых факультета ВМК МГУ, № 10, Управление портфелем пенсионного фонда обязательств.

2 Описание модели рынка и постановка задачи Рассматривается эффективный рынок, состоящий из безри скового актива, акций одного эмитента и безрисковых беску понных облигаций со всевозможными сроками до погашения из отрезка [0, T ]. Оптимальная стратегия подразумевает зна чение оптимальной доли средств, вкладываемых в безриско вый актив, акции и облигации с датой погашения T. Динами ка мгновенной процентной ставки rt, безрискового банковского счета Bt и цены акций St задается следующими стохастически ми дифференциальными уравнениями:

drt = µr (t, rt )dt + r (t, rt )dwt, r(0) 0, r (1) (2) dBt = rt Bt dt, B(0) = 1, dSt = µS (t, rt )St dt + S (t, rt )St dwt, S(0) 0. (3) r двумерный винеровский процесс, заданный на ве (wt, wt ) роятностном пространстве (, F, P) и порождающий фильтра цию Ft. Везде далее Et обозначает E(· | Ft ), E = E0. Предпо лагается, что коэффициенты уравнений удовлетворяют усло виям, гарантирующим существование сильного решения, при этом r, S 0. Для упрощения записи индексы и зависимо сти от переменных будут опускаться там, где их наличие не влияет на ход решения.

Предполагаем, что динамика rt дает аффинную временную структуру, т.е. цена облигации с датой погашения T имеет вид Pt = ea(t)b(t)rt для некоторых детерминированных функций a(t), b(t), определяемых уравнением (1). С помощью формулы Сборник статей молодых ученых факультета ВМК МГУ, № 10, 64 Н. А. Андреев, А. В. Дружинина Ито несложно получить, что динамика цены бескупонной об лигации Pt описывается уравнением dPt = µP (t, rt )Pt dt + P (t, rt )Pt dwt, P (0) 0.

r (4) Обозначим как Wt стоимость в момент t портфеля пен сионного фонда, состоящего из трех рассматриваемых акти вов. Стратегия управления характеризуется с помощью долей P S P S (1 t t, t, t ) = t от общей стоимости, вложенных в данные активы. Предполагается отсутствие торговых ограни чений, т.е. оптимальная стратегия ищется в классе A = {t, t [0, T ] : t Ft, Wt t квадратично интегрируем, 0 п.н.}. (5) WT LT При этом предполагается, что класс допустимых стратегий непуст. LT обозначает величину выплаты по обязательствам в терминальный момент T, которая в общем случае является случайной величиной. В таком случае считается, что случай ный процесс динамики популяции определен на собственном вероятностном пространстве (, F, P ) и порождает фильтра цию Ft, при этом LT m(FT ). Несложно видеть, что динамика стоимости портфеля описывается СДУ dBt dPt dSt dWt = (1 P S )Wt + P Wt + S Wt = Bt Pt St = (r + (µP r) P + (µS r) S )dt + P P dwt + S S dwt.

r (6) В качестве критерия оптимальности рассматривается мак симизация ожидаемой полезности от капитала, оставшегося по сле выплаты по обязательствам, для агента с постоянной отно сительной несклонностью к риску:

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Управление портфелем пенсионного фонда x (7) U (x) =, 0 1.

Задача оптимизации имеет следующий вид:

(8) EU (WT LT ) max.

A 3 Решение задачи Для случая двух активов (безрискового и акций) при отсут ствии ограничений задача была решена Мертоном [6], а позже рядом других авторов для случаев более сложной динамики цен и при наличии транзакционных издержек. Специфика за дачи управления для пенсионного фонда состоит в условии га рантированного исполнения обязательств перед участниками.

В рассматриваемом случае, гарантированное покрытие выплат общей суммой LT в конце периода. В данном разделе будет приведено решение задачи (8), основанное на подходе, предло женном в [1] для частного случая динамики (1) (4).

Рассмотрим вспомогательный процесс дисконтирования t t H P,H (u, ru ) dwu r Ht = exp µ (u, ru ) du 0 t S,H (u, ru ) dwu, удовлетворяющий СДУ dHt = µH (t, rt )Ht dt P,H (t, rt )Ht dwt S,H (t, rt )Ht dwt (9) r в предположении, что у уравнения существует сильное реше ние. Рассмотрим процесс Vt = Wt Gt, где Gt = Ht Et HT LT.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 66 Н. А. Андреев, А. В. Дружинина 1. Существуют µH, P,H, S,H такие, что для Теорема 1.

любого t A существует t A такой, что dV = (r +(µP r) P +(µS r) S )dt+ P P dwt + S S dwt.

r (10) 2. Пусть A = {t, t [0, T ] : t Ft, Vt t квадратично интегрируем, VT 0 п.н.}. Тогда max EU (WT LT ) = max EU (VT ).

A A Доказательство. 1) обозначим Vt = Ht Vt, Wt = Ht Wt, Gt = t = dWt dGt. Соглас = Ht Gt. Несложно убедиться, что dV но (6), dW = W dH + HdW + d[H, W ] = = W HµH dt P,H W Hdwt S,H W Hdwt + r +W H(r + (µP r) P + (µS r) S ) + P P W Hdwt + r + S S W Hdwt P,H P P W Hdt S,H S S W Hdt = = W H(r µH + (µP r) P P,H P P + (µS r) S S,H S S ))dt+W H( P P P,H )dwt +W H( S S S,H )dwt.

r Выберем µH, P,H, S,H следующим образом:

µH (t, rt ) = rt, (11) µP (t,r µS (t,rt )rt t )rt P,H (t, rt ) =, S,H (t, rt ) =, P (t,rt ) S (t,rt ) получая следующее уравнение для Wt :

dW = ( P P P,H )W dwt + ( S S S,H )W dwt.

r (12) 2) Gt = Ht Gt = Et HT LT M(Ft ). Согласно теореме о пред ставлении мартингала, существует единственный квадратично PS интегрируемый процесс t = (t, t ), такой, что Сборник статей молодых ученых факультета ВМК МГУ, № 10, Управление портфелем пенсионного фонда dG = P dwt + S dwt.

r (13) 3) Имеем два представления для Vt :

dV = dW dG = = ( P P P,H )W P dwt + ( S S S,H )W S dwt ;

r dV = V dH + HdV + d[H, V ] = = rV Hdt P,H V Hdwt S,H V Hdwt + HdV + d[H, V ].

r Приравнивая правые части, получаем dV = rV + ( P,H V + P P W P,H W H 1 P )dwt + r + ( S,H V + S S W S,H W H 1 S )dwt H 1 d[H, V ] (14) Таким образом, имеется два представления для Vt (10) и (14).

1 d[H, V ] :

Используем (10), (11) для вычисления H t H 1 d[H, V ] = H 1 ( P,H P P HV dt S,H S S HV dt) = = ( µP r) P + (µS r) S V dt. (15) Подставляя (10), (15) в (14), получаем ( P,H V + P P W P,H W H 1 P )dwt + r + ( S,H V + S S W S,H W H 1 S )dwt = = P P V dwt + S S V dwt, r откуда находим Сборник статей молодых ученых факультета ВМК МГУ, № 10, 68 Н. А. Андреев, А. В. Дружинина P P V = P,H (V W ) + P P W H 1 P, S S V = S,H (V W ) + S S W H 1 S, P P V = P P W P,H G H 1 P, S S V = S S W S,H G H 1 S, P P V +P,H G+H 1 P P =, P W (16) S S V +S,H G+H 1 S S =.

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

Доказанное утверждение позволяет решить сначала зада чу управления для Vt и найти t, оптимальную стратегию t исходной задачи можно вычислить потом с помощью отобра жения (16). Задача нахождения t известна и была решена в работе [6] для случая двух активов и постоянной мгновенной ставки. Случай трех активов и случайной ставки рассматри вается совершенно аналогично для функции цены I(t, rt, Vt ) = Et U (VT ). Уравнение Гамильтона–Якоби–Беллмана в этом слу чае имеет вид:

I I r I (r + (µP r) P + (µS r) S )V + max + µ+ t r V A 1 2 I r2 1 2 I P 2 P 2 2 ( + S S )V 2 = 0, (17) + + 2 2 r 2 V VT I(T, rT, VT ) =.

Замена I(t, rt, Vt ) = J(t, rt ) V обеспечивает наличие внут реннего максимума и позволяет записать (17) в виде Сборник статей молодых ученых факультета ВМК МГУ, № 10, Управление портфелем пенсионного фонда 2 µP r S r + µS P 1 2 J r J + J J r r + + = 0, r µ 2 r t 2(1) J(T, rT ) = 1.

(18) Уравнение (18) может быть решено численно с помощью фор мулы Фейнмана–Каца. Оптимальная стратегия находится ана литически из уравнения Гамильтона–Якоби–Беллмана:

I (µP 1 P µP r P = V I r) = JV (µ r) 2 =, 2 (1)P 2 J(1)V 2 V P P 2V V (19) S I = V (µS r) = JV 1 (µS r) µS r =.

2 I J(1)V 2 V S 2 (1)S S2 V V Подставляя (19) в (16), находим оптимальную стратегию исходной задачи с граничным условием WT LT :

µP (t,rt )rt Wt Gt P · 1 +Ht P = P (t,rt ), P t (t,rt )Wt µS (t,rt )rt Wt Gt (20) · 1 +Ht S S (t,rt ) S =, S t (t,rt )Wt P S P S t = (1 t t, t, t ).

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 70 Н. А. Андреев, А. В. Дружинина T T T µH (u,ru ) du P,H (u,ru ) dwu S,H (u,ru ) dwu r Et HT = Et e = 0 0 = Ht F (t, rt ) в силу марковости процесса мгновенной ставки. Поэтому G(t, rt, Ht ) = LT Ht F (t, rt ) = LT F (t, rt, Ht ). Предположим, что (t, r, H) непрерывно дифференцируема по t и дважды непре F рывно дифференцируема по r. Применяя формулу Ито к G, получаем F r F P,H F S,H r dG = LT H dwt + LT H dwt = r H H P = F r (t, rt ) µP (t,rt )rt Ht LT, t P (t,rt ) r (21) = S S = µ (t,rt )rt H L.

tT t S (t,rt ) Формула (21) позволяет моделировать процесс t на основе чис ленной оценки производной F.

r 4 Заключение В работе приводится решение задачи оптимального управ ления портфелем пенсионного фонда, т.е. агента, имеющего обязательства перед участниками. Отличительной особенно стью данной постановки является рассмотрение достаточно об щей стохастической совместной динамики цен активов. При этом учитывается такая важная характеристика как коррели рованность цен, которая вводится в модели посредством зави симости цен от общего фактора мгновенной ставки. Основ ными предположениями модели являются Сборник статей молодых ученых факультета ВМК МГУ, № 10, Управление портфелем пенсионного фонда 1. Одинаковый возраст участников на момент 0, что поз воляет говорить, что выплаты по обязательствам будут производиться в один и тот же момент для всех участни ков.

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

Заметим, что данное предположение не существенно для получения решения (20). Оно было использовано лишь для описания метода моделирования процесса t.

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

Полученное решение имеет интерпретируемый вид: доля вложений в рисковые активы (облигации и акции) растет с уве личением доходности актива относительно безрисковой (µP r и µS r соответственно) и уменьшается с увеличением вола тильности активов ( P S соответственно). С увеличением ве личины обязательств стратегия становится все более консерва тивной доля вложений в рисковые активы уменьшается. Вло жения в рисковые инструменты также снижаются при умень шении стоимости портфеля.

Предлагаемая постановка может быть обобщена с учетом смертности популяции и постоянного прихода новых участни ков. Основное отличие будет заключаться в введении дополни тельного семимартингального процесса накопленных дивиден дов. Вторым отличием будет введение последовательности вы плат по обязательствам. Оба процесса логичнее всего интегри ровать в существующую модель в качестве семимартингально го процесса потребления, в качестве еще одного управляюще го параметра. Нахождение оптимального потребления эквива Сборник статей молодых ученых факультета ВМК МГУ, № 10, 72 Н. А. Андреев, А. В. Дружинина лентно нахождению оптимальной величины взносов участни ков. Данная задача выходит за рамки работы и оставлена для дальнейшего исследования.

Список литературы [1] Deelstra G., Grasselli M., Koehl P. -F. Optimal investment strategies in the presence of a minimum guarantee // Insurance: Mathematics and Economics. 2003. №33. P.

189 207.

[2] Hakansson N. H. Optimal investment and consumption strategies under risk for a class of utility functions // Econometrica. –– 1970. –– №38. –– P. 587 607.

[3] Hakansson N. H. On optimal myopic portfolio policies, with and without serial correlation of yields // Journal of Business. –– 1971. –– №44. –– P. 324 334.

[4] Maranas C. D., Androulakis I. P., Floudas C. A., Berger A. J., Mulvey J. M. Solving long-term nancial planning problems via global optimization // Journal of Economic Dynamics and Control. –– 1997. –– №21. –– P.

1405 ––1425.

[5] Markovitz H. M. Portfolio selection // Journal of Finance. –– 1952. –– №7. –– P. 77 ––91.

[6] Merton R. Lifetime portfolio selection under uncertainty:

The continuous-time case // The Review of Economics and Statistics. 1969. №51(3). P. 247 257.

[7] Mossin J. Optimal multiperiod portfolio policies // Journal of Business. –– 1968. –– №41. –– P. 215 229.

[8] Samuelson P. A. Lifetime portfolio selection by dynamic stochastic programming // Review of Economics and Statistics. –– 1969. –– №51. –– P. 239 246.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Сборник статей молодых ученых факультета ВМК МГУ, № 10, 2013 УДК 519. Построение искусственных краевых условий для одномерной задачи распространения лазерного импульса в линейной среде с дисперсией третьего порядка © 2013 г. А. Д. Денисов antondnsv@gmail.com Лаборатория математического моделирования в физике 1 Введение В данной работе рассматривается проблема построения искусственных краевых условий для одномерного линейного уравнения Шрёдингера, учитывающего дисперсию третьего порядка. Точных неотражающих краевых условий для данного уравнения ранее предложено не было. Однако данная пробле ма является актуальной из-за широкого применения фемтосе кундных импульсов.

2 Постановка задачи Начально-краевая задача для безразмерного уравнения Шрёдингера с нулевыми краевыми условиями на левой гра нице и искусственным краевым условием на правой границе с начальным распределением в виде гауссова пучка записывает ся как:

2u 3u u u + + iD2 2 + D3 3 = 0, 0 z Lz, 0 t Lt, z t t t u(z0, t) = exp (t Lt )2, z0 (1) 0, [0, 1], 74 А. Д. Денисов u u u ut=0 = = 0, + aR + ibR u = 0, t z t t=0 t=Lt aR = + 2D2 R 3D3 2, bR = D2 2 2D3 3.

R R R Здесь u – медленно изменяющаяся амплитуда, R = (Lt ) – мгновенная частота светового импульса.

Известно аналитическое решение рассматриваемой задачи в случае нулевых краевых условий [1]:

U (z, t) = 2 1/2 |B|1/3 Ai 1 AB C 2 + i2C B 4/ 2 3AB 6C 2 3AB + 2C 2, (2) exp iC 3B 2 3B A = 2(t Lt z), B = 24D3 z, C = 4D2 z.

Здесь Ai - функция Эйри [2] от комплексного аргумента. Из-за специфики функций Эйри, решение вблизи z = 0 необходимо заменять асимптотическими формулами. В рамках данной ра боты ограничимся тем, что будем сравнивать аналитическое и численное решения начиная c z0 = 0.2. Хорошо известно [3], что функции Эйри осциллируют с уменьшающимся периодом и ин тенсивностью, но имеют бесконечную протяженность, поэтому их необходимо ограничивать при численном моделировании.

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

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Построение искусственных краевых условий I 1.0 z= z=0. z=0. 0. z= 0. 0. 0. 0. 20 25 30 35 40 t Рис. 1. Форма импульса в зависимости от z при = D2 = 0, D3 = 1.

3 Запись искусственного граничного усло вия Для записи искусственного краевого условия вводится но вая переменная u = ueit, которая учитывает распростране ние импульса направо. После этого получается уравнение от носительно u, в котором опускаются производные второго и третьего порядка. Далее выполняется обратная замена, откуда и получается искусственное краевое условие для правой гра ницы.

4 Преобразованное уравнение Из задачи (1) можно убрать слагаемые с D2, преобразовав её. Полученная система в некоторых случаях упрощает и уско ряет нахождение численного решения задачи.

Представляя u в виде u = ueit, частные производные пре Сборник статей молодых ученых факультета ВМК МГУ, № 10, 76 А. Д. Денисов образуются следующим образом:

u u u u eit, eit + i, u z z t t 2u 2u u iD2 2 iD2 eit 2 u, + 2i t t t 3u 3u 2u u D3 3 D3 eit + 3i 2 32 i3 u.

t t t t Откуда, группируя слагаемые, получаем:

2u u u 2D2 3D3 2 + i 2 (D2 + 3D3 ) + + z t t 3u + D3 3 + i D2 2 D3 3 = 0, u t где для обнуления второй производной считаем = D2 /3D3.

Для условия на правой границе получится условие, аналогич ное (1), но без слагаемых с D2. Окончательно получаем:

3u u u +a + D3 3 + ib = 0, u z t t D2 2D D a = + 2,b = 2, 3D3 3D3 27D u0 = u0 (t)eit, = D2 /3D3, (3) u u u ut=0 = = 0, + aR + ibR u = 0, t z t t=0 t=Lt aR = a 3D3 2, bR = 2D3 3.

R R Положив u = uibz для v легко найти аналитическое решение, т.к. получаем уравнение 3u u u (4) +a + D3 3 = 0, z t t решение которого хорошо известно (2).

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Построение искусственных краевых условий 5 Инварианты задачи Приведём инварианты задачи (1). В литературе они отсут ствуют.

Lt z u |u|2 + 2iD2 Im u I1 = |u| dt + + t 0 u0 2u u +D3 2Re u + D3 d, t2 t t Lt Lt 2 u u u u u I3 = iD2 D3 dt+ t2 t t t z u u u u + + 2iD2 Re + t Lt 2 u u u 2 u (5) +2iD3 Im D3 d.

t2 t t Lt 6 Разностная схема Введём равномерную сетку = z t :

z = {zj = jh, n = 0, 1,..., N, hN = Lz }, t = {tk = k, k = 0, 1,..., M, M = Lt }.

Определим на сеточную функцию A. Также будем использо вать следующие обозначения:

0. A = Aj,k = A(zj, tk ), A = A(zj+1, tk ), A = 0.5 A + A, 0.5 0.5 0. 0 A = 0.5( A k+1 A k1 )/, t Сборник статей молодых ученых факультета ВМК МГУ, № 10, 78 А. Д. Денисов 0.5 0.5 0.5 0. A k+1 2 A k + A k1 / 2, tt A = 0.5 0.5 0.5 0.5 0. 2 3, 0 A= 2 A k+1 + 2 A k1 A k A k+ t tt ttt A = (Ak+1 3Ak + 3Ak1 Ak2 ) 3.

В результате уравнение (1) во внутренних узлах сетки сво дится к следующему разностному виду:

AA 0. (6) + 0 + iD2 tt + D3 0 A = 0.

h t tt t Условие на правой границе записывается в следующем виде:

0.5 0.5 0. AM AM E A M + B A M 1 + C A M 2 0. + ibR A M = µ, (7) + aR h которое в случае (E, B, C) = (2, 2, 0) имеет порядок O(h2 + ) относительно точки (zj + h/2, tk ) при µ = µ0 = 0. Вто рой порядок достигается, если µ = µ1 = 0.5aR 2 u/t2 + + (1/6)aR 2 3 u/t3.

В случае (E, B, C) = (3, 4, 1) порядок равен O(h2 + 2 ), если µ = µ2 = (1/3)aR 2 3 u/t3.

Чтобы повысить порядок разностной аппроксимации гра ничного условия (7), необходимо аппроксимировать третью производную, используя уравнение для внутренних узлов сет ки: 3 u/t3 = (u/z + u/t + iD2 2 u/t2 )/D3. Запишем аппроксимацию µi :

0.5 0.5 0. D A M + G A M 1 + H A M 2 AM AM µ1 = a1 b1 + 2 h 0.5 0.5 0.5 0.5 0. A M 1 + G A M 1 + H A M DAM AM + + iD2, Сборник статей молодых ученых факультета ВМК МГУ, № 10, Построение искусственных краевых условий 0.5 0.5 0. AM AM 3 A M 4 A M 1 + A M µ 2 = b2 + + h 0.5 0. 0.5 0. D A M + G AM 1 A M 1 + H A M + iD2, a1 = 0.5aR, b1 = aR 2 /(6D3 ), b2 = aR 2 /(3D3 ).

Здесь важно отметить, что аппроксимация второй производной при (D, G, H) = (1, 2, 1) в точке k = M имеет первый поря док точности по. Но краевое условие (7) при этом учитывает вклад третьей производной. Схема порядка O(h2 + 2 ) для пра вого краевого условия (1) в данной работе не рассматривалась.

Тогда, учитывая последние выражения, запишем уравнение (7) в следующем виде (индекс i = 0 соответствует схеме (7) при µ = 0):

0.5 0.5 0. AM AM E A M + B A M 1 + C A M + qi + pi h 0.5 0.5 0. D A M + G A M 1 + H A M 2 0. +ri + ibR A M = 0, i {0, 1, 2} ;

p0 = 1, q0 = aR, r0 = 0, (E, B, C) {(2, 2, 0), (3, 4, 1)} ;

p1 = 1 + b1, q1 = aR + b1, r1 = a1 + iD2 b1, (E, B, C) = (2, 2, 0);

p2 = 1 b2, q2 = aR b2, r2 = iD2 b2, (E, B, C) = (3, 4, 1);

Запишем полученную задачу через систему разностных уравнений:

c0 A0 d0 A1 + e0 A2 = f0 ;

b1 A0 + c1 A1 d1 A2 + e1 A3 = f1 ;

aAk2 bAk1 + cAk dAk+1 + eAk+2 = fk, k = 2... M 2;

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 80 А. Д. Денисов aM 1 AM 3 bM 1 AM 2 + cM 1 ANt 1 dM 1 AM = fM 1 ;

(8) aM AM 2 bM AM 1 + cM AM = fM.

Учитывая граничные условия (1) и вводя обозначения: = 2 = 0.5D2 h/ 2, D3 = 0.5D3 h/ 3, легко вычислить 0.25h/, D коэффициенты полученной системы:

c0 = 1, d0 = e0 = f0 = 0, b1 = c1 = 1, d1 = e1 = f1 = 0, a = 0.5D3, b = D3 iD2 +, c = 1 i2D2, d = D3 iD2,e = 0.5D3, hz fk = A 0 + iD2 tt + D3 0 A, k = 2... M 2, 2 t tt t aM 1 = D3, bM 1 = 3D3 iD2 +, cM 1 = 1 i2D2 3D3, (9) dM 1 = D3 iD2, h fM 1 = A 0 + iD2 tt + D3 ttt A, 2 t aM = qi hC/(4 ) + ri hH/(2 ), bM = qi hB/(4 ) ri hG/(2 2 ), cM = pi + qi hE/(4 ) + ri hD/(2 2 ) + 0.5ihbR, i {0, 1, 2} ;

fM = pi AM qi h (EAM + BAM 1 + CAM 2 ) /(4 ) ri h (DAM + GAM 1 + HAM 2 ) /(2 2 ) 0.5ihbR AM.

Для решения системы (8), (9) будем использовать метод прямой пятиточечной прогонки [4], который может быть запи сан в виде Ak = F Ak1, Ak2. В работе [5] нами предложен консервативный разностный метод решения этой задачи с нуле выми граничными условиями. Также была проиллюстрирована связь между знаком дисперсии третьего порядка и разновид ностью применяемого в вычислениях метода прогонки.

Для корректности и устойчивости метода прогонки доста точно выполнение условий:

|c0 | |d0 | + |e0 |, |c1 | |b1 | + |d1 | + |e1 |, Сборник статей молодых ученых факультета ВМК МГУ, № 10, Построение искусственных краевых условий |c| |a| + |b| + |d| + |e|, c0, c1, c, cM 1, cM = 0, (10) |cM 1 | |aM 1 | + |bM 1 | + |dM 1 |, |cM | |aM | + |bM |.

Оценка в общем виде (i = 0, 1, 2) является трудоёмкой зада чей. По этой причине на практике используется непосредствен ная проверка выполнения этих условий, а также выбор ша гов, которые удовлетворяют условию h/ 3 1. Отметим, что это ограничение является очень обременительным с точки зре ния времени расчётов. В частности по этой причине построение искусственных краевых условий является актуальной задачей для рассматриваемой проблемы.

7 Вычисление мгновенной частоты волно вого пакета Выражение для частоты решения (2) в литературе отсут ствует, но её можно легко получить численно из определения комплексной функции u = uR + iuI. Так как arctan(uI /uR ) имеет разрывы, то необходимо осуществлять контроль этих то чек при вычислении разностной производной, учитывая это, нами был предложен способ вычисления мгновенной частоты импульса. Приведём её изменение вдоль z (Рисунок 2).

8 Результаты компьютерного моделирова ния Как видно из таблицы и рисунка 3, величина невязки = max |u|2 |A|2 по интенсивности комплексной амплитуды z,t уменьшается с уменьшением шагов сетки. Рассмотренные ва рианты для µi не дали ожидаемого повышения точности схемы.

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 82 А. Д. Денисов R z 0.2 0.4 0.6 0.8 1. Рис. 2. Зависимость мгновенной частоты аналитического реше ния (2) от z в точке t = Lt + 5 предполагаемой границы.

Таблица 1. Значения невязки в зависимости от шагов и µi при следующих значениях параметров: = 0, D2 = 1, D3 1, Lt = 195, Lt = 200, Lz = 1. Обозначения: a : µ0, (E, B, C) = (2, 2, 0), b :

µ0, (E, B, C) = (3, 4, 1) a b µ1 µ hz = 104, = 0.1 0.0473 0.0677 0.0673 0. · 105, 0.0391 0.0567 0.0565 0. hz = 5 = 0. Сборник статей молодых ученых факультета ВМК МГУ, № 10, Построение искусственных краевых условий z=0.4, =0.1 z=0.4, =0. 0.6 0. I I Analytical Analytical Numerical Numerical 0.4 0. 0.2 0. 0.0 0. 192 194 196 198 200 192 194 196 198 t t I I 0.6 0. z=0.6, =0.1 z=0.6, =0. Analytical Analytical Numerical Numerical 0.4 0. 0.2 0. 0.0 0. 192 194 196 198 200 192 194 196 198 t t 0.6 0. I I z=0.8, =0.1 z=0.8, =0. Analytical Analytical Numerical Numerical 0.4 0. 0.2 0. 0.0 0. 192 194 196 198 200 192 194 196 198 t t I I 0.6 0. Analytical Analytical z=1, =0.1 z=1, =0. Numerical Numerical 0.4 0. 0.2 0. 0.0 0. 192 194 196 198 200 192 194 196 198 t t Рис. 3. Распределение интенсивности в зависимости от z и величины шага по времени Сборник статей молодых ученых факультета ВМК МГУ, № 10, 84 А. Д. Денисов 9 Заключение В заключении отметим, что предложенное искусственное краевое условие может быть эффективно применено для рас сматриваемой задачи на одной из границ.

Настоящая работа выполнена при частичной финансовой поддержке гранта РФФИ № 12–01–91156_ГФЕН_а.

Список литературы [1] Miyagi M., Nishida S. Pulse spreading in a single-mode ber due to third-order dispersion. Appl. Opt., 18, p. 678, 1979.

[2] Abramowitz M., Stegun I. A. Handbook of Mathematical Functions. Dover, New York, 1970.

[3] Хонина С. Н., Волотовский С. Г. Зеркальные лазерные пуч ки Эйри. Компьютерная оптика, том 34, 203–213, 2010.

[4] Самарский А. А., Николаев Е. С. Методы решения сеточных уравнений. М.: Наука, 1978.

[5] Tromov V. A., Denisov A. D. Conservative Finite-dierence Scheme for the Problem of Laser Pulse Propagation in a Medium with Third-order Dispersion. Proceedings of IEEE East–West Design and Test Symposium (EWDTS’12). 225– 229, Kharkov, Ukraine, 2012.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Сборник статей молодых ученых факультета ВМК МГУ, № 10, 2013 УДК 519. Разностная схема для расчета 2D генерации полупроводниковой плазмы, индуцированной лазерным импульсом © 2013 г. В. А. Егоренков egorenkov-v-a@mail.ru Кафедра вычислительных методов 1 Введение Как известно, воздействие фемптосекундного лазерного импульса на среду, ввиду его высокой интенсивности, сопро вождается различными нелинейными эффектами [1], изучение которых имеет важное значение для различных прикладных задач. В частности, интерес представляет возможность постро ения оптически бистабильного элемента [2,3] на основе зависи мости коэффициента поглощения от концентрации свободных зарядов [4,5]. В ряде случаев в системе оптическое излучение – полупроводник развиваются сложные колебательные процес сы изменения характеристик полупроводника: концентраций свободных электронов и ионизированных доноров, потенциа ла электрического поля. Эти осцилляции характеристик сре ды подтверждены также аналитически для одномерной моде ли в рамках линейного анализа устойчивости решений соответ ствующей системы уравнений [6]. Для таких режимов взаимо действия возникает необходимость построения эффективных разностных схем. В [7] показано, что для одномерного случая известные в литературе разностные схемы в ряде случаев ока зываются не эффективными: итерационный процесс сходится слишком медленно или вообще не сходится. Ввиду этого по строение эффективных разностных схем для 2D задачи взаи 86 В. А. Егоренков модействия лазерного импульса с полупроводником является актуальной задачей.

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

2 + 2 = (n N ), x2 y n n n = Dx µn + Dy µn + t x x x y y y (1) + G (N, n, ) R (N, n), N = G (N, n, ) R (N, n), t I + 0 (N, n, )I = 0, 0 x Lx, 0 y Ly, t 0.

y Граничные и начальные условия задаются следующим обра зом:

= 0, 0 x Lx, y y=0,Ly = E0, 0 y Ly, x x=0,Lx n µn = 0, 0 x Lx, (2) Dy y y y=0,Ly n µn = 0, 0 y Ly, Dx x x x=0,Lx x Lx / 1 e10t, 0 x Lx, I|y=0 = exp 0.1Lx Сборник статей молодых ученых факультета ВМК МГУ, № 10, Разностная схема для 2D генерации плазмы. n|t=0 = N |t=0 = n0 eµ, 0 x Lx, (3) |t=0 = E0 x, 0 y Ly, I|t=0 = 0, что соответствует случаю нахождения полупроводника во внешнем электрическом поле.

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

N n n (4) G = q0 I (N, n, ), R =, p где (N, n, ) представляет собой нелинейный коэффициент по глощения среды.

В системе уравнений (1) – (4) введены следующие обозна чения: x – безразмерная поперечная координата, нормирован ная на радиус падающего оптического пучка, y – безразмерная продольная координата, нормированная на длину среды, Lx – безразмерная толщина поперечного сечения полупроводника, Ly – безразмерная длина полупроводника, t – время, измеряе мое в единицах времени рекомбинации свободных носителей за рядов, n(x, y, t) и N (x, y, t) – концентрации свободных электро нов в зоне проводимости полупроводника и ионизированных доноров, нормированные на максимально возможное в данных условиях их значение. Функция (x, y, t) – безразмерный по тенциал электрического поля, µ – коэффициент подвижности электронов, Dx, Dy - коэффициенты их диффузии по осям x и y соответственно. Параметр зависит, в частности, от мак симально возможной концентрации свободных носителей заря да, n0 – равновесное значение концентрации свободных элек тронов и ионизированных доноров полупроводника, p – время рекомбинации свободных носителей заряда. Заметим, что для выбранной нормировки времени этот параметр должен быть Сборник статей молодых ученых факультета ВМК МГУ, № 10, 88 В. А. Егоренков равен 1, однако он оставлен в (4) для удобства проведения рас четов. Функция I(x, y, t) описывает профиль интенсивности и временную форму оптического импульса, q0 – максимальное значение его интенсивности.

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

Коэффициент поглощения (N, n, ) аппроксимируем функцией (5) (N, n, ) = (1 N ) exp((1 n)), которая близка к одной из экспериментальных зависимо стей [8].

Для взаимодействия оптического излучения с полупровод ником справедлив закон сохранения заряда, который с учетом поставленных начальных условий записывается в виде (6) Q(t) = (n(x, y, t) N (x, y, t)) dxdy = 0.

3 Разностная схема Для построения разностной схемы введем в области = {0 x Lx } {0 y Ly } {0 t Lt } равномер G ные сетки = x y t, = x y t, = x y t :

x = {xi = ihx, i = 0, Nx, hx = Lx /Nx }, y = {yj = jhy, j = 0, Ny, hy = Ly /Ny }, y = {yj = (j 0.5)hy, j = 0, Ny + 1, hy = Ly /Ny }, t = {tk = k, k = 0, Nt, = Lt /Nt }, t = {t = (k + 0.5), k = 0, Nt 1, = Lt /Nt }.

k Определим на сеточные функции nh, Nh, h следующим об разом: nk = n(xi, yj, tk ), Nij = N (xi, yj, tk ), k = (xi, yj, tk ).

k ij ij Сборник статей молодых ученых факультета ВМК МГУ, № 10, Разностная схема для 2D генерации плазмы. Функцию nh определим на сдвинутой по времени сетке :

k = n(x, y, t ). Функцию I определим на сдвинутой по y nij ijk h : I k = I(x, y, t ). При этом будем использовать сле сетке ijk ij дующие обозначения:

k k k f = fij, fi±1 = fi±1j, fj±1 = fij±1, k k k k fi±0.5 = 0.5(fij + fi±1j ), fj±0.5 = 0.5(fij + fij±1 ), 0. k+ f = fij, f = 0.5(f + f ), 0. k k k I = Iij, Ii±1 = Ii±1j, I = Iij±1, I = 0.5(I + I), где f – одна из функций n, N,.

Для численного решения системы (1) – (5) построим кон сервативную разностную схему. Систему (1) аппроксимируем следующим образом [9]:

N N 0.5 0. =G R, nn = Dx nxx + Dy nyy + G R 0. i+1 i Dx µ ni+0.5 ni0.5 h2 h x x j+1 j Dy µ nj+0.5 nj0.5, 2 h hy y nn = Dx nxx + Dy nyy + G R (7) 0. i+1 i Dx µ ni+0. ni0. 2 h hx x j+1 j Dy µ nj+0. nj0., 2 h hy y xx + yy = ( N ), n I I 0.50. + 0 I = 0.

hy Сборник статей молодых ученых факультета ВМК МГУ, № 10, 90 В. А. Егоренков Построенная схема имеет второй порядок аппроксимации по пространству и времени на сетке =x y t.

В (7) дополнительно введены обозначения:

nN n2 nN n2 0.5 0 R=,R=, R = 0.5(R + R), p p 0.5 0.5 0. G = q0 I, G = q0 I, G = 0.5(G + G), 0. n = (1 N )e(1n), = (1 N )e(1 ), = 0.5( + ).

Граничные условия аппроксимируем с первым порядком по пространству, что обучловленно требованиям к консервативно сти:

y,i0 = 0, y,iNy = 0, i = 1, Nx 1, x,0j = E0, x,N j = E0, j = 1, Ny 1, x ( y,i0 ni0.5 y,i0 ) = 0, i = 1, Nx 1, n ( y,iNy niNy 0.5 y,iNy ) = 0, i = 1, Nx 1, n µ j = 1, Ny 1, (8) nx,0j (n0.5j x,0j + n0.5j x,0j ) = 0, µ nNx 0.5j x,Nx j + nx,Nx j = 0, j = 1, Ny 1, 2 + nN 0.5j x,N j x x i Nx / k (1 e(10k ) ), Ii0 = exp i = 1, Nx 1.

0.1Nx 4 Консервативность схемы Определим условия, при которых предложенная система разностных уравнений (7) удовлетворяет закону сохранения за ряда (6) и, соответственно, является консервативной. Для это го продифференцируем выражение для инварианта Q(t) (6), и Сборник статей молодых ученых факультета ВМК МГУ, № 10, Разностная схема для 2D генерации плазмы. аппроксимируем результат следующим образом:

Nx 1 Ny 1 nn N N (9) Q (t) = = 0.

hx hy i=1 j= Подставим в (9) первоее уравнение и сумму второго и третьего уравнений системы (7):

Nx 1 Ny (nt Nt ) = i=1 j= Dy Dy Dx nxx + 2 nyy + 2 nyy Dx µ (ni+0.5 x ni0.5 x ) 2hx Nx 1 Ny 1 Dy µ (nj+0.5 y nj0.5 y ). (10) = 2hy i=1 j=1 x ni0.5 x Dx µ ni+0. 2hx Dy µ y nj0.5 y nj+0. 2hy Вычислим суммы, входящие в правую часть равенства (10):

Nx nxx = ( x,Nx j nx,0j ), n hx i= Nx (ni+0.5 x ni0.5 x ) = n0.5j x,0j + nNx 0.5j x,Nx j.

i= Остальные суммы вычисляются аналогично.

В результате получим:

µ Ny 1 nx,0j (n0.5j x,0j + n0.5j x,0j ) + Dx + µ hx + nx,Nx j (nNx 0.5j x,Nx j + nNx 0.5j x,Nx j ) j=1 Сборник статей молодых ученых факультета ВМК МГУ, № 10, 92 В. А. Егоренков (ny,i0 µni0.5 y,i0 ) + ( y,i0 ni0.5 y,i0 ) + n Nx Dy.

+ (ny,iNy µniNy 0.5 y,iNy )+ 2hy + i= + ( y,iN niN 0.5 y,iN ) n y y y Теперь не трудно заметить, что краевые условия (8) при под становке в (10) обращают сумму в 0, а это значит, что схема (7) является консервативной.

5 Итерационный процесс Так как предложенная разностная схема (7) нелинейна, то для ее решения необходимо использовать итерационный про цесс:

s+1 s s N N 0.5 0. =G R, nn = Dx nxx + Dy nyy + G R 0. i+1 i Dx µ ni+0.5 ni0.5 h2 h x x j+1 j Dy µ nj+0.5 nj0.5, 2 h hy y s+ s s n n s+ (11) = Dx nxx + Dy n yy + G R 0. s s s s i+1 s i s Dx µ ni+0. ni0. h2 h x x s s s s j+1 s j s Dy µ nj+0. nj0., h2 h y y Сборник статей молодых ученых факультета ВМК МГУ, № 10, Разностная схема для 2D генерации плазмы. s+1 s+1 s+ s+ + yy = n N, xx s+1 s+1s+ I I 0.5 0. + 0 I = 0.

hy Для решения уравнения относительно использовался метод расщепления:

p+0.5 p p+0.5 p s+ s+ xx + yy = n N, p+1 p+0. p+0.5 p+ s+ s+ = + n N, xx yy где p - номер внутренних итераций.

Краевые условия (8) на итерациях принимают вид:

s+1 s+ y,i0 = 0, y,iNy = 0, i = 1, Nx 1, s+1 s+ x,0j = E0, x,Nx j = E0, j = 1, Ny 1, s+1 s+ n y,i0 = 0, n y,iNy = 0, i = 1, Nx 1, s µE nx,0j + (n0.5j + n0.5j ) = 0, j = 1, Ny 1, s µE + (nNx 0.5j + nNx 0.5j ) = 0, j = 1, Ny 1, nx,Nx j i Nx / k (1 e(10k ) ), = exp i = 1, Nx 1.

Ii 0.1Nx В (5) введены следующие обозначения:

ss s nN n2 s nN n2 0.5 s 0 R=, R = 0.5(R+ R),, R= p p Сборник статей молодых ученых факультета ВМК МГУ, № 10, 94 В. А. Егоренков s s 0.5s 0. s s 0. G = q0 I, G= q0 I, G = 0.5(G+ G), s s s s s 0. (1 n) (1n) = (1 N )e, = (1 N )e, = 0.5(+ ), s номер итерации.

Условия окончания итерацинного процесса:

s s | f f | 1 | f | + 2, где f – одна из функций n, N,, I.

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

При определенных параметрах системы в полупроводни ке устанавливался осциллирующий режим концентрации n. На рисунке 1 показано начало процесса, а на рисунке 2 устанав ливается периодичность. Расчет проводился с шагами hx = = 0.01, hy = 0.01, = 103 до момета времени T = 1500, при этом конечное занчение инварианта Q составило 4 104.

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Разностная схема для 2D генерации плазмы. a) b) c) d) Распределения концентраций свобод Рис. 1.

ных электронов, реализуемые при параметрах 0 = 2, Dx = Dy = 105, = 103, n0 = 0.01, µ = 1, = 2.553, = 3, q0 = 1.5 в моменты времени t = 80(a), 210(b), 320(c), 400(d).

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 96 В. А. Егоренков a) b) c) d) Распределения концентраций свобод Рис. 2.

ных электронов, реализуемые при параметрах 0 = 2, Dx = Dy = 105, = 103, n0 = 0.01, µ = 1, = 2.553, = 3, q0 = 1.5 в моменты времени t = 480(a), 500(b), 520(c), 850(d).

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Разностная схема для 2D генерации плазмы. Разработанная схема позволяет эффективно рассчитывать задачи взаимодействия ТГц импульса с полупроводником, ко торые все более актуальны в настоящее время [10, 11].

8 Благодарности Работа выполнена при частичной финансовой поддержке РФФИ(грант 12–01–31496 мол_а).

Список литературы [1] Делоне Н. Б., Крайнов В. П. Нелинейная ионизация атомов лазерным излучением. М.: Физматлит, 2001. 310 c.

[2] Гиббс Х. Оптическая бистабильность. Управление светом с помощью света. М.: Мир, 1988. 518 с, (Пер. с англ.

N. Y.: Academic Press, 1985).

[3] Розанов Н. Н. Оптическая бистабильность и гистерезис в распределенных нелинейных системах. М.: Наука, 1997. 334 с.

[4] Гасников А. И., Карамзин Ю. Н., Трофимов В. А. // Пись ма в ЖТФ. 1992. Т. 18, В. 6, с. 76–80.

[5] Бондаренко О. С., Трофимов В. А. // Письма в ЖТФ.

1996. Т. 22, В. 19, с. 6–9.

[6] Логинова М. М., Трофимов В. А. Осцилляции концентра ции свободных электронов при воздействии фемтосекунд ного светового импульса на полупроводник. // Оптика и спектроскопия. 2005. Т. 98, № 19, с. 1001–1007.

[7] Логинова М. М., Трофимов В. А.Разностная схема для за дачи воздействия фемтосекундного импульса на полупро водник при нелинейной подвижности. // ЖВМ и МФ.

2005. Т. 45, № 12, с. 2185–2196.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 98 В. А. Егоренков [8] Gibbs H. M., Olbright G. R., Peyghambarian N. et. al.

Kinks: Longitudinal excitation discontinuities in increasing absorption optical bistability. // Phys. Review. 1985.

Vol. 32, №. 1. P. 692–694.

[9] Марчук Г. И. Методы вычислительной математики. М.:

Наука, 1997. 456 с.

[10] Kersting R., Unterrainer K., StrasserG., Kaumann H. F., Gornik E. Few-cycle THz emission from cold plasma oscillations. // Phys. Review. 1997. Vol. 79 P. 3038– 3041.

[11] Arikawa T., Wang X., Belyanin A., Kono J. Giant tunable Faraday eect in a semiconductor magneto-plasma for broadband terahertz polarization optics. // Optics express.

2012. Vol. 20, №. 17. P. 19485–19492.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Сборник статей молодых ученых факультета ВМК МГУ, № 10, 2013 УДК 519.718. Некоторые задачи восстановления бесповторных функций © 2013 г. Д. В. Кафтан blond.programmist@gmail.com Кафедра математической кибернетики 1 Введение Проблема восстановления бесповторных функций состоит из задач, различающихся по трем основным параметрам: ба зис, цель тестирования [1] и используемый метод тестирова ния. Принципиальными являются отличия между четыремя типами базисов: базис {&, } (сокращенно MRO monotone read-once), базис B0 = {&,, ¬} (RO read-once), базисы, со стоящие из монотонных функций, повторных в базисе {&, } (MON monotone), и базисы, содержащие отрицание и функ цию, повторную в базисе B0 (BIG). Базисы считаются наслед ственными содержащими вместе с функцией все ее под фукции. Например, базис {&,, ¬} содержит также и 0, и 1.

Выделяют следующие цели: тестирование относительно беспо вторной альтернативы построение для бесповторной функ ции, существенно зависящей от всех своих переменных, про веряющего теста на множестве всех бесповторных функций тех же переменных (ALT read-once alternative [5]), класси ческая задача проверяющего тестирования на множестве всех бесповторных функций (TEST), диагностическое тестирование на множестве бесповторных функций, существенно зависящих от всех своих переменных (ESS essential, актуальность дан ной постановки обоснована в [6]), и классическая задача диа гностического тестирования на множестве всех бесповторных 100 Д.


В. Кафтан функций (DIAG diagnostic). Рассматриваемые диагностиче ские задачи являются условными (тест представляет собой ал горитм). Мы рассмотрим четыре оракула: точечные запросы (F (x)), оракул счетчик четности (PARITY), запросы тож дественности (IDENTITY) и оракул арифметическая сумма (SUM). Последние три оракула обращаются к произвольным подкубам. Под сложностью задачи подразумевается функция Шеннона [4] числа запросов. В силу того, что один из базисов содержит другой, справедливы следующие соотношения для задач, отличающихся только базисами: задачи с базисом MRO имеют не большую сложность, чем задачи с базисами RO и MON, а задачи с базисами RO и MON не большую, чем с базисом BIG. Аналогичные соотношения справедливы и для за дач, различающихся только целью тестирования или оракулом:

оракул F (x) имеет наименьшую мощность, поэтому сложность задач с его использованием будет наибольшей, а оракул SUM наибольшую мощность, и сложность задач с его использовани ем будет наименьшей;

задачи типа ALT имеют наименьшую сложность, а задачи типа DIAG наибольшую.

Один из основных вопросов проблемы восстановления бес повторных функций имеет задача полиномиальную или экспоненциальную сложность. Известно, что полиномиальную сложность имеют все задачи с базисом типа MRO [7], все за дачи с целью тестирования ALT [2] и все задачи с оракулом IDENTITY [2,3] (и, следовательно, с оракулом SUM). При этом из 64 = 4·4·4 задач остается неизвестной полиномиальность во семнадцати. О других задачах и связанных с ними результатах будет рассказано в заключении.

2 Результаты работы Рассмотрим задачу расшифровки булевой функции f, бес повторной в некотором монотонном базисе B и зависящей от переменных x1,..., xn (не обязательно существенно). Обозна Сборник статей молодых ученых факультета ВМК МГУ, № 10, Задачи восстановления бесповторных функций чим через T (n) функцию Шеннона числа точечных запросов, а через t(n) функцию Шеннона [4] сложности расшифровки f при помощи оракула IDENTITY. Известно, что t(n) = O(nl+2 ), где l наибольшая арность функции из базиса B [3].

Справедливо следующее утверждение:

Лемма 1. T (n) 2t(n).

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

Следствием леммы и результатов работы [3] является:

Теорема 1. T (n) = O(nl+2 ).

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

(1, 5)n.

Лемма 2. t(n) Доказательство. Как область, на которой конъюнкция обра щается в единицу, так и область запроса представляют собой подкубы. Пересечение подкубов есть подкуб. Ответ 1 на запрос возможен только лишь если они имеют единственную общую точку. Стандартное кодирование подкуба (грани) предусмат ривает подстановку соответствующей константы, либо прочер ка (в случае если обе константы возможны). Нетрудно убедить ся, что каков бы ни был подкуб, для каждой компоненты су ществуют два варианта координаты подкуба, пересекающегося Сборник статей молодых ученых факультета ВМК МГУ, № 10, 102 Д. В. Кафтан с ним в единственной точке (для прочерка константы, для константы она сама или прочерк). Рассмотрим пространство элементарных исходов – равновероятностное (1/3n ) простран ство элементарных конъюнкций. Для отличия их от нуля тре буется хотя бы один ответ 1. Событие Ai на i-й запрос получит ответ 1. Вероятность объединения событий не превосходит сум мы вероятностей: P (A1... At(n) ) P (A1 ) +... + P (At(n) ) = = 2n · t(n) · 1/3n. Так как P (A1... At(n) ) = 1, то t(n) (1, 5)n.

Следствием из леммы 2 является:

Теорема 2. Задача тестирования бесповторной функции в любом наследственном базисе, содержащем немонотонную и нелинейную функцию, оракулом счетчиком четности тре бует не менее (1, 5)n запросов в худшем случае.

3 Заключение Из оставшихся 18 задач экспоненциальную сложность, оче видно, имеют задачи с базисом RO и оракулом F (x). В рабо те [6] установлено, что задача ESS с оракулом PARITY для базиса RO полиномиальна, а для базисов BIG экспоненци альна. В данной работе показано, что задача TEST имеет экс поненциальную сложность для базиса RO и оракула PARITY (теорема 2). Также установлена полиномиальность всех задач с монотонными базисами (теорема 1). Таким образом, вопрос о полиномиальности всех 64 задач решен.

Список литературы [1] Чегис И. А., Яблонский С. В. Логические способы контро ля электрических схем // Тр. матем. ин-та им. В. А. Стек лова. 51. 1958. С. 270–360.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Задачи восстановления бесповторных функций [2] Chistikov, Dmitry V. Checking tests for read-once functions over arbitrary bases. // Lecture Notes in Computer Sci ence. —2012. —V. 7353. —P. 52–63.

[3] Чистиков Д. В. О связи задач диагностического и про веряющего тестирования бесповторных функций. // Дис кретная математика. 2011. Т. 23. № 1. C. 46–50.

[4] Shannon C. E. A symbolic analysis of relay and switching circuits // Trans. AIEE 57. —1938. —P. 713–723.

[5] Вороненко А. А. О проверяющих тестах для бесповтор ных функций // Математические вопросы кибернетики.

Физматлит. Вып. 11. 2002. С. 163–176.

[6] Вороненко А. А., Чистиков Д. В. Расшифровка бесповтор ных функций оракулом – счетчиком четности // Приклад ная математика и информатика. 2010. 34. С. 93–106.

[7] Гурвич В. А. О бесповторных булевых функциях // Успехи математическиз наук. 1977. 32. № 1. С. 183–184.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 104 Сборник статей молодых ученых факультета ВМК МГУ, № 10, УДК 517. Метод внутренней аппроксимации для поиска областей устойчивости аффинных полиномов © 2013 г. А. В. Мальцева amaltseva91@gmail.com Кафедра нелинейных динамических систем и процессов управления Введение Работа посвящена разработке алгоритма поиска областей устойчивости аффинных полиномов с использованием методов интервального анализа [3, 4]. Данная проблема возникает при решении одной из актуальных задач современной теории авто матического управления задачи об одновременной стабили зации линейных динамических объектов [2, 7].

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

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

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

1 Задача о поиске областей устойчивости для аффинных полиномов 1.1 Параметрические полиномы Параметрическим полиномом n-й степени на множестве Q Rk будем называть семейство полиномов, определяемое следующим образом:

a (s;

Q) = a0 (q) + a1 (q) s +... + an (q) sn, где ai (q) непрерывные функции на Q, q = (q1,..., qk ) Q вектор параметров, an (q) = 0 для всех q Q.

Аффинным полиномом n-й степени называют параметриче ский полином, для которого ai (q) = ai1 q1 +... + aik qk + ai0.

При этом, если ai0 = 0, то будем называть его однородным.

Интервальным полиномом n-й степени ([a] (s, Q) = [a0 ] + + [a1 ] s +... + [an ] sn ) называют параметрический полином, для которого 1. Q параллелотоп в Rk ;

2. ai (q) = qi, (q = (q1,..., qk )).

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 106 А. В. Мальцева Здесь через [ai ] обозначены так называемые интервальные чис ла, о которых пойдет речь в следующем параграфе.

1.2 Устойчивость параметрических полиномов Для параметрических полиномов можно ввести понятия устойчивости и неустойчивости, а именно, будем говорить, что параметрический полином a (s, Q) робастно устойчив, если для всех q Q полином a (s;

q ) = a0 (q ) + a1 (q ) s + +...+an (q ) sn устойчив (для любого фиксированного q поли ном a(s;

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

параметрический полином a (s, Q) робастно неустойчив, если для всех q Q полином a (s;

q ) = a0 (q ) + a1 (q ) s + +...+an (q ) sn неустойчив (т. е. один или несколько его корней лежат в правой полуплоскости комплексной плоскости или на мнимой оси);


параметрический полином a (s, Q) частично устойчив, если найдутся такие наборы параметров q, q Q, что полином a (s, q) устойчив, а a s, q неустойчив.

Далее, точку (q 1,..., q k ) Rk будем называть устойчи вой, если устойчив полином a (s;

q 1,..., q k ) и неустойчивой в противном случае.

Область G Q назовем областью устойчивости для по линома a (s;

Q), если все ее точки являются устойчивыми.

1.3 Общая постановка задачи Рассмотрим аффинный однородный полином третьей сте пени на Q Rk a (s;

Q) = an (q) sn +... + a1 (q) s + a0 (q) (1) Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации где q.

ai (q) = ai (q1,..., qk ) = (ai1,..., aik ).,.

qk q = (q1,..., qk ) Rk.

параллелотоп в Rk.

Пусть Q = q1, q1... qk, qk При решении различных задач теории управления (одно временная стабилизация, поиск минимального стабилизатора и т.д.) возникает следующая задача:

Существует ли область устойчивости G Q для задан ного аффинного полинома a (s;

Q)?

В общем случае сформулированная задача является нере шенной. Один из возможных подходов к решению изложен в работах Коровина С. К., Кудрицкого А. В., Фурсова А. С. [1,2], где для поиска областей устойчивости аффинных полиномов используются следующие методы и алгоритмы:

1. Методы интервального анализа;

2. Методы и алгоритмы исследования устойчивости пара метрических полиномов;

3. Алгоритм SIVIA (Set Invertor Via Interval Analysis) ал горитм обращения множеств (1).

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

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 108 А. В. Мальцева 2 Устойчивость параметрических полино мов 2.1 Некоторые понятия интервального анализа Основным объектом интервального анализа является так называемое вещественное интервальное число [3]. Веществен это односвязное подмножество ное интервальное число [x] из R, для простоты называемое интервалом. Множество ин тервальных чисел обозначают через IR.

Вещественным интервальным вектором [x] размерности n (или параллелотопом) называют подмножество из Rn, которое определяется как декартово произведение n замкнутых интер валов [x] = [x1 ]... [xn ]. Множество интервальных векторов размерности n обозначают через IRn.

Операции, введенные для интервалов [4], переносятся на па раллелотопы. Ширина непустого интервального вектора [x] определяется как w ([x]) max w ([xi ]).

1in это выпуклая комбинация образов вершинных Зонотоп точек параллелотопа при аффинном преобразовании.

Говорят [4], что интервальная функция [f ] : IRn IRm, действующая из пространства IRn в пространство IRn является функцией включения для функции f : Rn Rm, если для всех [x] IRn выполнено условие f ([x]) [f ] ([x]).

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации 2.2 Условия устойчивости параметрических поли номов.

Вначале рассмотрим проблему исследования на робастную устойчивость для интервального полинома степени 3 вида [a] (s;

Q) = [a0 ] + [a1 ] s + [a2 ] s2 + [a3 ] s3, (2) где [ai ] 0 для i = 0, 1, 2, 3.

Как известно, обычный полином 3-й степени с положитель ными коэффициентами u (s) = u0 + u1 s + u2 s2 + u3 s3, (3) ui 0, можно исследовать на устойчивость с помощью критерия Гур вица. Для этого рассмотрим функцию H (u) = u1 u2 u0 u3. То гда, в соответствии с критерием Гурвица, полином u (s) устой чив тогда и только тогда, когда H (u) 0.

Обобщим критерий Гурвица для исследования устойчиво сти интервального полинома (2). Для этого введем интерваль ную функцию [H] : IR4 IR [H] ([a]) = [a1 ] [a2 ] [a0 ] [a3 ], ([a] = ([a0 ], [a1 ], [a2 ], [a3 ])).

Фактически, H([a]) является функцией включения для функ ции (3).

Теорема 1. Интервальный полином (2) робастно устойчив тогда и только тогда, когда [H] ([a]) 0.

Теорема 1, по сути, представляет собой интервальный ана лог критерия Гурвица.

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

Пусть [H] ([a]) 0, т. е. [H] ([a]) (0, +), тогда [a1 ] [a2 ] [a0 ] [a3 ] (0, +), откуда следует (a1 a2 a0 a3 ) 0 для всех ai [ai ], i = 0, 3. Учитывая последнее неравенство получаем, Сборник статей молодых ученых факультета ВМК МГУ, № 10, 110 А. В. Мальцева что полином a0 +a1 s+a2 s2 +a3 s3 устойчив для всех ai [ai ], i = = 0, 3, а, следовательно, (2) робастно устойчив, т. к. он устой чив при любом наборе параметров, лежащих в указанных диа пазонах.

Докажем теперь необходимость. Пусть полином (2) робаст но устойчив, тогда полином a0 + a1 s + a2 s2 + a3 s3 устойчив для всех ai [ai ], i = 0, 3, то есть, a1 a2 a0 a3 0 для всех ai [ai ], i = 0, 3, из чего следует, что [a1 ] [a2 ] [a0 ] [a3 ] (0;

+), а, значит, [H] ([a]) (0, +), следовательно, [H] ([a]) 0.

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

Теорема 2. Интервальный полином робастно устойчив то гда и только тогда, когда [Hi ] ([a]) 0 для всех i.

Здесь Hi (a) i-й диагональный минор матрицы Гурви ца [6].

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

Чтобы воспользоваться теоремой 2 для получения условий робастной устойчивости параметрических полиномов рассмот рим преобразование произвольного параметрического полино ма a(s;

Q) = a0 (q) + a1 (q) s +... + an (q) sn, (4) где Q Rk параллелотоп, в интервальный вида [a] (s;

Q) = [a0 ] + [a1 ] s +... + [an ] sn, (5) которое будем называть V -преобразованием V : a (s;

Q) [a] (s;

P ) ([a] (s;

P ) = V {a (s;

Q)}).

По определению, при V -преобразовании коэффициенты поли нома (4) и интервального полинома (5) связаны следующим Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации образом a0 (q) [a0 ] a1 (q) [a1 ]..,..

..

an (q) [an ] где [ai ] = [ai ] ([q]), [q] = [q1 ]... [qk ] = Q, [ai ]([q]) функции включения для ai (q), P = [a0 ] [a1 ]... [a3 ], P IRn. Факти чески, функции включения [a]i ([q]) преобразуют параллелотоп Q IRk в интервалы [ai ] IR1 соответственно.

Теорема 3. Пусть полином [a] (s;

P ) = V {a (s;

Q)} робастно устойчив, тогда полином a (s;

Q) робастно устойчив.

Доказательство. Доказательство теоремы 3 следует из приве денной ниже цепочки импликаций: для всех q Q ai (q ) [ai ], i = 0, 1, 2,..., n a(s;

q ) устойчивый (поскольку ин тервальный полином [a] (s;

P ) робастно устойчив) a(s;

Q) робастно устойчив.

Рассмотрим подробно процедуру V -преобразования аф финного однородного полинома на примере полинома третьей степени на множестве параметров Q R a (s;

Q) = a0 (q) + a1 (q) s + a2 (q) s2 + a3 (q) s3, (6) q = (q1, q2, q3 ), Q = [q1 ] [q2 ] [q3 ], ai = ai (q1, q2, q3 ) = ai1 q1 + ai2 q2 + ai3 q в интервальный полином [a](s;

P ) = [a0 ] + [a1 ]s + [a2 ]s2 + [a3 ]s3, P = [a0 ] [a1 ] [a2 ] [a3 ] R4.

Линейное преобразование A : R3 R4, задаваемое равен ствами ai = ai (q1, q2, q3 ) = ai1 q1 + ai2 q2 + ai3 q3 q Q, i = 0, преобразует трехмерный параллелотоп Q в зонотоп Q, ле жащий в пространстве коэффициентов аффинного полинома.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 112 А. В. Мальцева Фактически, каждая точка зонотопа Q соответствует (как на бор коэффициентов) некоторому полиному семейства, задан ного аффинным полиномом.

Далее, преобразование V : a (s;

Q) [a] (s;

P ) с функциями включения (7) [ai ] = ai1 [q1 ] + ai2 [q2 ] + ai3 [q3 ] можно интерпретировать как преобразование в пространстве коэффициентов полиномов третьей степени, а именно его мож но понимать как преобразование параллелотопа Q в прямой параллелотоп P (при этом каждая точка параллелотопа P со ответствует некоторому полиному семейства, заданного интер вальным полиномом [a](s;

Q)), либо можно понимать как пре образование множества параметров Q во множество коэффи циентов P интервального полинома.

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

Перейдем теперь к вопросу о робастной неустойчивости па раметрических полиномов. Но прежде, сформулируем доста точное условие робастной неустойчивости интервальных поли номов.

Теорема 4. Пусть для полинома (5) найдется индекс i, для которого [Hi ] ([a]) 0, тогда [a] (s;

P ) робастно неустойчив.

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

Используя тот факт, что Q P можно сформулировать теорему о достаточном условии робастной неустойчивости па раметрических полиномов.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации Теорема 5. Пусть интервальный полином [a] (s;

P ) = = V {a (s;

Q)} робастно неустойчив, тогда параметрический полином a (s;

Q) робастно неустойчив.

Заметим, что критерий робастной устойчивости интерваль ных полиномов дает также известная теорема Харитонова [5], однако достаточное условие робастной неустойчивости (а его роль существенна в нахождении областей устойчивости аф финных полиномов) интервальных полиномов возможно сфор мулировать на основе критерия Гурвица (теорема 4), в то вре мя как критерий Харитонова не позволяет это сделать.

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

Введем следующие обозначения:

область устойчивости в пространстве параметров SQ Rk ;

область неустойчивости в пространстве парамет N Q ров Rk ;

область неопределенности в пространстве пара UQ метров Rk.

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

Основные шаги метода внешней аппроксимации:

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 114 А. В. Мальцева 1. К аффинному полиному a(s;

Q) = a0 (q) + a1 (q) s +... + + an (q) sn, q = (q1,..., qk ) Rk, Q = q1, q1... qk, qk Rk параллелотоп, применяется линейное преобразование A : Rk Rn+1, задаваемое равенствами ai = ai (q1,..., qk ) = ai1 q1 +... + aik qk q Q, i = 0, n, которое преобразует параллелотоп Q в зонотоп Q.

Q = {q : q = Aq, q Q} 2. Строится интервальный полином [a] (s;

P ) такой, что Q P. Для этого аффинный полином переводится в ин тервальный (преобразование параллелотопа Q в паралле лотоп P ) с помощью V -преобразования:

V : a (s;

Q) [a] (s;

P ), [ai ] = [ai ] (q1,..., qk ) = ai1 [q1 ] +... + aik [qk ].

В зависимости от устойчивости [a] (s;

P ) в соответствии с теоремами 2 – 5 делаются выводы о робастной устойчи вости a (s;

Q).

3. Если [a] (s;

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

вход: a (s;

Q), [H] ([a]), Q, ;

выход: S, N, U.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации Таблица 1. Алгоритм SIVIA.

1. Если [H] ([a]) 0 S = Q, N =, U = ;

end;

2. Если [H] ([a]) 0 N = Q, S =, U = ;

end;

3. Если w (Q) U = Q;

end;

4. SIV IA([H] ([a]), QL, );

SIV IA([H] ([a]), QR, );

где QL и QR получаются путем деления пополам наибольшего ребра параллелотопа Q.

Здесь точность работы алгоритма (фактически, деле ние параллелотопов происходит до тех пор, пока их ши рина не станет меньше ).

Результатом работы алгоритма явялется разбиение исход ного множества параметров Q на три подмножества Q = Sn Nn Un, где Sn, Nn, Un множества устойчивости, неустойчивости и неопределенности соответственно, сформированные после n-го шага работы алгоритма SIVIA.

3.1 Задача об областях неопределенности для аф финных полиномов.

Пусть Un область неопределенности после n-го шага ал горитма. µ (Un ) представляет собой сумму объемов неопре деленных параллелотопов (мера множества), 0 µ (Un ) V (Q).

Рассмотрим функцию (n) = µ (Un ). Будем говорить, что метод внешней аппроксимации сходится, если (n) 0. В связи с этим возникает задача об исследовании поведения функции (n) в процессе работы алгоритма SIVIA.

Можно привести пример классов аффинных полиномов, для которых (n) = µ (Q), т. е. (n) 0.

Сборник статей молодых ученых факультета ВМК МГУ, № 10, 116 А. В. Мальцева Для простоты рассмотрим аффинный полином 3-ей степени вида (6). Пусть [a] (s;

P ) = V {a (s;

Q)}. Тогда для интерваль ного полинома [a] (s;

P ) = [a0 ] + [a1 ] s + [a2 ] s2 + [a3 ] s3 матрица Гурвица имеет вид [a1 ] [a0 ] [a3 ] [a2 ] [a1 ] 0 0 [a3 ] с главными минорами [M1 ] = [a1 ];

[M2 ] = [a1 ] [a2 ] [a0 ] [a3 ];

[M3 ] = [M2 ] [a3 ].

Тогда при [ai ] = [a], i = 0, 3 или [a1 ] [a2 ] = [a0 ] [a3 ], учи тывая, что [a] [a] = [a a, a a], получим M1 0 и M2, M3 с границами, симметричными относительно нуля.

Так как ai вычисляются по формуле (7), то при использо вании алгоритма SIVIA (n) 0 и (n) = µ (Q) для всех n.

Пусть теперь для аффинного полинома a (s;

Q) выполняет ся равенство a1 (q) a2 (q) = a0 (q) a3 (q). Тогда при использова нии алгоритма SIVIA (n) 0 и (n) = µ (Q) для всех n.

Пример численного моделирования метода внешней ап проксимации [1, 2]:

Зададим точность алгоритма = 0.1. Рассмотрим поли ном a (s;

Q) = a3 (q1, q2, q3 ) s3 + a2 (q1, q2, q3 ) s2 + a1 (q1, q2, q3 ) s + + a0 (q1, q2, q3 ), q1 = [1, 2], q2 = [1, 3], q3 = [1, 2], a0 = q1 + 3q2, a1 = 2q1 + 4q2, a2 = 2q1 + q2 + q3, a3 = q1 + 3q2 + q3.

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

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации Graphic 2. 2 1. 1. 1. 1. measures 1. 1. 1. 1. 2.5 1. 2 0. 1. 1. 1. 1.2 2 1.8 1.6 1.4 1.2 1 0.8 0.6 0.4 0. 1 epsilons Рис. 1. Результат работы программы: разбиение исходного парал лелотопа на три множества (от темного к светлому N, U, S) и график изменения меры множества неопределенности.

1. Вскрытие механизма стремления (а также его отсут ствия) к нулю функции (n).

2. Разработка алгоритма поиска областей устойчивости с га рантированной сходимостью ((n) 0).

3.2 Причины отсутствия сходимости метода внешней аппроксимации.

Введем несколько определений.

Параллелотоп [q1 ]... [qn+1 ] Rn+1 называем вырож денным, если существует i {1,..., n + 1} такое, что qi = qi.

Порядок вырождения это количество i {1,..., n + 1} та ких, что qi = qi.

Интервальный полином [a](s;

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

При максимальном порядке вырождения n+1 паралле лотоп превращается в точку, а соответствующий интервальный полином – в обычный полином.

Рассмотрим пространство обычных полиномов q (s) = q1 + q2 s +... + qn+1 sn (q1,..., qn+1 ). Его можно представить в Сборник статей молодых ученых факультета ВМК МГУ, № 10, 118 А. В. Мальцева виде объединения двух множеств:

множество устойчивых полиномов;

•S •N множество неустойчивых полиномов.

Обозначим через S границу между этими множествами.

Для невырожденных интервальных полиномов верна сле дующая теорема.

Теорема 6. Для любого невырожденного интервального полинома (n) 0 при n.

Доказательство. Пусть [a] (s, Q) Q невырожденный па раллелотоп, тогда возможны следующие варианты:

1. Q S [H] (a) 0 Q множество робастной устой чивости;

2. Q N [H] (a) 0 Q множество робастной неустойчивости;

3. Q S = Q область неопределенности.

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

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

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

В связи с этим возникает вопрос: как быть в случае, если (n) 0?

Сборник статей молодых ученых факультета ВМК МГУ, № 10, Метод внутренней аппроксимации 4 Метод внутренней аппроксимации Для построения алгоритма, гарантирующего сходимость, предлагается следующая идея:

1. Переводим параллелотоп Q в зонотоп Q (Q = A (Q)).

2. Приближаем Q изнутри некоторым множеством P = параллелотопы в Rn+1 ).

= Pi (Pi 3. Исследуем на устойчивость P (описано далее).

Возможны два случая:

1. Q невырожденное множество в Rn+1 (содержит внутри себя некоторый невырожденный параллелотоп).

вырожденное множество в Rn+1.

2. Q Рассмотрим задачу: как для невырожденного Q построить P Q (фактически это вопрос о построении внутренней ап проксимации множества)?

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

1. Задаем 0.

2. Строим внутреннюю аппроксимацию множества Q. Q = (Pi ) B : µ (B ).

3. Пусть = min {µ (Pi )}. Задаем 0 такое, что ( i точность работы алгоритма SIVIA).

4. SIV IA = SIV IA (Pi ) (применение алгоритма SIVIA к i каждому параллелотопу Pi ).

Тогда, в результате работы SIVIA, общая мера множества неопределенности µ (U ) = µ (B ) + i i.

Поскольку, в силу теоремы 6 i (n) 0, n, мера множества неопределенности µ (U ) может быть сделана сколь угодно малой.



Pages:     | 1 || 3 | 4 |   ...   | 6 |
 





 
© 2013 www.libed.ru - «Бесплатная библиотека научно-практических конференций»

Материалы этого сайта размещены для ознакомления, все права принадлежат их авторам.
Если Вы не согласны с тем, что Ваш материал размещён на этом сайте, пожалуйста, напишите нам, мы в течении 1-2 рабочих дней удалим его.