Содержание к диссертации
Введение
Глава 1. Итерационный проекционный метод С. Качмажа 20
1.1. Метод алгебраической реконструкции 20
1.1.1. Метод Качмажа как итерационный метод уменьшения ошибки 26
1.1.2. Постановка задачи реконструктивной компьютерной томографии 30
1.2. Влияние последовательности выбора строк в алгоритме алгебраи
ческой реконструкции 37
1.2.1. Алгоритм упорядочивания строк с использованием четной перестановки 39
1.2.2. Рандомизированный алгоритм Качмажа 43
1.2.3. Классификация и обобщение способов выбора проекционной последовательности 46
1.3. Программные реализации модификации классического алгоритма
Качмажа и вычислительный эксперимент 48
1.3.1. Классический алгоритм Качмажа с циклическим правилом 50
1.3.2. Классический алгоритм Качмажа с квазиоптимальным правилом 52
1.3.3. Классический алгоритм Качмажа с использованием правила четной перестановки 53
1.3.4. Классический алгоритм Качмажа с использованием рандомизированного правила 54
1.3.5. Результаты вычислительных экспериментов для задачи компьютерной томографии с параллельной схемой сканирования 56
1.4. Заключение
Глава 2. STRONG Проекционный метод для решения задачи регуляризации решений
линейных систем STRONG 61
2.1. Регуляризованный итерационный проекционный метод 61
2.2. Рандомизированный регуляризованный метод алгебраической реконструкции 67
2.3. Квазиоптимальный регуляризованный метод алгебраической реконструкции 73
2.4. Результаты вычислительных экспериментов для задачи компьютерной томографии с параллельной схемой сканирования в присутствии возмущений 75
2.5. Заключение 78
Глава 3. Решение интегральных уравнений Фредгольма 1-го рода на основе расширенных регуляризованных нормальных уравнений 80
3.1. Понятие корректно поставленной задачи 80
3.2. Корректность постановки вычислительной задачи 84
3.3. Метод регуляризации А.Н. Тихонова для задачи нормального псевдорешения С.Л.А.У 88
3.3.1. Выбор параметра регуляризации при решении с.л.а.у. со множе ством правых частей 91
3.4. Решение дискретизированного интегрального уравнения Фредгольма с использованием регуляризованного алгоритма Качмажа. 93
3.4.1. О решении задачи Филлипса 93
3.4.2. О решении еще одной тестовой задачи реконструкции изображения 97
3.5. Заключение 100
Глава 4. Регуляризирующие алгоритмы вычисления оценок двумерных им пульсных характеристик искажающих систем 102
4.1. Об эффективной индексации двухуровневых теплицевых матриц 105
4.2. Постановка задачи идентификации 107
4.2.1. Приведение обобщенной задачи регуляризации к каноническому виду 111
4.2.2. Метод расширенных регуляризованных систем для обобщенной задачи регуляризации 112
4.3. Правило останова итерационных алгоритмов как параметр регуляризации 114
4.4. Решение задачи идентификации с использованием останова итерационного алгоритма по специальным правилам 118
4.5. Заключение 122
Заключение 123
Список литературы
- Постановка задачи реконструктивной компьютерной томографии
- Рандомизированный регуляризованный метод алгебраической реконструкции
- Метод регуляризации А.Н. Тихонова для задачи нормального псевдорешения С.Л.А.У
- Приведение обобщенной задачи регуляризации к каноническому виду
Постановка задачи реконструктивной компьютерной томографии
Метод алгебраической реконструкции [4, 5, 6] (Kaczmarz method, ART l, PCS2) - полностью последовательный строковый3 метод с длинной историей, богато описанный в литературе (библиография [112]). В классическом виде метод был впервые сформулирован и предложен польским математиком С. Качма-жем (S. Kaczmarz) для решения систем линейных алгебраических уравнений с квадратными и невырожденными матрицами [7, 56], а впоследствии независимо был использован Р. Гордоном и др. (R. Gordon; R. Bender; Herman [4, 6]) для реконструкции изображений по их проекционным данным в задаче рентгеновской компьютерной томографии.
В более общем случае последовательность, полученная по методу итераций С. Качмажа, также сходится к решению и для переопределенных, но обязательно совместных с.л.а.у Аи = f . Итерации в данном методе принято где к — номер итерации, и0 — начальное приближение, а Л — параметр4 релаксации [24]. В [56] показана фундаментальная связь данного метода с методами Гаусса-Зейделя и методом, предложенным de la Garza в [100], данная тема также исследуется автором в работах [62, 63]. Возможна также и другая форма записи таких итераций. Для этого рассмотрим отображение gi (и) : Rn — Rn, определяемое, как в [54], где / — номер макроитераций. Важно, что вычисления с использованием формулы (1.5) могут быть преобразованы к виду стационарного итерационного построен для случая, когда Afc = 1, У к как и везде, далее, если не оговорено иное. процесса первого рода, как это показано в [28]. Там же указывается, что последовательность {и1}5, каждый элемент которой определяется через (1.5), сходится для с.л.а.у. произвольного вида. Это доказывается следующей важной теоремой о сходимости. [28] Для любой6 матрицы А є Rmxn не содержащей нулевых строк, и для любого вектора f Є Кт итерационная процедура (1.5) генерирует последовательность {и1} векторов в Rn, которая сходится к
Если же не производить разделение на макро-и микроитерации, а напрямую исследовать итерационную формулу (1.1) для решения систем произвольного вида, то вопрос относительно сходимости последовательности {м } не имеет простого ответа. Упрощенный пример с несовместной с.л.а.у. (рассматривается ниже) является достаточно типичным для задач, в которых решение получается на основе серии некоторых измерений, каждое из которых может включать как устранимую, так и неустранимую погрешности.
Важно отметить, что для случая, когда с.л.а.у. совместна7 и матрица коэффициентов квадратная, в [27, 28] указывается, что алгоритм (1.1) также сходится, если
На практике установить совместность с.л.а.у. часто оказывается сложнее, чем найти само решение. Более того, исходная точная совместная задача после ее задания в памяти вычислительного устройства приобретает неустранимую погрешность, которая может приводить к тому, что система становится несовместной, показателен в этом отношении пример из Главы 3 — система (3.8). Этот результат обобщается в [29] для процедуры (1.5): доказывается, что в несовместном случае при использовании параметра релаксации, удовлетворяющего условию цикличности Xk = Xi-i и ограничениям (1.6), также присутствует сходимость последовательности {и1}.
Порядок обхода строк матрицы с.л.а.у (иначе — проекционная последовательность) в общем случае определяется сюрьективным отображением где т — количество уравнений в системе, Щ = N U {0}, а при применении формулы (1.2) является последовательным — от первой строки до последней, а с последней — на первую8. Безусловной является всю важность выбора параметра релаксации Л для сходимости, а также скорости сходимости итераций по методу Качмажа, при этом в настоящей работе основной акцент делается на различные модификации алгоритма (1.2) за счет управления способом задания последовательности j (к). Выбор последовательности строк при решении системы уравнений является важной частью алгоритма и предметом для построения различных модификаций формулы (1.2) с целью ускорения сходимости алгоритма.
Метод алгебраической реконструкции достаточно просто формулируется с использованием его геометрической интерпретации. Рассмотрим совместную с.л.а.у. с двумя неизвестными
Рандомизированный регуляризованный метод алгебраической реконструкции
Как отмечалось ранее, под проекционной последовательностью понимается порядок обхода уравнений системы в процессе итераций алгоритма Качма-жа. Рассматриваемая последовательность задается с использованием правила j (к), где к Є Щ, которое позволяет вычислить номер решаемого уравнения на итерации с номером к. Правило j (к) можно рассматривать как сюрьективное отображение множества Щ на множество {1,2,..., т}\ где т — количество уравнений в решаемой с.л.а.у. Аи = /, А є Rmxni и Є ІГ, f є Rm. Отображения вида (1.41) могут быть как детерминированными, так и случайными, рандомизированными. В последнем случае j (к) является дискретной случайной величиной с заданным законом распределения. Например, к детерминированным отображениям необходимо отнести:
В свою очередь детерминированные отображения могут задаваться рекуррентными правилами, например, квазиоптимальное правило является, рекуррентным так как формально для вычисления j (к) необходимо знать все предшествующие значения последовательности j (к — 1), j (к — 2),..., j (0). В то же время при использовании правила циклического обхода предыдущие значения Известные правила
Существенно важной является классификация рассматриваемых правил по их зависимости от данных. Например, результат (1.36), полученный в [37], качественно отличается от всех известных рандомизированных способов наличием зависимости от матрицы с.л.а.у — от исходных данных задачи. А квазиоптимальное правило неявно зависит как от матрицы системы, так и от вектора правой части, в то время как правило циклического обхода зависит только от размерности задачи.
Правила задания проекционной последовательности часто имеют эвристическое обоснование, а их применение оправдано только для узкого класса задач. Показательно в этом отношении правило четной перестановки (1.34), которое было специально разработано для решения задачи восстановления данных по измерениям, выполненным с использованием компьютерного томографа с параллельной схемой сканирования. При этом говорить о зависимости правила от данных в этом случае нельзя, правило построено авторами из интуитивных соображений о структуре самих данных, а не их значений.
В результате проведенного анализа становится ясно, что попытки преодолеть медленную скорость сходимости исходного алгоритма Качмажа с цикли 48 ческим правилом обхода начинают давать результат только тогда, когда вводится зависимость от данных для конкретной решаемой задачи:
Для того чтобы сформулировать программные реализации рассмотренных модифицированных методов Качмажа, необходимо определить правило, по которому вычисления будут прекращены, а полученное приближение будет принято за искомое решение задачи. В данном разделе предполагается, что метод Качмажа применяется только к совместным системам, для которых он гарантированно сходится к единственному решению. Для того чтобы минимизировать влияние регуляризующих свойств правил останова итераций, в данном разделе каждый алгоритм может быть остановлен с использованием двух критериев: где и — точное решение, которое считается известным из процедуры моделирования тестовой с.л.а.у. В качестве модельной задачи для проведения вычислительного эксперимента используется задача компьютерной томографии с параллельной схемой сканирования (Рисунок 1.6) для известного фантома Шеппа-Догана (Рисунок 1.3, справа). Генерация тестовых данных производится с использованием расширения к математическому пакету MatLab из [52]. Формальная схема программной реализации различных модификаций алгоритма Качмажа представлена на Рисунке 1.10.
Общая структура программной реализации алгоритма Качмажа. Красным выделены блоки для которых оценивается количество операций.
Необходимо отметить, что далее при оцешгї сложности исследуемых и предлагаемых модификаций алгоритма Качмажа не учитываются затраты на определение момента прекращения вычислений, но являются важным и сложность микроитерации и сложность вычисления проекционной последовательности.
Как уже отмечалось выше, циклическое правило выбора проекционной последовательности является исторически самым первым и наиболее простым для реализации. Существует два подхода к реализации периодического цикла: 1. Макро-и микроитерации. В алгоритм вводят два дополнительных понятия: макроитерации и микроитерации, при этом одна макроитерация считается выполненной, если в результате т микроитераций каждая строка с.л.а.у. последовательно участвовала в расчетах (как в процедуре (1.4)). Данная реализация приводит к необходимости реализации двух вложенных циклов, внешний из которых ограничивается максимальным количеством макроитераций, а внутренний цикл — последовательный проходит от первой строки с.л.а.у к последней. При этом, что крайне важно, ошибку текущего приближения считают только на макроитерациях; 2. Кольцевой буфер. Данная структура данных часто используется для буферизации данных при обработке сигналов, например, в [53]. Наиболее просто данный подход реализуется через прямое вычисление формулы (1.2). При этом, в отличие от предыдущего подхода, принято рассчитывать ошибку текущего приближения на каждой итерации, так как сама используемая структура данных не предполагает возможности определить начало или конец блока микроитераций. В данной работе используется именно этот подход, так как разделение на макро-и микроитерации в части других модификаций алгоритма Качмажа - идейно невозможно. При оценке количества операций, вычисления, которые необходимо произвести для нормировки строк матриц, не будем принимать в расчет, так как они постоянны, не зависят от количества итераций и являются этапом алгоритма, который преобразует исходные данные задачи в вид, удобный для дальнейших вычислений.
Метод регуляризации А.Н. Тихонова для задачи нормального псевдорешения С.Л.А.У
Задача решения интегрального уравнения Фредгольма первого рода где А [] — оператор Фредгольма, а функция К (s,t) — ядро уравнения, в общем, относится к классу некорректных задач, ее суть состоит в том, чтобы при фиксированных непрерывной функции ядра К (s,t) и функции / (s) найти и (t).
Ж. Адамаром были сформулированы три условия, выполнение которых необходимо, для того чтобы считать задачу корректной (или корректно поставленной)[44]: решение задачи существует — естественное требование, не выполнение которого может свидетельствовать о непригодности математической модели или о неправильной постановки задачи [50]; решение задачи единственно; в некоторых случаях [50], проблема невы 81 полнения данного условия снимается тем, что необходимо найти набор всех допустимых решений, и тогда за решение задачи принимается этот набор; решение задачи непрерывно зависит от входных данных, суть — решение устойчиво к малым возмущениям наблюдаемых данных. Формально данные правила определяются следующим образом в [43]. Рассмотрим уравнение в операторной записи общего вида
Определение 3.1 [43] Задача (3.2) называется корректно поставленной, если выполняются три условия: условие разрешимости: область значений QA оператора А совпадает с F; условие единственности: равенство А [щ] = А [г ] для некоторых щ, U i Є DA влечет равенство щ = и% то есть условие единственности обеспечивается тогда и только тогда, когда оператор А является инъективным [49] — если два образа совпадают, то совпадают и их прообразы; если выполнены первые два условия, то есть существует обратный к А оператор А 1[1], то необходимо выполнение условия устойчивости: обратный оператор А 1 непрерывен на F.
Наиболее часто на практике нарушается последнее условие — условие устойчивости обратного оператора, что приводит к достаточно парадоксальной ситуации [46]: решение задачи невозможно получить обычными численными методами. Ж. Адамар предполагал, что если поставленная задача является некорректной, то она не может встречаться на практике и иметь прикладное значение [44]: «Аналитическая задача всегда корректно поставлена в смысле, указанном ранее, когда есть механическое или физическое истолкование вопроса». С самого начала XX века некорректные задачи рассматривались как математический курьез, не имеющий отношения к приложениям математики [47]. Однако в результате перехода к постиндустриальному способу организации общества появилось преобладающее количество задач, оказавшихся некорректно поставленными, по Ж. Адамару, и связанными с обработкой и интерпретацией данных, полученных в результате различных приближенных измерений, физических экспериментов.
Интенсивное развитие методов решения некорректных задач связано с именами великих советских математиков А.Н. Тихонова, Г.И. Марчука, М.М. Лаврентьева, В.К. Иванова; созданная ими математическая школа позволила решать огромное [43] количество прикладных, задач поставляемых геофизикой, спектроскопией, электронной микроскопией, автоматическим регулированием, теплофизикой, гравиметрией, электродинамикой, оптикой, ядерной физикой, теорией плазмы и другими областями науки и техники.
Первые же работы А.Н. Тихонова в области решения некорректных задач связаны с задачей определения исторического климата Земли и вопросами мерзлоотведения, в результате которых впервые были получены условия, при которых решение задачи восстановления температурного режима планеты оказывается единственным [47]. А.Н. Тихонов в корне изменил подход к решению задач, связанных с обработкой данных неточных измерений, которые до него пытались решать точным образом. Наиболее просто основной тезис А.Н. Тихонова можно сформулировать следующим образом [47]: «Если исходные данные известны приближенно, то и оператор, описывающий процесс, должен быть заменен приближенным, чтобы преобразованная задача оказалась корректной».
Таким образом, А.Н. Тихонову принадлежит обобщение определения корректности задачи по Ж. Адамару, в котором вводится понятие условно кор 83 ректной задачи (или корректной по Тихонову). Основная [43] идея которого заключается в том, чтобы сузить область определения исходного оператора.
Определение 3.2 [43] Задача (3.2) называется корректной по Тихонову или условно корректной, если априори известно, что решение задачи (3.2) существует для некоторого класса данных из F и принадлежит априорно заданному множеству корректности М с DA; решение единственно во множестве М; бесконечно малым изменениям в f (таким, что не выводят решение задачи из множества М), соответствуют бесконечно малые изменения в решении.
Существенному прорыву в теории решения некорректных задач способствовал метод регуляризации, разработанный под руководством А.Н. Тихонова, который успешно был применен для решения различных как фундаментальных, так и прикладных задач [48]: решение интегрального уравнения первого рода, различные обратные задачи теории потенциала и теплопроводности, задача об аналитическом продолжении функции, неустойчивые задачи линейной алгебры, задачи математической экономики и теории оптимального управления, ряд важных обратных задач геофизики, астрофизики, оптической и нейтронной спектроскопии, компьютерной томографии, задачи теории распознавания образов и многие другие. Суть метода регуляризации заключается в стабилизации минимума отклонений Аи от / при помощи сглаживающего функционала Q(u), в результате чего решение задачи (3.2) сводится к минимизации параметрического функционала А.Н. Тихонова.
Приведение обобщенной задачи регуляризации к каноническому виду
Итерации по методу Качмажа рассмотрены неотрывно от исторического контекста развития данного метода. Историческая справка позволяет отводить исследуемому проекционному методу ведущее место при решении прикладных задач из различных областей знаний, в том числе и для анализа данных.
Основной результат, описанный в первой главе, — доказанная автором Теорема 1.3 о геометрической скорости сходимости итераций по методу Качмажа для случая совместных с.л.а.у Необходимо отметить, что, помимо прикладного значения при решении задач большой размерности, доказанная теорема несет и методологическое значение — показывается связь проекционного метода с методом наискорейшего координатного спуска и методом Гаусса-Зейделя. Предложенный алгоритм сравнивается с известными модификациями на примере решения модельной задачи компьютерной томографии и демонстрируется, что, несмотря на высокую сложность одной итерации (Таблица 1.1), количество таких итераций до сходимости уменьшается на порядки (Таблица 1.2).
Отдельное место в данной главе занимает предложенная автором классификация правил задания проекционной последовательности в различных модификациях итераций по методу Качмажа, она предлагается впервые и представляется автору наиболее полной. В первой главе описана историческая хронология развития проекционных методов и указывается место в ней предложенной квазиоптимальной модификации.
Основные результаты второй главы диссертации составляют три новых созданных метода решения задачи регуляризации с использованием проекционного метода: 1. Регуляризованный итерационный проекционный метод — используется правило циклического обхода. 2. Рандомизированный регуляризованный итерационный проекционный метод — используется рандомизированное правило. 3. Квазиоптимальный регуляризованный итерационный проекционный метод — используется квазиоптимальное правило, разработанное в первой главе.
Предложенный в главе квазиоптимальный регуляризованный алгоритм Качмажа обладает крайне быстрой скоростью сходимости, хотя и не удается получить конструктивной оценки скорости сходимости в Теореме 2.4. Очевидно, что параметр q зависит от матрицы системы и от начального приближения, однако записать эту зависимость в явном виде не удалось. Тем не менее для исследуемой задачи компьютерной томографии с параллельной схемой сканирования в присутствии возмущений в векторе правых частей алгоритм сходится быстрее чем рандомизированная модификация сконструированная по методу из [37]. Все алгоритмы из данного раздела опубликованы в работах (в том числе и в совместных) [60, 62, 63, 90].
Важно отметить, что при использовании результатов диссертации и исследований из [10, 11] теоретически решается проблема отыскания нормального псевдорешения для линейных алгебраических систем общего вида, как показано в [10], параметр регуляризации при этом должен удовлетворять условиям UJ2 smach иш mach, где smach — машинное эпсилон. Предложенные в третьей главе алгоритмы опубликованы в виде программного комплекса и его описания в работе [94].
В третьей главе рассмотрена задача решения интегрального уравнения Фредгольма первого рода на основе метода расширенных регуляризованных систем. Для модельной задачи интегральное уравнение приводится к задаче решения с.л.а.у. с использованием квадратурных формул по методу Бубнова-Галеркина. В работе не рассматривается влияние выбранного метода дискретизации на обусловленность получаемой с.л.а.у. В первых разделах приводятся понятия из теории некорректных задач, необходимые для последующего целостного изложения. Приведено важное понятие о корректности вычислительной задачи. Делается вывод о том, что если возмущения в матрице правой части с.л.а.у. отсутсвуют, то задача вычисления нормального псевдорешения всегда корректна — это следует из неравенства (3.9). В некоторых случаях неравенство (3.9) может оказываться слишком пессимистичным, тем не менее возмущения в решении всегда ограничены (для случая 6 (А) = 0), но могут быть крайне велики.
Для большого количества задач диагностики, по результатам косвенных наблюдений, эксперимент в приближенно одних и тех же условиях может быть повторен больше одного раза. Когда это не связано со временем протекания исследуемого явления или критичной важностью объема лучевой нагрузки — тем более. В результате подобные задачи приводят к необходимости решения с.л.а.у. со множеством правых частей. Для этого в Теореме 3.1 предложена оценка возмущений вектора правой части в неравенстве (3.11). Предлагаемый способ тестируется в последнем разделе главы с использованием модельной задачи.