VPF::[C++] Решение СЛАУ меотодом Зейделя

Ardour

Дата 13.6.2007, 05:24 (ссылка)

(нет голосов)

Загрузка ...

Новичок

Профиль
Группа: Участник
Сообщений: 5
Регистрация: 10.6.2007

Репутация: нет
Всего: нет

Дана система уравнений:
0.07596x1-0.02041x5=6;
0.03741x2-0.01176x3-10x4=-256;
-0.01176x2+5.07059x3-0.05882x4=-425;
-10x2-0.05882x3+10.14216x4=-241;
-0.02041+0.13152x5=-9.

При это после каждой итерации нужно выводить значение на экран. Погрешность 0.001.

Kuvaldis

Дата 13.6.2007, 08:23 (ссылка)

(нет голосов)

Загрузка ...

механик-вредитель

Профиль
Группа: Участник Клуба
Сообщений: 1189
Регистрация: 16.6.2006
Где: Минск

Репутация: 32
Всего: 61

Ardour,
Вот кусок из моей лабы. Разбирайся, вбивай свои данные

Код


//-----------------------------------------------------------------------
#include <iostream>
#include <stdlib.h>
#include <conio.h>
#include "Iter.h" 
//-----------------------------------------------------------------------
using namespace std;
//-----------------------------------------------------------------------
void Ziedel()
{
    fType** A;            // матрица системы
    fType*  b;            // правая часть
    fType*  r;            // вектор невязки
    fType*  r_old;        // старый вектор невязки
    fType*  x_old;        // старое решение 
    fType*    x_new;        // новое (уточненное) решение
    fType    sum1;        
    fType    sum2;
    fType    eps;        // требуемая погрешность
    fType    curr_eps;    // текущая погрешность
    FILE*    out;        // выходной файл
    int        iter;        // количество итераций
    int        n;            // порядок матрицы системы
    int     i, j;       // переменные цикла  

    // бесконечный цикл

    while (1)
    {

        system("cls");
        out = fopen("ZiedelOut.txt","w");

        cout << "Matrix order = ";
        cin >> n;
        fprintf(out,"Matrix order = %d\n", n);
        
        cout << "Precision = "; 
        cin >> eps;
        cout << endl;
        fprintf(out, "Precision = %0.*f \n\n",SHARP, eps);

        // инизиализация переменных для вычисления
        A     = CreateMatrix(n);
        b     = CreateVector(n);
        x_old = CreateVector(n);
        x_new = CreateVector(n);
        r     = CreateVector(n);
        r_old = CreateVector(n);
        iter = 1;
        
        
        while (iter <= MAX_ITER)
        { 

            // вычисление вектора x_new

            for (i = 0; i < n; i++)
            { 
                sum1 = 0;
                for (j = 0; j < i; j++)            
                {                                        
                    sum1 += A[i][j] * x_new[j];    
                }
                sum2 = 0;
                for (j = i + 1; j < n; j++)            
                {
                    sum2 += A[i][j] * x_old[j];    
                }
                x_new[i] = (b[i] - sum1 - sum2) / A[i][i];        
            }


            // вычисление вектора невязки r

            for (i = 0; i < n; i++)                        
            {                                        
                sum1 = 0;                                
                for (j = 0; j < n; j++)    
                {
                    sum1 += A[i][j] * x_new[j];        
                }
                r[i] = sum1 - b[i];                        
            }

            // получение и вывод текущей невязки
            sum1 = VectorNorm(r, n);
            printf("Iteration %*d: \nNevyazka norm = %0.*f \n", WIDTH_NUM, iter, OUT_SHARP, sum1);
            fprintf(out, "Iteration %*d: \nNevyazka norm = %0.*f \n", WIDTH_NUM, iter, OUT_SHARP, sum1);
            
            // определение отношения невязок
            sum2 = VectorNorm(r_old, n);
            sum2 = sum1 / sum2;
            printf("Nevyazka attitude = %0.*f \n", OUT_SHARP, sum2);
            fprintf(out, "Nevyazka attitude = %0.*f \n", OUT_SHARP, sum2);

            // вычисление текущей погрешности решения
            sum1 = VectorNorm(x_new, n);
            sum2 = VectorNorm(x_old, n);
            
            // выбираем максимальную норму для оценки: sum1 = max(sum1, sum2)
            if ( sum2 > sum1 )                                    
                sum1 = sum2;                                        
                                                                
            sum2 = GetDeltaNorm(x_new, x_old, n); 
            curr_eps = sum2 / sum1;               
            printf( "eps = %0.*f \n\n", OUT_SHARP, curr_eps );
            fprintf( out, "eps = %0.*f \n\n", OUT_SHARP, curr_eps );

            // выйти, если достигнута точность
            if (curr_eps <= eps)    
                break;            

            // переход к следующей итерации
            iter++;
            CopyVector(x_old, x_new, n);
            CopyVector(r_old, r, n);
        }

        // вывод результатов
        if (curr_eps <= eps)
        {
            // вывод результатов
            printf( "\nSolving vector: \n" );        
            fprintf( out, "\nSolving vector: \n" );

            for (i = 0; i < n; i++)
            { 
                printf( "X[%d] = %*.5f \n", i, OUT_SHARP, x_new[i] );        
                fprintf( out,"%0.*f \n", OUT_SHARP, x_new[i] );
            }
            fprintf(out,"\n");
        }
        else
        {
            cout << "Unable to find desicion" << endl;
        }

        // освобождение ресурсов
        FreeMatrix(A, n);
        delete [] b;
        delete [] x_old;
        delete [] x_new;
        delete [] r;
        fclose(out);
        
        cout << "ESC to quit" << endl;
        char sym = getch();
        if ( sym == 0x1B )            // выйти по ESC 
            break;                    
    }

    return;
}
//-----------------------------------------------------------------------

--------------------

Помни - когда ты спишь, враг не дремлет
Спи чаще и дольше, изматывай врага бессоницей

Weightlessness

Дата 13.9.2007, 15:40 (ссылка)

(нет голосов)

Загрузка ...

Новичок

Профиль
Группа: Участник
Сообщений: 7
Регистрация: 9.9.2007
Где: Antarctica

Репутация: нет
Всего: нет

Не могу найти ошибку в коде. Кажется неправильно преобразуется матрица. Компиллятор - Turbo C++ 3.0. Синтаксис нужен примерно такой же. Утром сдавать!

Код


#include <stdio.h>
#include <conio.h>
#include <math.h>

void main()
{

   float a[5][5]={                                                 //матрица коэффициентов
                 {0.78823,0,0,0,-0.21177},
                 {0,0.50067,-0.000586588,-0.4988,0},
                 {0,-0.0022874,0.98628,-0.011441,0},
                 {0,-0.495,-0.0029116,0.50204,0},
                 {-0.13434,0,0,0,0.86566},
                 };
   float d[5]={62.2536,12.17072,82.66675,-11.9295,-59.23782};  //правая часть

   float x[5]={0};                         //результат
   float b[5][5]={0};
   float g[5]={0};
   float x1[5];

   int i,j;
   int k;
   int h=0;

   float v=0; 

   float xk[5]={0};

   float e=0.001;    //погрешность

   k=0;

   clrscr();

   printf("Исходная матрица:\n\n");
   for(i=0;i<5;i++)
   {
      for(j=0;j<5;j++)
     printf("%5.5f   ", a[i][j]);
      printf("\n");
   }
   printf("\nПравая часть:\n\n");
   for(i=0; i<5; i++)
      printf("%5.5f  ", d[i]);
   printf("\n");

   for(j=0;j<5;j++)
   {
      g[j]=d[j];
      x[j]=g[j];
      x1[j]=x[j];
      for(i=0;i<5;i++)
      {
     if(i==j)   b[j][i]=1-a[j][i];

     else       b[j][i]=-a[j][i];
      }
   }

   printf("\nПреобразованная матрица:\n\n");
   for(i=0; i<5; i++)
   {
       for(j=0; j<5; j++)
      printf("%5.5f  ", b[i][j]);
       printf("\n");
   }
   printf("\n");

   do
   {
      k++;

      for(j=0;j<5;j++)
      {
     xk[j]=g[j];
     for(i=0;i<5;i++)
     {
        xk[j]+=b[j][i]*x1[i];
     }
     x1[j]=xk[j];
      }

      for(j=0;j<5;j++)
      {
          v=x[j]-xk[j];
          v=fabs(v);
          if(v<=e)
          {
               h++;
          }
      }

      if(h<5)
      {
     h=0;
     for(j=0;j<5;j++)
        x[j]=xk[j];
      }
   }
   while(h<5);

   printf("\Число итераций:\t");
   printf("%d\n\n", k);

   //Результат
   printf("Результат:\n\n");
   for(i=0;i<5;i++)
   {
       printf("%5.5f   ", x[i]);
   }


   getch();

}

Weightlessness

Дата 13.9.2007, 17:07 (ссылка)

(нет голосов)

Загрузка ...

Новичок

Профиль
Группа: Участник
Сообщений: 7
Регистрация: 9.9.2007
Где: Antarctica

Репутация: нет
Всего: нет

И еще вопрос: к какому виду нужно привести матрицу, чтобы добиться сходимости? Насколько я понял, сумма модулей элементов в правой части должна быть равна 1. Т.е. Нужно в каждой строке в правой части сложить модули элементов (e.g. sum), затем найти некоторую величину vel=1/sum, а потом умножить элементы каждой строки (слева и справо) на vel?

Kuvaldis

Дата 14.9.2007, 03:23 (ссылка)

(нет голосов)

Загрузка ...

механик-вредитель

Профиль
Группа: Участник Клуба
Сообщений: 1189
Регистрация: 16.6.2006
Где: Минск

Репутация: 32
Всего: 61

Weightlessness,

Цитата

И еще вопрос: к какому виду нужно привести матрицу, чтобы добиться сходимости?

Критерий сходимости достаточно сложен для анализа, поэтому на практике пользуются достаточным условием: диагональное преобладение матрицы, т.е.
модуль каждого диагонального эл-та больше суммы модулей ВСЕХ эл-тов ряда, в котором находится данный диагональный эл-т

--------------------

Помни - когда ты спишь, враг не дремлет
Спи чаще и дольше, изматывай врага бессоницей

Weightlessness

Дата 20.9.2007, 10:15 (ссылка)

(нет голосов)

Загрузка ...

Новичок

Профиль
Группа: Участник
Сообщений: 7
Регистрация: 9.9.2007
Где: Antarctica

Репутация: нет
Всего: нет

Решение задачи:

Код


#include <stdio.h>
#include <conio.h>
#include <math.h>

void main()
{

   float a[5][5]={                             //левая часть СЛАУ
         {0.945,-0.5882,0,0,0},
         {-0.017,0.9806,-0.0114,0,0},
         {0,-0.0153,0.9876,-0.0107,0},
         {0,0,-0.1111,0.9188,-0.5},
         {0,0,0,-0.75,0.9378},
         };
   float d[5]={50,-15.2371,-14.2308,0,-48};    //правая часть СЛАУ
   float x[5]={0};                             //массив неизвестных значений
   float b[5][5]={0};
   float g[5]={0};
   float x1[5]={0};
   int i,j;
   int k;
   int y=0;
   float z;
   float xk[5]={0};
   float e;                                  //погрешность
   clrscr();

   k=0;

   printf("Исходная матрица, являющаяся правой частью СЛАУ:\n\n");
   for(i=0;i<5;i++)
   {
      for(j=0;j<5;j++)
     printf("%.4f   ", a[i][j]);
      printf("\n");
   }

   printf("\nПравая часть СЛАУ:\n\n");
   for(i=0; i<5; i++)
      printf("%.4f  ", d[i]);
   printf("\n\n");

   printf("Введите требуемую погрешность нахождения неизвестных:");
   scanf("%f", &e);
   printf("\t%.5f\n\n", e);

   //Первоначальные преобразования
   for(j=0;j<5;j++)
   {
      g[j]=d[j];
      x[j]=g[j];
      x1[j]=x[j];
      for(i=0;i<5;i++)
      {
     if(i==j)   b[j][i]=1-a[j][i];

     else       b[j][i]=-a[j][i];
      }
   }

   printf("Преобразованная матрица:\n\n");
   for(i=0; i<5; i++)
   {
      for(j=0; j<5; j++)
         printf("%.4f   ", b[i][j]);
      printf("\n");
   }
   printf("\n");

   //Организация итерационного процесса (метод Зейделя)
   do
   {
      y=0;
      k=k+1;
      printf("Итерация номер %d\t\n", k);

      for(j=0;j<5;j++)
      {
     xk[j]=g[j];
     for(i=0;i<5;i++)
     {
        xk[j]+=b[j][i]*x1[i];
     }
     x1[j]=xk[j];
      }

      //проверка соответствия значений неизвестных заданной точности
      for(j=0; j<5; j++)
      {
         
     z=fabs(xk[j]-x[j]);
     printf("x[%d]: %.4f, xk[%d]: %.4f\n", j, x[j], j, xk[j]);
         printf("|x[%d]-xk[%d]|=%.4f", j, j, z);
     x[j]=xk[j];
     if(z<e)
     {
        y=y+1;
     }
         printf("\n");
      }
      printf("\n");
   }
   while(y!=5);

   printf("Число итераций:\t");
   printf("%d\n\n", k);

   //Printing the result
   printf("Результат, найденный с погрешностью %.5f:\n\n", e);
   for(i=0;i<5;i++)
   {
       printf("%.4f   ", x[i]);
   }


   getch();
}

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

0 Пользователей читают эту тему (0 Гостей и 0 Скрытых Пользователей)
0 Пользователей:
« Предыдущая тема \| Центр помощи \| Следующая тема »