Решение систем линейных алгебраических уравнений методом Гаусса

      Writeln; 

      Writeln('Введите расширенную матрицу  системы');

      ReadSystem(n, a, b);

      Writeln; 

      If Gauss(n, a, b, x) then begin

         Writeln('Результат вычислений по методу Гаусса');

         WriteX(n, x);

      end

      else

          Writeln('Данную систему невозможно  решить по методу Гаусса');

      Writeln;

End. 

Программа решения систем линейных уравнений по методу Гаусса 

Введите порядок  матрицы системы (макс. 10)

>4 

Введите расширенную  матрицу системы

A     1     2     3     4     b 

1     3.2   5.4   4.2   2.2   2.6

2     2.1   3.2   3.1   1.1   4.8

3     1.2   0.4   -0.8  -0.8  3.6

4     4.7   10.4  9.7   9.7   -8.4 

Результат вычислений по методу Гаусса

x1 =  5.0000000000E+00

x2 = -4.0000000000E+00

x3 =  3.0000000000E+00

x4 = -2.0000000000E+00 

Пример №2: 
Решение системы линейных уравнений (возможно переопределенной) методом Гаусса. Определяется ситуация, что система не имеет решений. Ситуация, когда система имеет более чем одно решение не рассматривается. В случае удачного завершения возвращает нуль. 
 
function LinGauss(M, N: Integer; Data: PExtended; X: PExtended): Cardinal; 
var 
PtrData: PExtended; 
PtrData1, PtrData2: PExtended; 
Temp: Extended; 
I, J, Row: Integer; 
Max: Extended; 
MaxR: Integer; 
begin 
Assert(M >= N, 'Invalid start data'); 
for I := 0 to N - 1 do // Для каждой переменной 
begin 
 
// 1. Поиск максимального элемента 
PChar(PtrData) := PChar(Data) + I * (N + 2) * SizeOf(Extended); 
MaxR := I; 
Max := PtrData^; 
for J := I + 1 to M - 1 do 
begin 
PChar(PtrData) := PChar(PtrData) + (N + 1) * SizeOf(Extended); 
if Abs(PtrData^) > Abs(Max) then 
begin 
Max := PtrData^; 
MaxR := J; 
end; 
end; 
// 2. А вдруг неразрешима? 
if Abs(Max) < 1.0E-10 then 
begin 
Result := $FFFFFFFF; 
Exit; 
end; 
// 3. Меняем местами строки 
if MaxR <> I then 
begin 
PChar(PtrData1) := PChar(Data) + MaxR * (N + 1) * SizeOf(Extended); 
PChar(PtrData2) := PChar(Data) + I * (N + 1) * SizeOf(Extended); 
for J := 0 to N do 
begin 
Temp := PtrData1^; 
PtrData1^ := PtrData2^; 
PtrData2^ := Temp; 
PChar(PtrData1) := PChar(PtrData1) + SizeOf(Extended); 
PChar(PtrData2) := PChar(PtrData2) + SizeOf(Extended); 
end; 
end; 
// 4. Пересчет направляющей строки 
PChar(PtrData) := PChar(Data) + I * (N + 1) * SizeOf(Extended); 
for J := 0 to N do 
begin 
PtrData^ := PtrData^ / Max; 
PChar(PtrData) := PChar(PtrData) + SizeOf(Extended); 
end; 
// 5. Пересчет всей оставшейся части таблицы 
PtrData1 := Data; 
for Row := 0 to M - 1 do 
begin 
if Row = I then 
begin 
PChar(PtrData1) := PChar(PtrData1) + (N + 1) * SizeOf(Extended); 
Continue; 
end; 
PChar(PtrData2) := PChar(Data) + I * (N + 1) * SizeOf(Extended); 
Temp := PExtended(PChar(PtrData1) + I * SizeOf(Extended))^; 
for J := 0 to N do 
begin 
PtrData1^ := PtrData1^ - Temp * PtrData2^; 
PChar(PtrData1) := PChar(PtrData1) + SizeOf(Extended); 
PChar(PtrData2) := PChar(PtrData2) + SizeOf(Extended); 
end; 
end; 
end; 
// 6. Проверка того, что система переопределена 
PChar(PtrData) := PChar(Data) + N * (N + 1) * SizeOf(Extended); 
for I := N to M - 1 do 
for J := 0 to N do 
begin 
if Abs(PtrData^) > 1.0E-10 then 
begin 
Result := $FFFFFFFF; 
Exit; 
end; 
PChar(PtrData) := PChar(PtrData) + SizeOf(Extended); 
end; 
// Все ОК 
PChar(PtrData) := PChar(Data) + N * SizeOf(Extended); 
for I := 0 to N - 1 do 
begin 
X^ := PtrData^; 
PChar(X) := PChar(X) + SizeOf(Extended); 
PChar(PtrData) := PChar(PtrData) + (N + 1) * SizeOf(Extended); 
end; 
Result := 0; 
end; 

Пример №3:

Решение системы  линейных алгебраических уравнений на VBA

 

Option Explicit

  Sub rewenie()

    Dim i As Integer

    Dim j As Integer

    Dim r() As Double

    Dim p As Double

    Dim x() As Double

    Dim k As Integer

    Dim n As Integer

    Dim b() As Double

    Dim file As Integer

    Dim y() As Double

     file = FreeFile

    Open "C:\data.txt" For Input As file

    Input #file, n

    ReDim x(0 To n * n - 1) As Double

    ReDim y(0 To n - 1) As Double

    ReDim r(0 To n - 1) As Double

    For i = 0 To n - 1

        For j = 0 To n - 1

            Input #file, x(i * n + j)

        Next j

        Input #file, y(i)

    Next i

    Close #file

    For i = 0 To n - 1

    p = x(i * n + i)

        For j = 1 To n - 1

            x(i * n + j) = x(i * n + j) / p

        Next j

        y(i) = y(i) / p

            For j = i + 1 To n - 1

            p = x(j * n + i)

                For k = i To n - 1

                    x(j * n + k) = x(j * n + k) - x(i * n + k) * p

                Next k

            y(j) = y(j) - y(i) * p

        Next j

Next i

' Верхнетреугольная матрица

For i = n - 1 To 0 Step -1

p = y(i)

    For j = i + 1 To n - 1

    p = p - x(i * n + j) * r(j)

    Next j

r(i) = p / x(i * n + i)

Next i

' Обратный ход

For i = 0 To n - 1

MsgBox r(i)

Next i

'<!--End of Code-->

End Sub 

Пример №4: 

//HOWTO solve system of linear equations

//ax=b

//where a[N][N] is square matrix,

//b[N] and x[N] are columns 

//data are stored in input text file: 

//input file format:

//a_11 a_12 a_13... a_1N b_1

//a_21 a_22 a_23... a_2N b_2

//...

//a_N1 a_N2 a_N3... a_NN b_N 

#include <stdio.h>

#include <process.h>

float **a, *b, *x;

int N;

char filename[256];

FILE* InFile=NULL;

void count_num_lines(){

   //count number of lines in input file - number of equations

   int nelf=0;       //non empty line flag

   do{

      nelf = 0;

       while(fgetc(InFile)!='\n' && !feof(InFile)) nelf=1;

       if(nelf) N++;

   }while(!feof(InFile));

} 

void freematrix(){

   //free memory for matrixes

   int i;

   for(i=0; i<N; i++){

       delete [] a[i];

   }

   delete [] a;

   delete [] b;

   delete [] x;

} 

void allocmatrix(){

   //allocate memory for matrixes

   int i,j;

   x = new float[N];

   b = new float[N];

   a = new float*[N];

   if(x==NULL || b==NULL || a==NULL){

       printf("\nNot enough memory to allocate for %d equations.\n", N);

       exit(-1);

   }

   for(i=0; i<N; i++){

       a[i] = new float[N];

       if(a[i]==NULL){

         printf("\nNot enough memory to allocate for %d equations.\n", N);

       }

   }

   for(i=0; i<N; i++){

       for(j=0; j<N; j++){

         a[i][j]=0;

       }

       b[i]=0;

       x[i]=0;

   }

} 

void readmatrix(){

   int i=0,j=0;

   //read matrixes a and b from input file

   for(i=0; i<N; i++){

       for(j=0; j<N; j++){

         fscanf(InFile, "%f", &a[i][j]);

       }

       fscanf(InFile, "%f", &b[i]);

   }

} 

void printmatrix(){

   //print matrix "a"

   int i=0,j=0;

   printf("\n");

   for(i=0; i<N; i++){

       for(j=0; j<N; j++){

         printf(" %+f*X%d", a[i][j], j);

       }

       printf(" =%f\n", b[i]);

   }

} 

void testsolve(){

   //test that ax=b

   int i=0,j=0;

   printf("\n");

   for(i=0; i<N; i++){

       float s = 0;

       for(j=0; j<N; j++){

         s += a[i][j]*x[j];

       }

       printf("%f\t%f\n", s, b[i]);

   }

}

void printresult(){

   int i=0;

   printf("\n");

   printf("Result\n");

   for(i=0; i<N; i++){

       printf("X%d = %f\n", i, x[i]);

   }

}

void diagonal(){

   int i, j, k;

   float temp=0;

   for(i=0; i<N; i++){

       if(a[i][i]==0){

         for(j=0; j<N; j++){

             if(j==i) continue;

             if(a[j][i] !=0 && a[i][j]!=0){

               for(k=0; k<N; k++){

                   temp = a[j][k];

                   a[j][k] = a[i][k];

                   a[i][k] = temp;

               }

               temp = b[j];

               b[j] = b[i];

               b[i] = temp;

               break;

             }

         }

       }

   }

}

void cls(){

   for(int i=0; i<25; i++) printf("\n");

}

void main(){

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

   cls();

   do{

       printf("\nInput filename: ");

       scanf("%s", filename);

       InFile = fopen(filename, "rt");

   }while(InFile==NULL);

   count_num_lines();

   allocmatrix();

   rewind(InFile);

   //read data from file

   readmatrix();

   //check if there are 0 on main diagonal and exchange rows in that case

   diagonal();

   fclose(InFile);

   printmatrix();

   //process rows

   for(k=0; k<N; k++){

       for(i=k+1; i<N; i++){

         if(a[k][k]==0){

             printf("\nSolution is not exist.\n");

             return;

         }

         float M = a[i][k] / a[k][k];

         for(j=k; j<N; j++){

             a[i][j] -= M * a[k][j];

         }

         b[i] -= M*b[k];

       }

   }

   printmatrix();

   for(i=N-1; i>=0; i--){

       float s = 0;

       for(j = i; j<N; j++){

         s = s+a[i][j]*x[j];

       }

       x[i] = (b[i] - s) / a[i][i];

   } 

   InFile = fopen(filename, "rt");

   readmatrix();

   fclose(InFile);

   printmatrix();

   testsolve();

   printresult();

   freematrix();

} 

Поскольку весь метод уже описан, коротко опишу, как эта программа получает уравнения  на вход, что выдает в результате решения.  

Чтобы воспользоваться  этой программой, нужно запустить  скомпилированный исполняемый файл. В первую очередь программа спросит, откуда ей брать коэффициенты для  уравнений. Создайте в любом текстовом  редакторе (но только не в Word-е а, например в notepad-е) файл, где напишите коэффициенты уравнений построчно через пробел, приблизительно так:

a11  a12  a13  ...  a1N  b1

a21  a22  a23  ...  a2N  b2

a31  a32  a33  ...  a3N  b3