99703

Программа реализующяа решение плохо обусловленных систем методом регуляции А.Н. Тихонова

Курсовая

Информатика, кибернетика и программирование

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

Русский

2016-10-08

117.26 KB

7 чел.

Содержание

  1.  Постановка задачи………………………………………………………….…..1
  2.  Теоретические сведения……………………………………………………….2
  3.  Реализация алгоритма Тихонова……………………………………………...8
  4.  Руководство пользователя……………………………………………………15
  5.  Вывод…………………………………………………………………………..15

Постановка задачи.

Написать программу реализующую решение плохо обусловленных систем методом регуляции А.Н.Тихонова

Теоретические сведения.

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

Пусть имеем систему линейных алгебраических уравнений  Ax=b (1.1). Предположим, что матрице А задана, является не вырожденной и плохо обусловленной.  В  этом  случае  решение  системы (1.1) существует, единственно  и  неустойчиво,  т.е.  задача  некорректна.  Некорректно поставленные  задачи  создают  значительные  трудности  при  их приближенном решении. Как правило, исходные данные прикладных задач известны с некоторой погрешностью. Отыскание неустойчивых решений прикладных  задач  приводит  к  тому,  что  классические  методы вычислительной  математики  могут дать  результат,  сколь  угодно  сильно отличающийся  от точного решения,  при условии, что  исходные данные отличаются сколь угодно мало.

Рассмотри простейший случай, когда матрица А - симметричная матрица. Пусть  –  ее собственные значения, упорядоченные в порядке убывания . Соответствующую ортогональную систему собственных векторов обозначим через   . Решением  системы Аx=b является вектор

X=, где .

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

.

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

Рассмотрим способы, иллюстрированные на примерах матриц с >0.

Первый способ.

Задавая некоторое α>0, находим решение xα системы

(A+ αΕ)xα =.

Оно записывается в виде xα = . Так как   , то наличие малого параметра  не существенно изменит слагаемое, соответствующие  большим . В тоже время при  имеем <<. это означает что введение параметра   приводит к существенному уменьшению роли слагаемых, соответствующим малым . Подбор оптимального значения  обычно осуществляется экспериментально, сравнивая результаты расчетов при различных .

Второй способ.

Заключается в следующем. Будем решать систему уравнений каким-то итерационным способом. Рассмотрим случай итерации по формуле

(*)

При некотором начальном приближении x0. Пусть

= ,  

Подставим эти выражения в (*) и, приравнивая коэффициенты  про ei , получим соотношение .

Последовательно подставляя  через предыдущие, имеем  =.

Если <1 , то при n 

 

Поэтому . Пусть  При больших значениях      быстро стремится к 0 с ростом n и   близко к своему предельному значению.

В других случаях решение задачи находят, минимизируя некоторый функционал, близкий к F(x)=(Ax,x)-2(b,x), например функционал F(x)+ (x,x)  с малым

Успешность  применения описанных приемов в случае несимметричных матриц А в существенной степени  зависит от  структуры жордановой формы и от ряда свойств матрицы. Здесь часто решение находят минимизируя функционал (Аx-b, Аx-b)+(x,x), при малых >0.

Метод А.Н.Тихонова применим к задаче (1.1)  в ее общей постановке (без условия самосопряженности).   

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

Во многих прикладных задачах таким условием является требование, чтобы длина вектора решения х была минимальной среди всех решений этой системы, т.е. чтобы скалярное произведение (х,х) было минимальным. Такое решение системы (1.1) называется нормальным и обозначается xn . Вычисление на ЭВМ минора матрицы А, равного нулю, обычно приводит к результату отличному, от нуля даже при малых ошибках округления. В такой ситуации совершенно не ясно, какой минор нужно считать базисным минором  матрицы А. Алгоритм становится некорректным.

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

x-b, Аx-b) = 0. (1.2)  

Регуляризованная  задача  формулируется  так:

Найти  решение  xα, зависящее от α , при котором функционал Fα(x) =(Ax-b, Ax-b)+α(x,x)  (1.3), где  α>0 , является  минимальным,  и  выбрать  параметр  α    так,  чтобы регуляризованное  значение хα  было  наиболее  близким  к  нормальному решению хn системы (1.1) . Регуляризирующий функционал Fα(x) запишем, вычисляя скалярное произведение в ортонормированном базисе:  

Fα(x) =     (1.4)  

Должно выполняться необходимое условие минимума функции:

Fα(x)=Fα(…, , поэтому имеем

=, j=   (1.5)

Равенство (1.5) в матричной форме имеет вид:

(ATA+ αΕ)xα =AT b,α>0     (1.6)

Таким  образом,  решение  неопределенной  системы (1.1) сведено  к решению однозначно разрешимой системы (1.6). Решение этой  системы можно  найти,  используя  прямые  или  итерационные  метода  решения корректных задач. Таким образом, простейший вариант метода А.Н. Тихонова заключается в отыскании точки минимума функционала  

Fα(x)= α

x      (1.7)

где α - малый положительный параметр.

Эта задача на минимум равносильно определению x из уравнения   

?????????αx+Α*Ax=A*b.

Пример.

Для реализации численного примера был выбран метод Тихонова решения плохо обусловленных СЛАУ.  В качестве исходной была взята СЛАУ Ax=b, имеющая в матричной записи вид:

Определитель матрицы коэффициентов этой системы близок к нулю – он равен 0.000125. Попробуем решить эту систему с помощью обратной матрицы:

x=A-1b

Получим    x1=316

                  x2=-990

                  x3=832

Теперь предположим, что правая часть нам известна приближенно, с погрешностью 0.1  Изменим, к примеру, третий элемент вектора-столбца с 1 на 1.1 :

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

               x1=348

       x2=-1090

       x3=916.

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

Будем искать решение методом Тихонова. В теоретической части было показано, что целесообразно использовать регуляризирующий оператор следующего вида:   (ATA+ αΕ)xα =AT b, где E – единичная матрица, xα  – приближенное нормальное решение, AT – транспонированная исходная матрица, α – параметр регуляризации,       b – правая часть, заданная неточно. Эту задачу можно решать стандартными методами, задав предварительно функцию α = α (b), удовлетворяющую условиям теоремы Тихонова.

Приближение к нормальному решению

X(1)= 3.47

X(2)=-1.08

X(3)= 9.15

Значение правой части при подстановке прибл. решения

B1(1)= 9.99

B1(2)= 1.00

B1(3)= 1.09

Реализация алгоритма Тихонова.

Создаем транспонированную матрицу _At матрицы _A[j][i], которая, в свою очередь является главной матрицей программы

void _TRANS_A(){                 

     for(int i=0;i<size;i++)

   for(int j=0;j<size;j++)

    _At[i][j]= _A[j][i];

}

//-----------------------------------------------------

Создаем единичную матрицу   _ONE[i][j]

for(int i=0;i<size - 1;i++)

 BB[i] = new double [size];

   for(int i=0;i<size;i++)

 for(int j=0;j<size;j++)

  if(i==j)

    _ONE[i][j] = 1;

   else

                _ONE[i][j] = 0;

//-----------------------------------------------------

Определяем знак коэффициента при малом определителе

double _COEFF(int _FPARAM,int _SPARAM){

  int _SUMM = _FPARAM + _SPARAM;

  double _TMP = 1.0;

     for(int k=0;k<_SUMM;k++)

   _TMP*=(-1);

  return _TMP;

}

//-----------------------------------------------------

считает определитель на 1 размер меньше исходной матрицы

double _ALL_ADD(int _FPARAM,int _SPARAM){

         double **TMP = new double*[size - 1];

    for(int i=0;i<size-1;i++)

   TMP[i] = new double [size - 1];

    for(int i=0;i<size - 1;i++)

     for(int j=0;j<size - 1;j++)

      TMP[i][j] = 0;

       for(int v=0;v<size;v++)

     for(int w=0;w<size;w++)

      {

        AA[v][w]._VALUE = _A2[v][w];

        AA[v][w]._FLAG = true;

      }

     for(int v=0;v<size;v++)

      {

              AA[_FPARAM][v]._FLAG = false;//фолсим строку 

              AA[v][_SPARAM]._FLAG = false;//фолсим столбец

      }

    k = 0;

    int x = 0,y=0;

     {

    for(int v=0;v<size;v++)

     for(int w=0;w<size;w++)

      {

      if(AA[v][w]._FLAG)

       {

          TMP[x][y] = AA[v][w]._VALUE;

          y++;

           if(y == size - 1)

         {

           x++;

           y = 0;

         }

        if(x == size - 1) {break;break;}

       }

      }

     }

   

    double res = DET(TMP,size - 1);

    delete [] TMP;

    return res;

    

}

//------------------------------------------------------------

определитель матрицы А2; где А2==АТ *А+α*Е (для разложения используем первую строчку)

double DET_COUNT(){

 double S1;

  for(int z=0;z<size;z++)

   S1+=_A2[0][z]*_COEFF(0,z)*_ALL_ADD(0,z);

 return S1;

 }

//-------------------------------------------------------------

Находим обратную матрицу А2==_Ar[j][i]

void _REV_MATR(){

 for(int i=0;i<size;i++)

 for(int j=0;j<size;j++)

  _Ar[j][i] = (double)(_COEFF(i,j)*_ALL_ADD(i,j)/DET_COUNT());

 }

//-------------------------------------------------------------

Определяем разницу между правой и левой частью уравнения по абсолютному значению и сравниваем с ε, т.е. сравниваем |АТ*В- (АТ *А+α*Е)*Хα|<ε, если условие выполниется, то заканчиваем цикл

bool _ALL_RIGHT(){   

 for(int i=0;i<size;i++)

 {

    if(abs(RIGHT[i] - _U1[i]) > precision)

    return false;

 }

 return true;

}

//-------------------------------------------------------------

Функция, результатом которой является _FPARAM в степени _SPARAM

double _POW(double _FPARAM,int _SPARAM){

  double S2 = 1;

  int z;

       if(_SPARAM==0)

  return 1;

 else

        for(z=0;z<_SPARAM;z++)

   S2*=_FPARAM;

   return S2;

}

//------------------------------------------------------------

void MAIN(){

 обозначаем, что _A[i][j] главная  матрица

       for(int i = 0;i<size; i++ )

    for(int j = 0;j<size;j++)

_A[i][j] = MAIN_MATRIX[i][j];

транспонируем главную матрицу

   _TRANS_A();      

                                            

транспонированная матрицу _At[i][j] умножаем на правый столбец (получаем правую часть уравнения Тихонова==_B[i])

   for(int i=0;i<size;i++){      

     _S = 0;                    

      for(int j=0;j<size;j++)     

    {

        _S+=_At[i][j]*RIGHT[j];

        _B[i] = _S;

     

       }

    }

Транспонированную главную матрицу умножаем на главную матрицу получем матрицу_A1[i][j], т.е. _A1[i][j]==АТ

        for(int i=0;i<size;i++)

   for(int j=0;j<size;j++)

    {

    _S = 0;

    for(int k=0;k<size;k++){

     _S+= _At[i][k]*_A[k][j];

         _A1[i][j] = _S;

     }

    }

     Сначала количесвто итераций равно 0

             q = 0;

    do{

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

     this->_ALPHA = (double)(q / _POW(4,q));

    

находим матрицу равную произведению α и единичной матрицы, т.е. _One1[i][j]==α*Е

        for(int i=0;i<size;i++)

     for(int j=0;j<size;j++)

      _One1[i][j] = _ONE[i][j]*_ALPHA;

находим матрицу  _A2[i][j]==(АТ *А+α*Е)

    for(int i=0;i<size;i++)

     for(int j=0;j<size;j++)

      _A2[i][j] = _One1[i][j]+_A1[i][j];

находим обратную матрицу _A2[i][j]== +=_Ar[i][j] ==(АТ *А+α*Е)-1

    this->_REV_MATR();

находим результирующий столбец _Z[i]== (АТ *А+α*Е)-1Т

    for(int i=0;i<size;i++){

     _S = 0;

     for(int j=0;j<size;j++){

        _S+=_Ar[i][j]*_B[j];

        _Z[i] = _S;

      }

     }

 

находим матрицу произведения главной матрицы и результирующего столбца _U1[i], т.е. _U1[i] == (АТ *А+α*Е)*Хα    

    for(int i=0;i<size;i++){

     _S = 0;

     for(int j=0;j<size;j++){

      _S+=_A[i][j]*_Z[j];

          _U1[i] = _S;

         }

     }

    q++;

     }

пока не выполнится условие |АТ*В- (АТ *А+α*Е)*Хα|<ε будем  продолжать цикл, увеличивая на каждом шаге k на единицу

    while(!_ALL_RIGHT());

 

}

Руководство пользователя

При запуске программы появится окно.

9

8

7

6

5

4

3

2

1

  1.  Исходная матрица(система линейных уравнений).
  2.  Решение системы линейных уравнений.
  3.  Векторный столбец правой части уравнения Ax=b.
  4.  Порядок матрицы.
  5.  Интерактивно задаваемый диапазон элементов.
  6.  Интерактивно задаваемый диапазон определителя.
  7.  Загрузка  матрицы из файла.
  8.  Данные о программе.
  9.  Интерактивно задаваемая точность.

Вывод.

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


 

А также другие работы, которые могут Вас заинтересовать

11773. Решение задачи целочисленного ЛП с помощью динамического программирования 481.5 KB
  Курсовая работа по дисциплине МАТЕМАТИЧЕСКИЕ МЕТОДЫ ИССЛЕДОВАНИЯ ОПЕРАЦИЙ на тему Решение задачи целочисленного ЛП с помощью динамического программирования АННОТАЦИЯ Курсовая работа содержит 40 страниц 8 формул 17 таблиц 10 литературных источников. В ...
11774. Розв’язання несумісних задач лінійної оптимізації в діалоговому режимі 199.06 KB
  Звіт до лабораторної роботи №4 на тему: Розвязання несумісних задач лінійної оптимізації в діалоговому режимі З курсу: Математичні методи дослідження операцій Мета: Знайти оптимальний розвязок розподілу ресурсів за умови первинної наявності несумісн
11775. Предмет та задачі дослідження операцій 165 KB
  Предмет та задачі дослідження операцій Предмет та історія виникнення дослідження операцій ДО. Основні поняття ДО та етапи операційного дослідження Пряма та обернена задачі ДО. Детерміновані задачі ДО. Проблема вибору розвязків в умовах невизначеності. ...
11776. Пошук оптимального розв’язку багатокритерійних лінійних задач 153.21 KB
  Звіт до лабораторної роботи №5 на тему: Пошук оптимального розвязку багатокритерійних лінійних задач З курсу: Математичні методи дослідження операцій Мета: Вивчити методологію розвязання багатокритерійних оптимізаційних задач на прикладі задачі розпо...
11777. Аналіз оптимізаційних задач. Аналіз оптимального розв’язку. Параметричний аналіз. Графічне представлення 139.27 KB
  З в і т про виконання лабораторної роботи 3 на тему: Аналіз оптимізаційних задач. Аналіз оптимального розвязку. Параметричний аналіз. Графічне представлення. З курсу: Математичні методи дослідження операцій Мета: Вивчити методи аналізу задач лінійного пр
11778. Основні етапи розв’язання задач на дослідження операцій 77.35 KB
  Звіт до лабораторної роботи №1 з дисципліни Математичні методи дослідження операцій ЛАБОРАТОРНА РОБОТА №1 Тема роботи: Основні етапи розвязання задач на дослідження операцій. Питання: Поняття про задачі дослідження операцій. Поняття моделі класифік...
11779. ТРАНСПОРТНАЯ ЗАДАЧА ЛИНЕЙНОГО ПРОГРАММИРОВАНИЯ (ТЗЛП) 256 KB
  ТЕМА 1 ТРАНСПОРТНАЯ ЗАДАЧА ЛИНЕЙНОГО ПРОГРАММИРОВАНИЯ ТЗЛП 1. Содержательная постановка и формальная модель ТЗЛП. 2. Условие существования решения ТЗЛП. Построение формальной модели ТЗЛП при нарушении усло...
11780. Математичні методи дослідження операцій. Лабораторні роботи 465 KB
  ЛАБОРАТОРНІ РОБОТИ 16 з дисципліни Математичні методи дослідження операцій ЛАБОРАТОРНА РОБОТА № 1 ТЕМА: Графічний метод розвязування задач лінійного програмування ЛП. МЕТА: Опанувати графічний метод побудови допустимої області ЗЛП та ліній рівн
11781. Первоначальные сведения и правила работы в операционной системы Windows’2000/2003 59.5 KB
  Лабораторная работа № 1. Тема: Первоначальные сведения и правила работы в операционной системы Windows2000/2003. Цель работы: Освоить основные приемы и методы работы с объектами операционной системы Windows2000/2003 XP. Содержание работы: Освоение основных приемов настройки и ра...