Секция “Краевые задачи механики сплошной среды, численные и численно-аналитические методы решения”
Вид материала | Задача |
- Краевые задачи механики сплошной среды, численные и численно-аналитические методы решения, 51.31kb.
- Секция “Численные и численно-аналитические методы решения краевых задач, 58.61kb.
- Модели и аналитические методы механики сплошной среды направление подготовки, 21.02kb.
- Научная программа конференции. На конференции предполагается работа следующих секций:, 29.46kb.
- Аннотация программы учебной дисциплины «Математические модели механики сплошных сред», 55.95kb.
- Аннотация рабочей программы дисциплины основы механики сплошной среды уровень основной, 22.79kb.
- Секция “Численные методы и пакеты прикладных программ механики конструкций”, 112.71kb.
- Численные методы газовой динамики и теплопереноса, 16.69kb.
- Секция “Краевые задачи в физике и химии твердого тела”, 31.36kb.
- Э. В. Прозорова «Вычислительные методы механики сплошной среды» СпбГУ, 1999, 119.9kb.
Секция “Краевые задачи механики сплошной среды, численные и численно-аналитические методы решения”
УДК 517.958:57, 519.632.8, 519.671
Михайлова А.Г.
Московский государственный университет приборостроения и информатики, г. Москва
О ТРЕБОВАНИЯХ И ВЫБОРЕ МЕТОДА РЕШЕНИЯ ОБРАТНОЙ ЗАДАЧИ ИМПЕДАНСНОЙ ТОМОГРАФИИ
Обсуждаются спектральные свойства матриц систем разрешающих уравнений, к которым приводится обратная задача импедансной томографии. Сформулированы рекомендации к выбору методов решения.
Метод электрической импедансной томографии. Электрическая импедансная томография – это метод получения трехмерной структуры биологического объекта. При приложении к биологическому объекту переменного тока, на поверхности регистрируются вызванные потенциалы. Зная величины приложенных токов и регистрируемых потенциалов, находят пространственное распределение электрической проводимости внутри биологического объекта. Поскольку различные участки биологического объекта, в зависимости от вида ткани и ее физиологического состояния, имеют различные свойства (проводимость, диэлектрическая проницаемость), знание пространственного распределения этих свойств позволяет восстановить внутреннюю структуру этого биологического объекта.
Задача восстановления трехмерной структуры по области Ω сводится к решению уравнения
,
где – точка в области Ω,
φ(x, w) – регистрируемый потенциал на поверхности объекта,
– комплексная электрическая проводимость,
- электрическая проводимость биологической ткани,
- диэлектрическая проницаемость,
– угловая частота приложенного тока.
Учет отсутствия источников и приемников электрического тока, дискретности электродов и сверхпроводимости материала электродов, а также электрохимических эффектов на границе «электрод – кожа» приводит к полной модели:
на ∂Ω;
;
на интервалах между электродами;
на el, ;
.
Система уравнений. Наличие достаточно большого количества экспериментов, в каждом из которых найдено распределение тока через границу Σ, соответствующее распределению потенциала на этой же границе, в общем случае, формально означает, что известно отображение (оператор): , определяемый электрическими свойствами биоткани. Таким образом, для численного решения задачи получаем систему уравнений вида , где - вектор свойств биоткани в каждом узле дискретизации, - данные измерений, - оператор системы, получаемый, например, методом конечных элементов (МКЭ).
Матрица системы. При моделировании биологических объектов сложной формы число узлов дискретизации и количество образуемых на них элементов может достигать больших значений (например, Байфорд с коллегами [1] для моделирования головы человека использует – 23219 элементов – для кожи головы, 63982 элементов – для черепа, 47140 элементов – для черепно-мозговой жидкости, 21574 элемента – для мозга. Итого, в сумме – более 155 тысяч элементов). Матрица, полученная МКЭ, является сильно разреженной, что дает возможность работать с большим количеством узлов дискретизации, используя экономный формат хранения разреженных матриц. Типичная матрица ЭИТ показана на рисунке 1.
Рисунок 1 - Типичный характер разреженности матрицы системы для задачи импедансной томографии
Использование классических прямых методов решения системы уравнений (методы Гаусса и Жордана-Гаусса, метод Крамера, матричный метод и др.) требует больших вычислительных ресурсов (особенно при решении трехмерных задач) и, таким образом, для случая импедансной томографии являются неподходящими. Также их нельзя использовать в случае сингулярности матрицы системы, что имеет место в импедансной томографии.
Методы подпространства Крылова. Вероятно, самым используемым семейством методов для решения больших и разреженных систем, относящихся к группе нестационарных методов, является семейство методов подпространства Крылова. Методы подпространства Крылова имеют два свойства, которые делают их особенно эффективными для решения систем больших размерностей. Первым достоинством является низкая требовательность к ресурсам памяти за счет использования для ортогонализации метода Ланцоша. Вторым достоинством является то, что исходная матрица , а иногда и ее сопряженная матрица , необходимы только в форме векторного произведения. Важность этого свойства заключается в том, что явное формирование матрицы не требуется, пока имеется возможность получения векторного произведения.
К основным методам подпространства Крылова относят метод сопряженных градиентов и метод двусопряженых градиентов, метод минимальных невязок, обобщенный метод минимальных невязок и метод квазиминимальных невязок, а также LSQR–метод и симметричный LQ–метод. Выбор того или иного метода решения зависит, прежде всего, от спектральных характеристик матрицы – от числа обусловленности, спектрального радиуса и собственно собственных значений матрицы.
По своей природе задача импедансной томографии является плохо обусловленной, и число обусловленности может составлять от нескольких десятков до нескольких тысяч в зависимости от неидеальности формы представляемого биологического объекта. При этом спектральный радиус всегда будет значительно больше единицы. Таким образом, матрица является сингулярной, хоть и положительно определенной, и требуется использовать предобуславливатели.
Например, для простого двумерного случая формы (рисунок 2) и его дискретизации (рисунок 3), состоящей из 41 узла, 68 элементов (треугольников), матрица системы является разреженной (рисунок 4) при моделировании объекта как однородного с вектором переменных =1 получено решение (рисунок 5) со значениями невязок (таблица 1).
Рисунок 2 - Простой двухмерный объект
Рисунок 3 - Дискретизация двухмерного объекта
Рисунок 4 - Визуализация разреженности матрицы. Степень заполненности - 0.15
Рисунок 5 - Решение со значениями невязок из таблицы 1
Таблица 1 – Невязки
Метод | Относительная невязка |
Метод сопряженных градиентов (CG) | 3.4595e-006 |
Метод двусопряженных градиентов (Bi-CG) | 4.2744e-004 |
Стабилизированный метод двусопряженных градиентов (Bi-CGStab) | 4.2744e-004 |
Метод минимальных невязок (MINRES) | 3.7045e-003 |
Обобщенный метод минимальных невязок (GMRES) | 4.4671e-004 |
Метод квазиминимальных невязок (QMR) | 3.5576e-004 |
LSQR–метод (LSQR) | 1.1750e-003 |
Симметричный LQ–метод (symmLQ) | 1.2352e-002 |
Свойства биологического объекта через свойства матрицы. Возможны два случая: 1) объект рассматривается как обладающий только резистивными свойствами – тогда матрица Якоби, а также вектор неизвестных(электрическая проводимость) являются действительными положительными величинами, 2) объект рассматривается как обладающий и резистивными, и диэлектрическими свойствами – матрица системы может иметь отрицательные и комплексные собственные значения, матрица Якоби может имеет действительные и мнимые компоненты, которые можно объединить по схеме , а вектор неизвестных принимает вид , где - частота прикладываемого воздействия.
Нелинейность системы уравнений. В классическом варианте методы подпространства Крылова используются для решения системы линейных уравнений, но они также могут быть использованы для решения и систем нелинейных уравнений, как в случае импедансной томографии, при предварительной процедуре линеаризации, заключающейся (например, для метода Ньютона) в том, что для нелинейной дифференцируемой системы уравнений вида , при начальном приближении , последующая итерация определяется как , . Матрица является матрицей Якоби и определяется для как .
Физическая природа процессов, описываемых конкретной матрицей Якоби, требует наличия отрицательных действительных частей для поддержания стабильности биологической системы, что также способствует ухудшению свойств матрицы системы (более подробно о матрице Якоби применительно к задаче импедансной томографии, например, [2]).
Выясненные свойства матрицы чувствительности позволяют обосновать выбор численного метода решения обратной задачи.
СПИСОК ЛИТЕРАТУРЫ
1. Bayford R.H., Gibson A., Tizzard A., Tidswell T., Holder D.S. “Solving the forward problem in electrical impedance tomography for the human head using IDEAS (integrated design engineering analysis software). A finite element model tool”, Physiol. Meas., #22, 2001. – pp. 55-64.
2. Михайлова А.Г. Физический смысл и построение матрицы чувствительности (Якобиан) для решения обратной задачи импедансной томографии. / А.Г. Михайлова, А.С. Кравчук // Труды XXI конференции «Математические методы в технике и технологии ММТТ-21», Саратов, май 2008. – с. 76-81.