Добавил решение Матвея

2025.04.04/dist/Krivoruchenko_SK/Makefile

@@ -34,7 +34,7 @@ gdb: CFLAGS = -mfpmath=sse -std=gnu99 -g -O0
 gdb: clean $(TARGET)
 
 # Профилировочная сборка (gprof)
-prof: CFLAGS = -mfpmath=sse -std=gnu99 -pg -O3 -fno-omit-frame-pointer
+prof: CFLAGS = -mfpmath=sse -std=gnu99 -pg -O3
 prof: clean $(TARGET)
 
 clean:

2025.04.04/dist/Krivoruchenko_SK/solve.c

@@ -4,118 +4,10 @@
 #include "matrix.h"
 #include <float.h>
 #include <math.h>
-#include <stdio.h>
 
 #define EPS 1.2e-16
 
-// c - changes in rows
-int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c)
-{
-	double norm = get_matrix_norm(n, A);
-	double eps = EPS*norm;
-
-//	printf("NORM = %lf EPS = %lf\n", norm, eps);
-
-	// Проходимся по главным минорам
-	for (int k = 0; k < n; ++k) {
-		double maximum = -1.;
-		int max_i = 0, max_j = 0;
-		
-		// Ищем максимальный элемент минора
-		for (int i = k; i < n; ++i)
-			for (int j = k; j < n; ++j)
-				if (fabs(A[i*n + j]) > maximum) {
-					maximum = fabs(A[i*n + j]);
-					max_i = i;
-					max_j = j;
-				}
-			
-//		printf("\n------- K = %d -------\n", k);
-//		printf("Maximum = %lf i = %d j = %d\n", maximum, max_i, max_j);
-
-		// Если максимальный по модулю элемент равен нулю, значит матрица вырождена
-		if (fabs(maximum) <= eps)
-			return SINGULAR;
-		
-		// Меняем строки местами, если максимум находится не в k строке
-		if (max_i != k)
-		{
-			for (int j = 0; j < k; ++j)
-			{
-				double swap = X[k*n + j];
-				X[k*n + j] = X[max_i*n + j];
-				X[max_i*n + j] = swap;
-			}
-	
-			for (int j = k; j < n; ++j)
-			{
-				double swap = X[k*n + j];
-				X[k*n + j] = X[max_i*n + j];
-				X[max_i*n + j] = swap;
-				
-				swap = A[k*n + j];
-				A[k*n + j] = A[max_i*n + j];
-				A[max_i*n + j] = swap;
-			}
-		}
-
-		// Меняем столбцы местами
-		if (max_j != k)
-		{
-			int swap_temp = c[max_j];
-			c[max_j] = c[k];
-			c[k] = swap_temp;
-
-			for (int i = 0; i < n; i++)
-			{
-				double swap = A[i*n + k];
-				A[i*n + k] = A[i*n + max_j];
-				A[i*n + max_j] = swap;
-			}
-		}
-
-//		printf("BEFORE GAUSS\n");
-//		printf("Original matrix:\n");
-//		print_matrix(A, n, n);
-//		printf("Inverse matrix:\n");
-//		print_matrix(X, n, n);
-		
-		gauss_inverse(n, k, A, X);
-		
-//		printf("AFTER GAUSS\n");
-//		printf("Original matrix:\n");
-//		print_matrix(A, n, n);
-//		printf("Inverse matrix:\n");
-//		print_matrix(X, n, n);
-	}
-
-	gauss_back_substitution(n, A, X);
-
-	for (int k = 0; k < n; ++k)
-	{
-		const int kn = k*n;	
-		int i = c[k];
-		
-		while (i != k)
-		{
-			const int in = i*n;
-			const int swap_int = c[i];
-			c[i] = i;
-			i = swap_int;
-
-			for (int j = 0; j < n; ++j)
-			{
-				double swap_temp = X[in+j];
-				X[in+j] = X[kn+j];
-				X[kn+j] = swap_temp;
-			}
-		}
-	}
-
-	return 0;
-}
-
-// Прямой ход Го ----- йда
+// Прямой ход Гауса
 void gauss_inverse(const int n, const int k, double * restrict A, double * restrict X)
 {
 	const double inv_akk = 1./A[k*n + k];
@@ -156,11 +48,97 @@ void gauss_back_substitution(const int n, double * restrict A, double * restrict
 			for (int j = 0; j < n; ++j)
 				X[i*n + j] -= X[k*n + j] * aik;	
 		}
-
-//		printf("\n------- K = %d -------\n", k);
-//		printf("Original matrix:\n");
-//		print_matrix(A, n, n);
-//		printf("Inverse matrix:\n");
-//		print_matrix(X, n, n);
 }
 
+// c - changes in rows
+int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c)
+{
+	double norm = get_matrix_norm(n, A);
+	double eps = EPS*norm;
+
+	// Проходимся по главным минорам
+	for (int k = 0; k < n; ++k) {
+		double maximum = -1.;
+		int max_i = 0, max_j = 0;
+		
+		// Ищем максимальный элемент минора
+		for (int i = k; i < n; ++i)
+			for (int j = k; j < n; ++j)
+				if (fabs(A[i*n + j]) > maximum) {
+					maximum = fabs(A[i*n + j]);
+					max_i = i;
+					max_j = j;
+				}
+			
+		// Если максимальный по модулю элемент равен нулю, значит матрица вырождена
+		if (fabs(maximum) <= eps)
+			return SINGULAR;
+		
+		// Меняем строки местами, если максимум находится не в k строке
+		if (max_i != k)
+		{
+			for (int j = 0; j < k; ++j)
+			{
+				double swap = X[k*n + j];
+				X[k*n + j] = X[max_i*n + j];
+				X[max_i*n + j] = swap;
+			}
+	
+			for (int j = k; j < n; ++j)
+			{
+				double swap = X[k*n + j];
+				X[k*n + j] = X[max_i*n + j];
+				X[max_i*n + j] = swap;
+				
+				swap = A[k*n + j];
+				A[k*n + j] = A[max_i*n + j];
+				A[max_i*n + j] = swap;
+			}
+		}
+
+		// Меняем столбцы местами
+		if (max_j != k)
+		{
+			int swap_temp = c[max_j];
+			c[max_j] = c[k];
+			c[k] = swap_temp;
+
+			for (int i = 0; i < n; i++)
+			{
+				double swap = A[i*n + k];
+				A[i*n + k] = A[i*n + max_j];
+				A[i*n + max_j] = swap;
+			}
+		}
+		
+		gauss_inverse(n, k, A, X);
+	}
+
+	gauss_back_substitution(n, A, X);
+
+	for (int k = 0; k < n; ++k)
+	{
+		const int kn = k*n;	
+		int i = c[k];
+		
+		while (i != k)
+		{
+			const int in = i*n;
+			const int swap_int = c[i];
+			c[i] = i;
+			i = swap_int;
+
+			for (int j = 0; j < n; ++j)
+			{
+				double swap_temp = X[in+j];
+				X[in+j] = X[kn+j];
+				X[kn+j] = swap_temp;
+			}
+		}
+	}
+
+	return 0;
+}
+
+
+

2025.04.04/dist/Krivoruchenko_SK/solve.h

@@ -1,8 +1,8 @@
 #ifndef SOLVE_H
 #define SOLVE_H
 
-int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c);
 void gauss_inverse(const int n, const int k, double * restrict A, double * restrict X);
 void gauss_back_substitution(const int n, double * restrict A, double * restrict X);
+int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c);
 
 #endif

2025.04.04/dist/Ulyanov_MT/Makefile

@@ -0,0 +1,12 @@
+FLAGS = -mfpmath=sse -fstack-protector-all -W -Wall -Wextra -Wunused -Wcast-align -Werror -pedantic -pedantic-errors -Wfloat-equal -Wpointer-arith -Wformat-security -Wmissing-format-attribute -Wformat=1 -Wwrite-strings -Wcast-align -Wno-long-long -std=gnu99 -Wstrict-prototypes -Wmissing-prototypes -Wmissing-declarations -Wold-style-definition -Wdeclaration-after-statement -Wbad-function-cast -Wnested-externs -Wmaybe-uninitialized -O3
+all: a.out
+a.out: a.o array.o matrix.o
+	gcc a.o array.o matrix.o -lm -o a.out
+array.o:
+	gcc -c $(FLAGS) -o array.o array.c
+matrix.o:
+	gcc -c $(FLAGS) -o matrix.o matrix.c
+a.o:
+	gcc -c $(FLAGS) -o a.o task20.c
+clean:
+	rm -f *.o *.out

2025.04.04/dist/Ulyanov_MT/array.c

@@ -0,0 +1,81 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include "io_status.h"
+#include "array.h"
+#include "matrix.h"
+
+io_status read_matrix(double* a, int n, int m, const char* name)
+  {
+  int i, j;
+  FILE* fp;
+  if (!(fp = fopen(name, "r")))
+      return ERROR_OPEN;
+  for (i = 0; i < n; i++)
+      {
+      for (j = 0; j < m; j++)
+          {
+          if (fscanf(fp, "%lf", a + i * m + j) != 1)
+              {
+              fclose(fp);
+              return ERROR_READ;
+              }
+          }
+      }
+  fclose(fp);
+  return SUCCESS;
+  }
+  
+void print_matrix(const double* a, int n, int m, int p)
+  {
+  int np = (n > p ? p : n);
+  int mp = (m > p ? p : m);
+  int i, j;
+  for (i = 0; i < np; i++)
+      {
+      for (j = 0; j < mp; j++)
+          printf(" %10.3e", a[i * m + j]);
+      printf("\n");
+      }
+  }
+  
+void init_matrix(double* a, int n, int m, int k)
+  {
+  int i, j;
+  for (i = 0; i < n; i++)
+      {
+      for (j = 0; j < m; j++)
+          a[i * m + j] = f(k, n, m, i+1, j+1);
+      }
+  }
+  
+double f(int k, int n, int m, int i, int j)
+  {
+  switch (k)
+      {
+      case 1: return (n >= m ? n : m) - (i >= j ? i : j) + 1;
+      case 2: return (i >= j ? i : j);
+      case 3: return (i - j >= 0 ? i - j : j - i);
+      case 4: return 1./(i + j - 1);
+      }
+  return -1e308;
+  }
+  
+void init_idmatrix(double* b, int n)
+  {
+  int i, j;
+  for (i = 0; i < n; i++)
+      {
+      for (j = 0; j < n; j++)
+          {
+          if (i == j) b[i * n + j] = 1;
+          else b[i * n + j] = 0;
+          }
+      }
+  }
+  
+void init_ind(int* ind, int n)
+  {
+  int i;
+  for (i = 0; i < n; i++) ind[i] = i;
+  }

2025.04.04/dist/Ulyanov_MT/array.h

@@ -0,0 +1,6 @@
+io_status read_matrix(double*, int, int, const char*);
+void print_matrix(const double*, int, int, int);
+void init_matrix(double*, int, int, int);
+double f(int, int, int, int, int);
+void init_idmatrix(double*, int);
+void init_ind(int*, int);

2025.04.04/dist/Ulyanov_MT/io_status.h

@@ -0,0 +1,7 @@
+typedef enum io_status_
+  {
+  SUCCESS,
+  ERROR_OPEN,
+  ERROR_READ,
+  ERROR_MATRIX,
+  } io_status;

2025.04.04/dist/Ulyanov_MT/matrix.c

@@ -0,0 +1,63 @@
+#include <math.h>
+#include "matrix.h"
+
+double count_r1(double* a, double* b, int n)
+  {
+  int i, j, u;
+  double sum, sumstl, max = 0;
+  if (n > 4000) return 0;
+  for (j = 0; j < n; j++)
+      {
+      sumstl = 0;
+      for (i = 0; i < n; i++)
+          {
+          if (i == j) sum = -1;
+          else sum = 0;
+          for (u = 0; u < n; u++)
+              sum += a[i * n + u] * b[u * n + j];
+          sumstl += fabs(sum);
+          }
+      if (sumstl > max) max = sumstl;
+      }
+  return max;
+  }
+
+double count_r2(double* a, double* b, int n)
+  {
+  int i, j, u;
+  double sum, sumstl, max = 0;
+  if (n > 4000) return 0;
+  for (j = 0; j < n; j++)
+      {
+      sumstl = 0;
+      for (i = 0; i < n; i++)
+          {
+          if (i == j) sum = -1;
+          else sum = 0;
+          for (u = 0; u < n; u++)
+              sum += b[i * n + u] * a[u * n + j];
+          sumstl += fabs(sum);
+          }
+      if (sumstl > max) max = sumstl;
+      }
+  return max;
+  }
+  
+double count_norm(double* a, int n)
+  {
+  int i, j;
+  double max = 0, sum;
+  for (i = 0; i < n; i++)
+      {
+      sum = 0;
+      for (j = 0; j < n; j++) sum += fabs(a[i * n + j]);
+      if (sum > max) max = sum;
+      }
+  return max;
+  }
+  
+int is_null(double a, double norm)
+  {
+  if (fabs(a) <= EPS * norm) return 1;
+  return 0;
+  }

2025.04.04/dist/Ulyanov_MT/matrix.h

@@ -0,0 +1,6 @@
+#define EPS 112e-18
+
+double count_r1(double*, double*, int);
+double count_r2(double*, double*, int);
+double count_norm(double*, int);
+int is_null(double, double);

2025.04.04/dist/Ulyanov_MT/task20.c

@@ -0,0 +1,176 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <time.h>
+#include <math.h>
+#include "io_status.h"
+#include "array.h"
+#include "matrix.h"
+
+io_status inverse(double*, double*, int*, int);
+
+int main(int argc, char* argv[])
+    {
+    io_status res;
+    int task = 20;
+    double* a;
+    double* b;
+    int* ind;
+    int n, p, k;
+    char* name = 0;
+    double r1, r2;
+    double t;
+     if (!((argc == 4 || argc == 5) && sscanf(argv[1], "%d", &n) == 1 && sscanf(argv[2], "%d", &p) == 1 && sscanf(argv[3], "%d", &k) == 1 && k >= 0 && k <= 4))
+        {
+        printf("Usage: %s n p k [file]\n", argv[0]);
+        return 1;
+        }
+    if (k == 0) name = argv[4];
+    a = (double*) malloc(n * n * sizeof(double));
+    if (!a)
+        {
+        printf("Not enough memory\n");
+        return 2;
+        }
+    if (name)
+        {
+        io_status ret;
+        ret = read_matrix(a, n, n, name);
+        if (ret != SUCCESS)
+            {
+            switch(ret)
+                {
+                case ERROR_OPEN:
+                    printf("Can not open file %s\n", name);
+                    break;
+                case ERROR_READ:
+                    printf("Can not read file %s\n", name);
+                    break;
+                default:
+                    printf("Unknown error %d in file %s\n", ret, name);
+                }
+            free(a);
+            return 3;
+            }
+        }
+    else
+        init_matrix(a, n, n, k);
+    b = (double*) malloc(n * n * sizeof(double));
+    if (!b)
+        {
+        printf("Not enough memory\n");
+        free(a);
+        return 2;
+        }
+    init_idmatrix(b, n);
+    ind = (int*) malloc(n * sizeof(int));
+    if (!ind)
+        {
+        printf("Not enough memory\n");
+        free(a);
+        free(b);
+        return 2;
+        }
+    init_ind(ind, n);
+    printf("Matrix:\n");
+    print_matrix(a, n, n, p);
+    t = clock();
+    res = inverse(a, b, ind, n);
+    t = (clock() - t) / CLOCKS_PER_SEC;
+    if (res == ERROR_MATRIX)
+        {
+        printf("Method is not applicable or inverse matrix does not exist\n");
+        free(a);
+        free(b);
+        free(ind);
+        return 0;
+        }
+    if (name) read_matrix(a, n, n, name);
+    else init_matrix(a, n, n, k);
+    r1 = count_r1(a, b, n);
+    r2 = count_r2(a, b, n);
+    printf("Inverse matrix:\n");
+    print_matrix(b, n, n, p);
+    printf ("%s : Task = %d Res1 = %e Res2 = %e Elapsed = %.2f K = %d N = %d\n", argv[0], task, r1, r2, t, k, n);
+    free(a);
+    free(b);
+    free(ind);
+    return 0;
+    }
+    
+io_status inverse(double* a, double* b, int* ind, int n)
+  {
+  int tmp, i, j, k, ind_max_i, ind_max_j;
+  double bag, mod_el_a, el_aik, el_akk, max, norm;
+  norm = count_norm(a, n);
+  for (k = 0; k < n; k++)
+      {
+      max = 0;
+      ind_max_i = k;
+      ind_max_j = k;
+      for (i = k; i < n; i++)
+          {
+          for (j = k; j < n; j++)
+              {
+              mod_el_a = fabs(a[i * n + j]);
+              if (mod_el_a > max)
+                  {
+                  max = mod_el_a;
+                  ind_max_i = i;
+                  ind_max_j = j;
+                  }
+              }
+          }
+      for (j = k; j < n; j++)
+          {
+          bag = a[k * n + j];
+          a[k * n + j] = a[ind_max_i * n + j];
+          a[ind_max_i * n + j] = bag;
+          }
+      for (j = 0; j < n; j++)
+          {
+          bag = b[k * n + j];
+          b[k * n + j] = b[ind_max_i * n + j];
+          b[ind_max_i * n + j] = bag;
+          }
+      for (i = 0; i < n; i++)
+          {
+          bag = a[i * n + k];
+          a[i * n + k] = a[i * n + ind_max_j];
+          a[i * n + ind_max_j] = bag;
+          }    
+      tmp = ind[k];
+      ind[k] = ind[ind_max_j];
+      ind[ind_max_j] = tmp;
+      el_akk = a[k * n + k];
+      if (is_null(el_akk, norm)) return ERROR_MATRIX;
+      for (j = k + 1; j < n; j++) a[k * n + j] = a[k * n + j] / el_akk;
+      for (j = 0; j < n; j++) b[k * n + j] = b[k * n + j] / el_akk;
+      for (i = 0; i < n; i++)
+          {
+          if (i != k)
+              {
+              el_aik = a[i * n + k];
+              for (j = k + 1; j < n; j++) a[i * n + j] = a[i * n + j] - el_aik * a[k * n + j];
+              for (j = 0; j < n; j++) b[i * n + j] = b[i * n + j] - el_aik * b[k * n + j];
+              }
+          }
+      }
+      tmp = -1;
+      for (i = 0; i < n; i++)
+          {
+          while (i != tmp)
+              {
+              tmp = ind[i];
+              for (j = 0; j < n; j++)
+                  {
+                  bag = b[tmp * n + j];
+                  b[tmp * n + j] = b[i * n + j];
+                  b[i * n + j] = bag;
+                  }
+              ind[i] = ind[tmp];
+              ind[tmp] = tmp;
+              }
+          }
+
+  return SUCCESS;
+  }