CC BY
10УДК 621.372.542
Некорректные задачи и многокритериальное программирование
© А. А. Грешилов МГТУ им. Н.Э. Баумана, Москва, 105005, Россия
Рассмотрено решение некорректных задач методами многокритериального математического программирования, позволяющими избежать введения параметров регуляризации. Используются одновременно метод сжатия области допустимых значений и целевое программирование, позволяющие учесть неотрицательность и ограниченность решения. Метод показан на примере определения параметров ядерного взрыва по изотопам криптона и ксенона. При регистрации малого числа изотопов применяют объединение двух видов мгновенного деления урана 235 и плутония 239 в один вид деления. Одновременно рассматривают варианты механизма ядерного взрыва. Находят точечные оценки вкладов различных видов деления в суммарную активность изотопов. Для определения момента сепарации вклады рЫу рассчитывают при разных значениях и выбирают значение , при котором отношение Г \2
п т т
^ Д ^) - ^ а*ше (¿ц, IХр^у ) / ^ (Р^у- )2 минимально. Для этой цели сформи-
i=l I У =1 ) У=1
рован функционал Е, из которого искомые вклады получают путем дифференцирования его по элементам системы линейных алгебраических уравнений (СЛАУ) при фиксированных значениях удельной активности на каждой итерации.
Ключевые слова: методы регуляризации, продукты деления, ядерный взрыв, сепарация изотопов, многокритериальное программирование, целевое программирование, метод сжатия, итерационный метод решения.
В процессе решения задач различной природы часто имеет место ситуация, когда небольшие изменения исходных данных приводят к значительным изменениям в решении. Такие задачи называются некорректными [1-5]. Приведем строгое определение некорректных задач.
Пусть х и у — некоторые соответственно искомые и наблюдаемые характеристики модели исследуемого процесса; х является элементом метрического пространства и, а у — элементом метрического пространства Е. Задан оператор А, действующий из и в Е и устанавливающий причинные связи между искомыми характеристиками х модели и входными (наблюдаемыми) данными у. Область определения оператора БА с и, область значений (А = А(БА) с Е. Тогда можно записать операторное уравнение
Ax = y. (1)
Задача (1) считается поставленной корректно, если ее решение удовлетворяет следующим требованиям (условиям Адамара):
• решение существует для любого y е Qa с F (условие разрешимости уравнения);
• решение единственно в U (условие однозначности);
• решение непрерывно зависит от y, т. е. если приращения Ay стремятся к нулю, то изменения решения Ax также стремятся к нулю (условие устойчивости).
Иначе говоря, задача решения уравнения (1) корректна, если существует однозначное и непрерывное отображение A"1, область определения которого DA1 = QA = F.
Если нарушается хотя бы одно из перечисленных требований, задача решения уравнения (1) становится некорректно поставленной.
В общем случае решение Х1 = Aне обладает свойством устойчивости к малым изменениям правой части y1(х). В качестве меры некорректности задачи выступают числа обусловленности матрицы системы.
Родоначальником методов решения некорректных задач является А.Н. Тихонов. Первые некорректные задачи решены в 1940-х годах. Большое развитие они получили после 1960-х годов.
Задача нахождения приближенного решения уравнения (1), устойчивого к малым изменениям его правой части, сводится:
• к нахождению регуляризирующих операторов;
• определению параметра регуляризации а [1] по дополнительной информации о задаче, например по значению погрешности 8, с которой задается правая часть уравнения y8 .
Регуляризирующие алгоритмы могут быть получены, если использовать идею стабилизации минимума уклонения х. Функционал Q(x) , определенный на непустом множестве UQ с U, так называемый стабилизатор, вводится, если:
• Q[ х] > 0 для всех x е Uq ;
• множество Qc ={х: хeUq, Q[х]< Cj является р -компактным при любом C = const > 0, т. е. из любой последовательности { xk j е Qc можно выбрать подпоследовательность {xku j, р -сходящуюся к
некоторой точке x е Qc ;
*
• множество Uq = Uqo U*, где U* — множество точек мини-
2
мума функции J(x) = pF (Ax, y) , не пусто.
Далее берется какая-либо положительная последовательность {а^ ], сходящаяся к нулю, например а к = а 04k, |q| < 1, и при каждом к = 1, 2,... на множестве Uq определяется функционал Тихонова [1]
Tk (x) = J (x) + а к Q[ x], x eUq .
Здесь J(x) — некий функционал, определяющий соответствие наблюдаемых значений и решения.
Минимум функционала Тихонова для различных значений к определяет минимизирующую последовательность {xk ] , сходящуюся к
регуляризованному решению.
Затем появился другой подход к решению некорректных задач, названный методом статистической регуляризации. В работе [5] впервые было введено понятие функции плотности вероятности для оценок вектора решения. В данном методе предполагается, что известны априорные функции плотности вероятности исходных данных и параметра регуляризации. Считают, что элементы вектора у не коррели-рованы, измеряются со среднеквадратическим отклонением (СКО), равным а (СКО шума), и распределены по нормальному закону. По формуле Байеса получают апостериорное распределение для x.
Целевая функция метода регуляризации посредством минимизации максимальной энтропии [5] отличается от целевой функции метода А.Н. Тихонова лишь стабилизатором и имеет следующий вид:
n
J(x, а) = 11 Ax - y|| + а V xt ln xt ^ min. (2)
2 ■ 1 x
г=1
Поскольку целевая функция (2) является нелинейной, нельзя записать простое аналитическое выражение для вектора x , доставляющего ей минимум, как это было сделано для метода Тихонова. В работе [2] доказана выпуклость целевой функции (2), поэтому для поиска ее минимума можно использовать любой метод локальной оптимизации нелинейных функций векторного аргумента. Следует отметить, что целевая функция (2) предполагает положительность всех элементов вектора x .
Метод регуляризации посредством ограничения числа итераций [5] основан на минимизации нелинейной целевой функции следующего вида:
m 1 n n
J (x) = V _ V aijxj- ytln V aijxj ^ min. (3)
i=1 yi j=1 j=1 x
Целевая функция (3) имеет единственную точку минимума [5]. Регуляризирующий эффект данного метода обусловлен:
а) наличием в (3) слагаемого, подобного энтропии,
aijXj;
Уi 1п
j=1
б) ограничением количества итераций процесса минимизации.
Верхняя граница количества итераций зависит от особенностей конкретной решаемой задачи и, вообще говоря, определяется эмпирически. Данный факт вносит в метод некоторую «нестрогость», однако в работе [5] показана его эффективность для решения ряда задач.
В 2005 г. А.И. Жданов опубликовал метод регуляризации несовместных СЛАУ (в том числе и неполного ранга), основанный на преобразовании несовместной СЛАУ к эквивалентной расширенной совместной СЛАУ с симметричной матрицей [4]. К полученной СЛАУ применяется теорема, доказанная В. А. Морозовым и С.Ф. Ги-лязовым, которая дает способ выбора параметра регуляризации для метода А.Н. Тихонова [2].
Исходная СЛАУ преобразуется к следующей:
(4)
где г = у - Ах .
Методы - и £р -регуляризации [3] являются модификацией
метода регуляризации А. Н. Тихонова, предназначенной для решения некорректных задач с разреженным вектором решения. Целевые функции для методов £1 - и £ р -регуляризации имеют вид
ц2
" I А" г у
=
_ А т 0 _ х 0
JA (х, а) = 11 Ах - у ||2 + а ||х|1 ^ mm;
J£р (х, а) = ||Ах - у||2 +а||х||p ^ гот, р < 1,
(5)
(6)
где х
г п , V/* k
Е1
V i=l
В работе [3] обоснована эффективность методов £1 - и £р-
регуляризации применительно к линейным задачам, в которых истинный вектор решения является разреженным. Под разреженным вектором решения будем понимать тот факт, что большинство его элементов являются нулевыми (имеют очень малые значения) и лишь несколько элементов имеют большие значения. Если разреженный вектор представить графически, он будет иметь несколько резких пиков.
Разделение методов на £1 и £ р с р < 1 связано с тем, что, как показано в [3], целевая функция (5) является выпуклой, чего нельзя ска-
к
зать о целевой функции (6). Следовательно, возникает вопрос о выборе начального приближения для £ p -метода, так как для минимизации целевых функций применяют методы локальной оптимизации (они просты в реализации и обладают гораздо меньшей вычислительной сложностью, чем £ p -метод). В работе [3] в качестве начального приближения для £ p -метода используют решение, полученное
более простым методом формирования пучка.
Для вычисления ковариационной матрицы решения воспользуемся методом, основанным на теореме Крамера — Рао [6]. Следует отметить, что в общем случае целевые функции (5), (6) не являются непрерывно дифференцируемыми из-за наличия в их выражениях знаков модуля.
Существуют и другие методы решения некорректных задач [1-5].
Необходимость выбора параметра регуляризации, показателя p в методе £ p -регуляризации, невыпуклость целевой функции заставляют искать другие подходы к решению некорректных задач.
Рассмотрим решение некорректной задачи с помощью многокритериального математического программирования [6]. Двух-критериальная задача математического программирования в данном случае записывается следующим образом:
J1 Az - u||2 ^ min,
* 2 (7)
J2 = V zt ^ min, zt > 0
i=1 z
при ограничениях
Az = u.
Для аналитического решения многокритериальная задача сводится к однокритериальной с помощью пороговой оптимизации (метод e -ограничений), используется целевое программирование (архимедова модель и модель с приоритетами) [6].
Метод пороговой оптимизации (метод e -ограничений) приводит к различным возможным комбинациям целевых функций и ограничений:
min || Az - u||2 при ||z||p <8 ; (8)
min II Az - u|| 2 при ||z||p <8 ; (9)
при ||Az - u||2 <ß; (10)
• II IIP N л l|2 min z 11 A~
ii "p
min||при ||Az - u||2 <Vß. (11)
В настоящей статье основное внимание уделено задачам (8) и (10). Задача (8) при р = 1 может быть решена методами квадратичного программирования, задача (10) решается методами нелинейного программирования [6].
Для решения двухкритериальной задачи (7) применено целевое программирование (архимедова модель и модель с приоритетами) [6].
Во всех рассматриваемых методах определяются точечные оценки решения. В реальных условиях параметры систем и наблюдаемые величины измерены или заданы с определенной погрешностью. Это приводит к тому, что от реализации к реализации решение будет меняться. Чтобы судить об эффективности того или иного метода и учитывать возможный разброс решений, необходимо получать не только точечные, но и интервальные оценки решения.
Дисперсии оценок амплитуд были получены методом статистических испытаний и аналитически. Число испытаний для задачи квадратичного программирования — 10 000, для задачи нелинейного программирования — 100, для задачи целевого программирования — 100 (архимедова модель) и 30 (модель с приоритетами), для задачи нелинейного программирования с учетом энтропии — 100. Аналитически ковариационная матрица оценок (дисперсии) определялась путем обращения матрицы вторых производных логарифма функции правдоподобия, полученной на основе исходной системы нелинейных уравнений, при найденных значениях оценок.
Продемонстрируем этот подход на примере задачи идентификации тайно проведенного ядерного взрыва [5, 7, 8], когда зарегистрировано только два изотопа ксенона.
В качестве практически оправданного допущения для предлагаемого способа измеренные активности изотопов рассматривают как детерминированные, подверженные аддитивной помехе, оценки параметров которых подлежат определению. Числа обусловленности в
этих задачах достигают 1026 .
При мгновенном делении / -й изотоп появляется в результате
различных видов деления, и его измеренная активность Д (() выражается следующим образом [5, 7, 8]:
т _
^ ау (9, л, А, 1, 1Ч= Д (), I = 1, п, (12)
j=1
где а^ (9, л, А, (, ) — активность I -го изотопа при j -м виде деления
для одного акта распада, вычисленная с учетом сепарации на момент времени I > ^, т. е. удельная активность; 9 — вектор параметров,
характеризующих сепарацию измеряемых изотопов от предшеству-
ющих им; ц — вектор независимых выходов изотопов (при - -м виде деления); X — вектор постоянных распада; Х — время наблюдения; ^ — предполагаемый момент сепарации изотопов криптона и
ксенона от предшествующих им изотопов по цепочкам радиоактивных превращений; р — доля 1 -го изотопа в образце (значение р обычно неизвестно); Nj — число делений - -го вида.
До момента сепарации удельная активность аП (9, Ц, X, ) определяется формулой [5]
а*- (9, ц, X, ^) = <
Ртх ПР
Пр -1
ц„х„ ехр(-Х„??) + Ц Цр Х1р
р=11Р =1
"--1 п ехр^д
П гр \ I
•р 'Р ^^ "
Зр =
1р П (Ч -Ч)
Яр =1р
Яр **Р
,(13)
где Ц1 — независимый выход 1 -го изотопа; пр — номер исследуемого
изотопа по р -й ветви; п — максимальный член из { п
{пр}; рп
число
ветвей цепочки; (пр -1) — число изотопов, предшествующих исследуемому по р -й ветви распада; у Гр — доля г -го члена цепочки, полученного из (г -1) -го члена по р -й ветви; Хг , X г , Xз , X ц — постоянные распада изотопов, имеющих соответственно номера 1р, г р, зр, цр по р -й ветви, причем г р < г р < п р -1, г р < з р < п р, г р < цр < п р и цр Ф зр;
— время, когда произошло мгновенное отделение исследуемого изотопа от предшественников, после чего распад изотопа идет по экспоненте с постоянной распада X п.
После момента сепарации изотопы распадаются по своим постоянным распада X,:
а- (9, ц, X, Х, ) = а- (9, ц, X, )ехр(-X!(Х -Х*)), 1 = Щ - = Гт (14)
где а- (9, Ц, X, Хц) — удельная активность, рассчитанная по формуле (13) на момент времени .
Уравнения вида (12) составляются для каждого измеряемого изотопа ксенона, в результате формируется СЛАУ
а11(*д, *) а12(*д, *) а21(*д, *) а22(*д, *)
ап1(*а 5 *) ап2 , *)
а1т (*д, *) " рМ " " Д(*)"
а2т (*д, *) р#2 = и)
апт (*д, *) _РМт _ _ А (*) _
(15)
в которой определению подлежат неизвестные вклады источников радиоактивности рМ в суммарную активность изотопов криптона и
ксенона.
Первый этап решения задачи идентификации источников радиоактивных благородных газов (РБГ) — определение времени сепарации изотопов криптона и ксенона. Временной отрезок, которому принадлежит момент сепарации, можно найти, «достроив» относительные активности изотопов в различных видах деления «в обратном времени» от момента измерения без учета влияния предшествующих им изотопов и определив точки пересечения линий, проведенных из экспериментальных точек, с относительными активностями, построенными с учетом влияния предшествующих изотопов по цепочке распада изотопов.
А (Хе133т) !А (Хе133)
10
-1
10
-2
10
10"
-з
10
-5
- 1
-// 2 / \
1 / = з ч
В
8 10 ч
Рис. 1. Графики относительной активности для
133т
двух изотопов ксенона (Хе
1 — и235; 2 — Ри239
Хе133):
На рис. 1 приведены графики относительной активности для
133 133
двух изотопов ксенона (Хе т, Хе ), экспериментальная точка В соответствует моменту измерения активностей * = 12 ч после события. Чтобы не усложнять рисунок, изображены только «граничные»
235
линии, соответствующие делению урана 235 и плутония 239 (и у и
239
Puf ) нейтронами спектра делений (вместо шести возможных видов
тт235 тт235 тт238 тт238 „ 239 т> 239 ч тт - л л
деления: иу , и14 , иу , и14 , Риу и Ри14 ). Нижний индекс 14
определяет нейтроны с энергией 14 МэВ. Как видно на рис. 1, момент сепарации принадлежит интервалу от ^ = 3 ч до ^ = 4 ч после события.
Задавая значения времени внутри отрезка [хь Х2 ] и решая систему (15) для моментов Хц, соответствующих точкам отрезка, в качестве момента сепарации принимают время Х™", для которого сумма
квадратов невязок системы (15) минимальна.
Второй этап решения задачи идентификации ядерного взрыва — определение для каждого фиксированного момента сепарации Хц
оценок числа делений рN-. При заданном Хц система (15) является
линейной относительно неизвестных рN -, - = 1, т . Поскольку элементы А = {а- (Хц, Х)} не могут быть точно рассчитаны (независимые
выходы известны с погрешностями) и активности изотопов А1 (Х) также измеряются с ошибками, будем считать, что элементы А и измеренные активности А1 (Х) — независимые случайные величины, распределенные по нормальному закону с математическими ожиданиями, равными соответственно а^^ц, Х) и А/ше(Х), и дисперсиями,
.2 / .Л „ „2 |
равными соответственно а2 (а- (Хц, Х)) и а2 (А1 (Х)):
(16)
а- (Хц, Х) = ауие(Хц, Х) ±е-,
А (Х) = АГ(Х) ±§1,
где а17гие(Хя, Х), А/ше(Х) — истинные значения соответственно удельных и измеренных активностей изотопов (которые нам неизвестны); е- — погрешности определения удельных активностей а- (Хц, Х); 8 г — ошибки измерения активностей А (Х) РБГ в атмосфере.
Для учета погрешностей как в измеренных активностях А (Х), так и в элементах {а-(Хц, Х)} используется определение ортогональной регрессии (конфлюэнтный анализ) [5], и в силу независимости случайных величин а- (Хц, Х) и А (Х) можно записать функционал следующего вида:
1 п =2 г
2 г=1
ау (г9, X) - , X)
1гие /
У ^
О'
\ау Од,Х))
417)
где рЫ у, у = 1, да, — подлежащие определению вклады источников
радиоактивности в суммарную активность; а^ , X) — неизвестные
точные значения удельных активностей, оценки которых определяются в процессе нахождения рЫу ; агу (Хд, X) — удельные активности,
рассчитанные по формулам (13), (14), по имеющим погрешности независимым и кумулятивным выходам элементов изобарных цепочек радиоактивных превращений; Д (X) — измеренные в пробе активности РБГ.
В точке минимума (17) должны выполняться следующие условия:
Э(рЫу)
Э^
= 0, У = 1, да; (18)
ры,=Ы
Эа<Т(^, X)
= 0, г = 1, п. (19)
а, / )=аГ(',, ' )
Несмотря на линейность систем уравнений (18), (19) при фиксированном , задача является вычислительно н е к о р р е к т н о й
в силу плохой обусловленности системы (18). Отношение максимального и минимального значений собственных чисел матрицы системы (18) достигает порядка 1026.
Сначала при аУ^^д, X) = ау (1Ц, X) решают СЛАУ (18) методами
многокритериального математического программирования (метод сжатия области допустимых значений, целевое программирование) и
находят первое приближение оценки (рЫу )1.
Далее для получения оценок истинных значений а*™^», X) при
заданном на каждом шаге вычисления оценок рЫу , у = 1, да, используется условие (19), что приводит к решению дополнительно п систем линейных уравнений с да неизвестными следующего вида:
I
(рЫг ) ^ , , г) _ а* , г) , (р^)А (г)
г=1 а
(4(г))
-а1г ,г) + -
а 2(а^)
^-;--+ -г,-
(а1у (г9, г)) а2 ( 41))
I _ 1, и; V _ 1, т.
Полученные оценки истинных значений а*™6^, г) должны принадлежать области неопределенности измеренных значений ау (гд, г) :
а
(гд, г) - аЦие (гя, г) < 3а (ау (гя, г)), I _ 1, п; у _ 1, т.
Если это условие не выполняется, то а*™6^, г), у _ 1, т, которые не
удовлетворяют данному неравенству, следует заменить на значения ближайших граничных точек. Из-за этого может происходить увеличение значений функционала / на новых точных значениях переменных по сравнению с предыдущим шагом итерационного процесса, что приводит к снижению скорости сходимости итерационного процесса или к возникновению колебаний. Чтобы значения функционала не увеличивались после пересчета оценок а*™6^, г), у _ 1, т,
элементы из набора оценок а*™6^, г), у _ 1, т, для которых произошло увеличение соответствующих слагаемых функционала ¥1 по сравнению с их значениями на предыдущей итерации, следует заменить на соответствующие значения а*™6^, г) для предыдущего шага.
После определения оценок истинных значений а*™6^, г) находят
очередное приближение рру к решению р/, у _ 1, т, методами
многокритериального математического программирования. Критерием останова алгоритма служит несущественное различие значений соответственно компонентов вектора р/ и функционала / на соседних итерациях, т. е. выполнение неравенств
Р>Ni) -(рЫ
I-1
рЫ;
<У 1,
Р\- Р\-1
<У 2,
(20)
где (рру) — очередное приближение к решению на I -й итерации;
у 1, у 2 — некоторые числа (малые десятичные дроби, например 0,001), определяющие точность вычисления оценок рру .
При решении методами многокритериального математического программирования:
1) формируют двухкритериальную задачу математического программирования:
( Л2
m
true
j=z a (t) - Z ar^, )
i=1 v j=1
m
J2 = Z PNj ^ mi
^ min,
pNj
j (21)
min
j=1 ^
при ограничениях pNj > 0, j = 1, m ;
2) используя метод пороговой оптимизации или целевое программирование, от двухкритериальной задачи математического программирования (21) переходят к однокритериальной задаче посредством перевода всех, кроме одного, из указанных выше функционалов в условия ограничений.
Метод пороговой оптимизации (или метод e -ограничений) приводит к различным возможным комбинациям целевых функций и ограничений. В алгоритме используют следующие их виды:
2m
n m m _
minZ A(t)-Z je(tq,t)(pNj) при ZPNj <8; pNj >0, j = 1,m, (22)
pNj i =11 j=1 J j=1
n ( m ^
minzpNj приz a(t)-Z4T(tq,t)(pNj)
pNjj=1 i=1 v j=1
<ß; pNj > 0, j = 1, m. (23)
Задача (22) является задачей квадратичного программирования, задача (23) — задачей нелинейного программирования.
Оценки правых частей ограничений 8 и ß могут быть получены при независимой минимизации функционалов J1 и J2 при ограничениях pNj > 0, j = 1, m . При этом может применяться любой метод
математического программирования.
В целевом программировании существует две модели решения — архимедова и модель с приоритетами.
При использовании архимедовой модели все целевые функции переводят в ограничения и осуществляют минимизацию взвешенной суммы меры их отклонений от ограничений:
min { -(w1d1 + w2d2 )} при J1 (pNy) + d1 < ß, J2 (pNy) + d2 <8, (24)
di, d2
где — весовые коэффициенты, ^ = 1; — отклонения от
г=1
ограничений.
В модели с приоритетами осуществляют последовательный перевод целевых функций в ограничения и минимизацию отклонения значений целевых функций от ограничений. При этом найденное на данном шаге значение отклонения ^ используют как оптимальное отклонение на следующем (г +1 )-м шаге:
шаг 1: тт(-ё) при ^ (рЫ-) + ё1 < Р;
шаг 2: тт(-ё2 ) при 31 (рЫ- ) + ё1опт < Р, 32 (рЫ- ) + ё2 < §.
ё2
При идентификации ядерного взрыва по малому числу изотопов ксенона используется тот факт, что выход этих изотопов слабо зависит от энергии нейтронов. Это иллюстрирует рис. 2, на котором представлены графики усредненных кумулятивных выходов изо-
133 135
топов Хе и Хе в зависимости от долей кумулятивных выходов, соответствующих нейтронам спектра деления и нейтронам с энергией 14 МэВ. Кумулятивные выходы известны с погрешностями до 5 %.
*1к. %
Рис. 2. Зависимость выходов изотопов ксенона от энергии нейтронов (в диапазоне от энергии нейтронов спектра деления до энергии 14 МэВ):
1 — изотоп Хе133, материал и235 ; 2 — изотоп Хе , материал Ри2; 3 — изотоп Хе ,
материал И ; 4 — изотоп Хе , материал Ри239
На рис. 2 видно, что усредненные значения выходов (для примера приведены значения, соответствующие равным долям нейтронов спектра деления и нейтронов с энергией 14 МэВ) укладываются в эти погрешности.
Для идентификации взрыва по двум измеренным изотопам применяют объединение двух видов деления урана 235 (нейтронами спектра деления и нейтронами с энергией 14 МэВ) в один вид деления и объединение двух видов деления плутония 239 (нейтронами спектра деления и нейтронами с энергией 14 МэВ) в один вид деления, что приводит к сокращению числа идентифицируемых видов деления (вместо четырех рассматривают два).
При этом удельную активность рассчитывают по формуле (13): • для делящегося материала урана 235 с вектором независимых выходов
,и235 „„и2/5 <Л \ и,235
ли = 1 +(1 - d )л
у235 и235
где л 1 и Л 14 — независимые выходы элементов изобарных цепочек при делении урана 235 нейтронами спектра деления и нейтронами с энергией 14 МэВ соответственно; c — параметр, учитываю-
U235 у235 р ,
щий доли независимых выходов л 1 и л 14 в их сумме, ci е [ 0, 1J;
• для делящегося материала плутония 239 с вектором независимых выходов
р,, 239 р, 239 р,239
Лри = С2ЛРи1 + (1 -С2)лр"14 ,
р,2139 р,12439
где л 1 и л 14 — независимые выходы элементов изобарных цепочек при делении плутония 239 нейтронами спектра деления и нейтронами с энергией 14 МэВ соответственно; С2 — параметр, учи-
ри2139 ри12439
тывающий доли независимых выходов л 1 ил 14 в их сумме,
С2 е[0, 1].
Задавая двумерную сетку по q и С2 с шагом Aq и Дс2 соответственно и находя минимум (17) для разных q и С2, за истинные
U235 ри239
вклады источников pN и pN в суммарную активность изотопов ксенона принимают те, при которых сумма квадратов невязок системы (13) м и н и м а л ь н а.
Проведено имитационное моделирование реализации предлагаемого способа на персональном компьютере с процессором Intel Celeron 2,40 ГГц с объемом оперативной памяти 768 Мбайт в математическом пакете Matlab 7.0.
Имитировали ситуацию отбора пробы через 6 дней после взрыва
ОГ 1 Л 1 1 лл
и измерения активностей пяти изотопов (Kr m, Xe m, Xe m,
133 135
Xe , Xe ). Результаты рассчитаны при условии, что момент сепарации предполагается известным: 3 ч после события. Значения изме-
ряемых активностей были аддитивно «зашумлены» гауссовым шумом с СКО, равным 5 % их «точного» значения.
Рассматривались следующие комбинации видов деления (возможные наборы переменных pNj):
235 133
1) Uh + фон по Xe (два неизвестных источника);
2) и|35 + Uf + U^5 (три неизвестных источника);
235 239 239
3) Uh + Pu f + Pu14 (три неизвестных источника);
235 235
4) U f + U14 (два неизвестных источника);
235 235 133
5) U f + Uf4 + фон по Xe1JJ (три неизвестных источника);
239 239
6) Pu^ + Puff (два неизвестных источника);
239 239 133
7) Pu f + Puff + фон по Xe1JJ (три неизвестных источника);
8) U^5 + U^5 + Puf + Pu^9 (четыре неизвестных источника);
m тт235 , тт235 , т> 239 , „ 239 , , v 133 ,
9) U f + U14 + Pu f + Pu14 + фон по Xe (пять неизвестных источников).
Здесь U^5 — реакторный выброс (данные по реакторам взяты из справочной литературы).
Истинным решением является четвертая комбинация, отно-
235 235
сительный вклад Uf равен 100, относительный вклад U14 равен
100. Результаты моделирования сведены в таблицу.
В графе «Метод решения» указаны четыре метода решения задачи идентификации ядерного взрыва (квадратичное программирование, нелинейное программирование, архимедова модель, модель с приоритетами), предлагаемые в данном способе, которые сравнивались с методом решения, используемым в аналоге (регуляризация Тихонова).
В графе «Номер комбинации видов деления» указан номер комбинации видов деления, обеспечившей для соответствующего метода решения (регуляризации Тихонова, методов многокритериального математического программирования) из всех девяти комбинаций наименьшую сумму квадратов невязок системы уравнений (15).
В графе «Порядок числа обусловленности матрицы системы» указаны порядки чисел обусловленности матрицы системы (15), соответствующей номеру комбинации видов деления, приведенному во второй графе таблицы.
В графе «Оценка решения» перечислены оценки вкладов видов деления, присутствующих в указанных во второй графе таблицы комбинациях видов деления. Например, для регуляризации Тихонова
наилучшей с точки зрения суммы квадратов невязок является комби-
235
нация 2. Этой комбинации соответствуют три вида деления: , рассчитанный относительный вклад которого 35,99 (деление тепло-
235
выми нейтронами); и у , рассчитанный относительный вклад кото-
235
рого 42,97, и И14 , рассчитанный относительный вклад которого 110,27. Аналогично для остальных методов решения.
В графе «Сумма квадратов невязок» приведены значения суммы квадратов невязок системы уравнений (15), рассчитанных для указанных в таблице комбинаций видов деления и оценок их вкладов в суммарную активность изотопов криптона и ксенона.
В графе «Время работы алгоритма, мин» указано в минутах время получения оценки вкладов соответствующим методом.
Результаты решения задачи идентификации разными методами (источником радиоактивности являются И235 и И2|5 , точное решение 100 и 100)
Метод решения Номер комбинации видов деления Порядок числа обусловленности матрицы системы Оценка решения Сумма квадратов невязок Время работы алгоритма, мин
Регуляризация Тихонова 2 106 35,99; 42,97; 110,27 751,39 1,46
Квадратичное про-граммиро-вание 9 1014 84,94; 108,18; 0; 0; 10,99 74,97 9,28
Нелинейное программирование 9 1014 43,53; 120,76; 0; 0; 0 5096,43 11,95
Архимедова модель 9 1014 84,94; 108,18; 0; 0; 10,98 74,97 13,40
Модель с приоритетами 9 1014 84,94; 108,18; 0; 0; 10,99 74,97 19,78
Из данных таблицы видно, что метод регуляризации Тихонова дал отрицательный результат — в решении присутствует значительный относительный вклад атомного реактора (которого нет в истинном решении). В методах многокритериального программирования (квадра-
тичное, нелинейное, целевое программирование (архимедова модель и модель с приоритетами)) при использовании дополнительного условия на неотрицательность переменных получено положительное решение. Оптимальный результат соответствует девятой комбинации видов деления. Это не противоречит истинному решению, так как вклады тех видов деления, которых не было в истинном решении, незначительны (большинство вкладов равны нулю).
Таким образом, преимуществом предлагаемого способа является повышение эффективности и достоверности решения некорректных задач, что принципиально важно в военном деле, космонавтике и авиации.
ЛИТЕРАТУРА
[1] Тихонов А.Н. О решении некорректных задач и методе регуляризации. Докл. АН СССР, 1963, т. 161, № 3, с. 501-504.
[2] Морозов В.А. Алгоритмические основы методов решения некорректных задач. Вычислительные методы и программирование, 2003, т. 45, с. 130-141.
[3] Malioutov D.M. A Sparse Signal Reconstruction Perspective for Source Localization with Sensor Arrays. Master of Science thesis. Massachusetts Institute of Technology, 2003.
[4] Жданов А.И. Регуляризация неустойчивых конечномерных линейных задач на основе расширенных систем. Журнал вычислительной математики и математической физики, 2005, т. 45, № 11, с. 1919-1927.
[5] Грешилов А.А. Некорректные задачи цифровой обработки информации и сигналов. 2-е изд. Москва, Логос, 2009, 360 с.
[6] Грешилов А.А. Математические методы принятия решений. 2-е изд. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2014, 645 с.
[7] Грешилов А.А., Лебедев А.Л. Способ идентификации ядерного взрыва по изотопам криптона и ксенона. Пат. № 2407039 Российская Федерация, 2010, бюл. № 35, 21 с.
[8] Грешилов А.А., Тетюхин А.А. Алгоритм идентификации источников радиоактивных благородных газов. Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки, 2003, № 2, с. 3-19.
Статья поступила в редакцию 10.12.2014
Ссылку на эту статью просим оформлять следующим образом:
Грешилов А.А. Некорректные задачи и многокритериальное программирование. Инженерный журнал: наука и инновации, 2015, вып. 2. URL: http://engjournal.ru/catalog/arse/itae/1367.html
Грешилов Анатолий Антонович родился в 1939 г., окончил Московский инженерно-физический институт в 1964 г. Д-р техн. наук, профессор кафедры «Высшая математика» МГТУ им. Н.Э. Баумана. Автор более 200 научных работ, в том числе более 30 монографий, 30 авторских свидетельств и патентов. e-mail: [email protected]
А.А. rpernmoB
Ill-posed problems and multicriteria programming
© A.A. Greshilov
Bauman Moscow State Technical University, Moscow, 105005, Russia
Solving ill-posed problems by methods of multicriteria mathematical programming has been considered. Several methods of multicriteria mathematical programming (method of compression of acceptance region and goal programming) are used simultaneously allowing considering additional types of restrictions (nonnegativity of the solution, bound-edness of solution) which must be met by evaluation of solution and which do not require definition of the regularization parameters necessary in the classical methods of regular-ization. When registering a small number of isotopes the merger of the two types of uranium- 235 instant fission into one kind of division and two types ofplutonium-239fission into one kind of division is used. Simultaneously different variants of the nuclear explosion mechanism are considered. Determination of contributions of different fission kinds into the total activity of isotopes of krypton and xenon is performed by formation of a functional for a given moment of separation t and time measurements t of functions Ft.
They are obtained from the functional by its differentiation on the elements pNj of the
system of linear algebraic equations (SLAE) at fixed values of the specific activity atjue (tq, t). Solving SLAE is performed by generating multiple objective functions and
using the 4 methods of multicriteria mathematical programming, reducing multicriteria problem to one-criterion problem with constraints. The solutions of the mentioned one-criterion problem with constraints are obtained by the iterative computational procedures with the given specific activity dtrue (tq, t) and when its assessment was specified
for each iteration. The point estimates of the contributions of different fission kinds into the total activity of isotopes are determined. For definition of the moment of separation tq contributions pNj are calculated for different values t and the value at which the
( \2 n m m 2
ratio £ A (t)- Y af* (tq, t )(pNj) / Y (pNj) is minimal is chosen.
i =11 J =1 ) J =1
Keywords: regularization methods, fission products, nuclear explosion, separation of isotopes, multicriteria programming, targeted programming, compression method, iterative method of solving
REFERENCES
[1] Tikhonov A.N. Doklady Akademii nauk SSSR - Reports of the USSR Academy of Sciences, 1963, vol. 161, no. 3, pp. 501-504.
[2] Morozov V.A. Vychislitelnye metody I programmirovanie: Novye vychislitelnye tekhnologii - Computational Methods and Programming: New Computing Technologies, (Electron. Edition). 2003, vol. 4, pp. 130-141.
[3] Malioutov D.M. A Sparse Signal Reconstruction Perspective for Source Localization with Sensor Arrays. IEEE Transactions on Signal Processing, 2005, vol. 53, no. 8, pp. 3010-3022.
[4] Zhdanov A.I. Zhurnal vychislitelnoy matematiki i matematicheskoi fiziki RAN — Journal of Computational Mathematics and Mathematical Physics RAS, 2005, vol. 45, no. 11, pp. 1919-1927.
[5] Greshilov A.A. Nekorrektnye zadachi tsifrovoy obrabotki informatsii i signalov [Some Ill-Posed Problems of Digital Information and Signal Processing]. Moscow, Universitetskaya kniga; Logos Publ., 2009, 360 p.
[6] Greshilov A.A. Matematicheskie metodyprinyatiya resheniy: Uchebnoe posobie dlya vuzov [Mathematical Methods of Decision-Making: study book for higher school]. Moscow, BMSTU Publ., 2014, 645 p.
[7] Greshilov A. A., Lebedev A.L. Sposob identificatsii jadernogo vzriva po isoto-pam kriptona I ksenona. [A Method for Determining the Concentration of Inert Gas Isotope in the Mixture of Fission Products]. Patent № 2407039 Russian Federation, 2010, bulletin № 35, 21 p.
[8] Greshilov A.A., Tetukhin A.A Vestnic MGTU im. N.E. Baumana. Seria Estestvennye nauki — Herald of the Bauman Moscow State Technical University. Series: Natural Sciences, 2003, no. 2, pp. 3-19.
Greshilov A. A. (b. 1939) graduated from Moscow Engineering Physics Institute, Department of Experimental and Theoretical Physics, in 1964. Doctor of Engineering Sciences, professor at the Department of Higher Mathematics at Bauman Moscow State Technical University. The author of more than 200 scientific papers, including more than 30 monographs, 30 patents and Certificates of Authorship in the field of the development of mathematical methods of considering uncertainty of the initial information in the problems of mathematical physics, pattern recognition, forecasting, and other technical applications. e-mail: [email protected]
