Метод Гаусса посредством использования valarray

@bobiko · Регистрация: 16.11.2015

Author24 — интернет-сервис помощи студентам

Здравствуйте, дорогие форумчане!
Задался идеей реализовать метод гаусса посредством использования valarray, ну уже пожалел об этом. Идея в чем: хотел использовать 1 vallarray для всех вычислений, но не вижу это возможным, поэтому использую 1 valarray для хранения всех матричных значений + ответы и 2 для хранения временных значений, чтобы при преобразовании их не потерять. Для простоты решил задать все значения сразу, прямой ход делал без хранения ответов в виде нулей и все получалось, а как только добавил их, сразу несуразные неопределенности в выводе. Как я понимаю, где-то оно дальше считает и не видит выхода, но хоть убей, не могу найти ошибку. Потом еще нужен обратный ход же сделать. Вопрос: я правильно начал делать его, или нет?
Исходник

Кликните здесь для просмотра всего текста

C++

#include "stdafx.h"
#include <iostream>
#include <stdlib.h>
#include <valarray>
 
using namespace std;
 
int main()
{
    int i, j, n, m, a, b, c, d, k, l;
    cout << "Number of unknown quantities: " << endl;
    cin >> m;
    n = (m + 1) *m;
    valarray<float> gauss = { 2, 1, 1, 1, 1, 2, 1, 1, 1, 1, 2, 1, 0, 0, 0 }; //там, где нули будут сохранятся ответы
    valarray<float> tmp1(1);
    valarray<float> tmp2(1);
    l = gauss.size();
 
    for (i = 0; i < m; i++) {
        for (j = 0; j < m - i; j++) {
                a = (m + 1) * i;            //значимый диагональный элемент строки
                b = a + ((m + 1)*(j + 1));  // индекс соответствующей ему неизвестной в след.строке
                d = a + ((m + 1)*(j + 1));  //индексы строки
                    tmp1 = gauss[slice(a + i, 1, 1)]; //сохраняем а, чтобы не потерять при изменениях
                    tmp2 = gauss[slice(b + i, 1, 1)]; //также с b
                    for (k = 0; k < (m + 1) - i; k++) { gauss[slice(a + i + k, 1, 1)] /= tmp1; }  //делим первую строку на а
 
                    for (k = 0; k < (m + 1) - i; k++) { gauss[slice(a + i + k, 1, 1)] *= tmp2; } //умножаем первую строку на b
 
                    for (k = 0; k < (m + 1) - i; k++) { gauss[slice(d + i + k, 1, 1)] -= gauss[slice(a + i + k, 1, 1)]; } // отнимает от второй строки преобразованную первую
 
                    for (k = 0; k < (m + 1) - i; k++) { gauss[slice(a + i + k, 1, 1)] /= tmp2; } // возвращаем вид первой в исходное
 
 
                    for (size_t i = 0; i < 4; i++) { cout << gauss[i] << ' '; }cout << endl;
                    for (size_t i = 4; i < 8; i++) { cout << gauss[i] << ' '; }cout << endl;
                    for (size_t i = 8; i < 12; i++) { cout << gauss[i] << ' '; } cout << endl << endl;
                }
 
    }
    //якобы обратный ход
    for (i = 0; i <= m; i++) {
        tmp1 = gauss[slice(n - 1, 1, 1)];
        for (j = 1; j <=i; j++) {
            a = (n - 1) - ((m + 1) * i);
            tmp1 = gauss[slice(a, 1, 1)];
            tmp2 = gauss[slice(a - j, 1, 1)];
            tmp2 *= gauss[slice(l - (j - 1), 1, 1)];
            tmp1 -= tmp2;
        }
        gauss[slice(l - i, 1, 1)] = tmp1;
    }
 
    for (size_t i = 0; i < 4; i++) { cout << gauss[i] << ' '; }cout << endl;
    for (size_t i = 4; i < 8; i++) { cout << gauss[i] << ' '; }cout << endl;
    for (size_t i = 8; i < 12; i++) { cout << gauss[i] << ' '; } cout << endl << endl;
    for (size_t i = l - m; i < l; i++) {
        cout << ' ' << gauss[i] << endl;
    }
    //check --- http://matematikam.ru/solve-equations/sistema-gaus.php
    return 0;
}

Добавлено через 3 часа 35 минут
Все, тему можно закрывать, уже все решил.) Волшебный форум какой-то, ведь как только написал, так сразу все замечать начал)

@AndrSlav · 15.12.2017, 22:22

Понимаю, что уже не надо, но я как раз пытаюсь изучать библиотеку

, так что попробовал решить тоже СЛАУ методом Гаусса. Может кто что подскажет насчет скрытых ошибок. Так вроде решается.

C++

#pragma hdrstop
#pragma argsused
 
#ifdef _WIN32
#include <tchar.h>
#else
  typedef char _TCHAR;
  #define _tmain main
#endif
 
#include <stdio.h>
#include <iostream>
#include <vector>
#include <algorithm>
#include <functional>
#include <string>
#include <memory>
#include <math.h>
#include <valarray>
#include <iterator>
   using std::vector;
   using std::valarray;
   using std::slice;
   using std::slice_array;
 
   int amax_val(const valarray<double>& v, vector<bool> trig)
   {int m= 0;
    for(int i=0; i<v.size(); ++i, fabs(v[i])>fabs(v[m]) && trig[i] ? m=i : m);
    return m;
   };
 
 int _tmain(int argc, _TCHAR* argv[])
{
  int dim= 3;
  vector<int> ind(3);
  vector<bool> trig;
  for(int i=0; i<dim; ++i, trig.push_back(true));
 
  std::cout<<"Initial matrix:";
  double mas[]= {1,5,4,2,9,15,7,2,4,3,8,2};
  valarray<double> matr(mas, 12);    //Матрица коэффициентов
  valarray<double> col(dim);         //Столбец
  valarray<double> rr(dim+1);        //Строка
  vector<double> rez;
  int m;
 
          for(int i=0; i<matr.size(); ++i)
          {
           if(!(i%(dim+1))) std::cout<<'\n';
           std::cout<<matr[i]<<"  ";
          };
 
  for(int i=0; i<dim; ++i)
  {
   col= matr[slice(i, dim, dim+1)];
   m= amax_val(col, trig);
   ind[i]= m;
   trig[m]= false;
   const slice_array<double>& row= matr[slice(m*(dim+1), dim+1, 1)];
   row/= (rr= col[m]);
 
        for(int j=0; j<dim; ++j)
        {
             if(j!=m)
             {
              const slice_array<double>& row2= matr[slice(j*(dim+1), dim+1, 1)];
              row2/= (rr= col[j]);
              row2-= matr[slice(m*(dim+1), dim+1, 1)];
             };
        };
  };
 
  for(int i=0; i<dim; ++i)
  {
   const slice_array<double>& row= matr[slice(ind[i]*(dim+1), dim+1, 1)];
   row/= (rr= matr[i+ ind[i]*(dim+1)]);
   rez.push_back(matr[(ind[i]+1)*(dim+1)- 1]);
  };
 
 
  std::cout<<'\n'<<'\n'<<"Rezult: ";
  for(int i=0; i<matr.size(); ++i)
  {
   if(!(i%(dim+1))) std::cout<<'\n';
   std::cout<<matr[i]<<"  ";
  };
 
  std::cout<<'\n'<<'\n';
  std::cout<<rez[0]<<" "<<rez[1]<<" "<<rez[2]<<'\n';
 
  char ch;
  std::cin>>ch;
 
    return 0;
}

	15.12.2017, 22:22