From 44ed843999b686c8677d7affda4d2599c3549a8b Mon Sep 17 00:00:00 2001 From: AZEN-SGG Date: Tue, 8 Apr 2025 12:54:42 +0300 Subject: [PATCH] =?UTF-8?q?=D0=9D=D0=B5=20=D1=80=D0=B0=D0=B1=D0=BE=D1=82?= =?UTF-8?q?=D0=B5=D1=82...?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- 2025.04.04/dist/Krivoruchenko_SK/Makefile | 3 ++ 2025.04.04/dist/Krivoruchenko_SK/main.c | 2 +- 2025.04.04/dist/Krivoruchenko_SK/matrix.c | 2 +- 2025.04.04/dist/Krivoruchenko_SK/solve.c | 55 ++++++++++++++--------- 4 files changed, 40 insertions(+), 22 deletions(-) diff --git a/2025.04.04/dist/Krivoruchenko_SK/Makefile b/2025.04.04/dist/Krivoruchenko_SK/Makefile index dd7fc07..f3a7bf4 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/Makefile +++ b/2025.04.04/dist/Krivoruchenko_SK/Makefile @@ -19,5 +19,8 @@ OBJ = main.o solve.o array_io.o init_f.o matrix.o $(TARGET): $(OBJ) gcc $^ -o $@ -lm +gdb: CFLAGS = -mfpmath=sse -std=gnu99 -g +gdb: clean $(TARGET) + clean: rm -f *.o *.out diff --git a/2025.04.04/dist/Krivoruchenko_SK/main.c b/2025.04.04/dist/Krivoruchenko_SK/main.c index 4c552dd..53ba17d 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/main.c +++ b/2025.04.04/dist/Krivoruchenko_SK/main.c @@ -67,7 +67,7 @@ int main(int argc, char *argv[]) t = clock(); res = t14_solve(n, a, x, c); - t = clock() - t; + t = (clock() - t) / CLOCKS_PER_SEC; if (res == SINGULAR) { diff --git a/2025.04.04/dist/Krivoruchenko_SK/matrix.c b/2025.04.04/dist/Krivoruchenko_SK/matrix.c index bd6cc82..ddea7e4 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/matrix.c +++ b/2025.04.04/dist/Krivoruchenko_SK/matrix.c @@ -93,6 +93,6 @@ double get_matrix_norm(const int n, const double * restrict A) if (sum > norm) norm = sum; } - + return norm; } diff --git a/2025.04.04/dist/Krivoruchenko_SK/solve.c b/2025.04.04/dist/Krivoruchenko_SK/solve.c index 35917b0..78d9994 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/solve.c +++ b/2025.04.04/dist/Krivoruchenko_SK/solve.c @@ -1,12 +1,20 @@ #include "solve.h" #include "io_status.h" #include "array_io.h" +#include "matrix.h" #include #include +#include // 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 = DBL_EPSILON*norm; + + if (norm < DBL_EPSILON) + return SINGULAR; + // Проходимся по главным минорам for (int k = 0; k < n; ++k) { double maximum = -1.; @@ -25,7 +33,7 @@ int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c) } // Если максимальный по модулю элемент равен нулю, значит матрица вырождена - if (fabs(maximum) < DBL_EPSILON) + if (fabs(maximum) < eps) return SINGULAR; // Меняем строки местами, если максимум находится не в k строке @@ -71,6 +79,12 @@ int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c) } gauss_inverse(n, k, A, X); + + printf("\n------- K = %d -------\n", k); + printf("Original matrix:\n"); + print_matrix(A, n, n); + printf("Inverse matrix:\n"); + print_matrix(X, n, n); } gauss_back_substitution(n, A, X); @@ -118,37 +132,32 @@ void gauss_inverse(const int n, const int k, double * restrict A, double * restr { const int kn = k*n; const int kk = kn + k; - const double inv_akk = 1./A[kn + k]; - A[kn + k] = 1.; + const double inv_akk = 1./A[kk]; + A[kk] = 1.; - for (int ij = kn; ij <= kn+k; ++ij) + for (int ij = kk+1; ij < kn+n; ++ij) { - double xij = X[ij]; - if (fabs(xij) > DBL_EPSILON) X[ij] = xij*inv_akk; - } - - for (int ij = kn + k+1; ij < kn+n; ++ij) - { - double aij = A[ij], xij = X[ij]; + double aij = A[ij]; if (fabs(aij) > DBL_EPSILON) A[ij] = aij*inv_akk; - if (fabs(xij) > DBL_EPSILON) X[ij] = xij*inv_akk; } + for (int ij = kn; ij < kk; ++ij) + X[ij] *= inv_akk; + + X[kk] = inv_akk; + for (int i = k+1; i < n; ++i) { const int in = i*n; const double aik = A[in + k]; A[in + k] = 0; - X[in + k] -= X[kk] * aik; + X[in + k] = -inv_akk * aik; - for (int j = 0; j < k; ++j) - X[in + j] -= X[kn + j] * aik; + for (int ij = in, kj = kn; kj < kk; ++ij, ++kj) + X[ij] -= X[kj] * aik; - for (int j = k+1; j < n; ++j) - { - A[in + j] -= A[kn + j] * aik; - X[in + j] -= X[kn + j] * aik; - } + for (int ij = in+k+1, kj = kk+1; ij < in+n; ++ij, ++kj) + A[ij] -= A[kj] * aik; } } @@ -169,6 +178,12 @@ void gauss_back_substitution(const int n, double * restrict A, double * restrict for (int j = 0; j < n; ++j) X[in + j] -= X[kn + 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); } }