diff --git a/2025.04.04/dist/Krivoruchenko_SK/Makefile b/2025.04.04/dist/Krivoruchenko_SK/Makefile index cb41b2f..1536e8e 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/Makefile +++ b/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: diff --git a/2025.04.04/dist/Krivoruchenko_SK/solve.c b/2025.04.04/dist/Krivoruchenko_SK/solve.c index 292a871..57fdaf3 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/solve.c +++ b/2025.04.04/dist/Krivoruchenko_SK/solve.c @@ -4,118 +4,10 @@ #include "matrix.h" #include #include -#include #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; +} + + + diff --git a/2025.04.04/dist/Krivoruchenko_SK/solve.h b/2025.04.04/dist/Krivoruchenko_SK/solve.h index 2e54ac9..ac8d030 100644 --- a/2025.04.04/dist/Krivoruchenko_SK/solve.h +++ b/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 diff --git a/2025.04.04/dist/Ulyanov_MT/Makefile b/2025.04.04/dist/Ulyanov_MT/Makefile new file mode 100644 index 0000000..701ed29 --- /dev/null +++ b/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 diff --git a/2025.04.04/dist/Ulyanov_MT/array.c b/2025.04.04/dist/Ulyanov_MT/array.c new file mode 100644 index 0000000..54aae42 --- /dev/null +++ b/2025.04.04/dist/Ulyanov_MT/array.c @@ -0,0 +1,81 @@ +#include +#include +#include +#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; + } diff --git a/2025.04.04/dist/Ulyanov_MT/array.h b/2025.04.04/dist/Ulyanov_MT/array.h new file mode 100644 index 0000000..39c4ed3 --- /dev/null +++ b/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); diff --git a/2025.04.04/dist/Ulyanov_MT/io_status.h b/2025.04.04/dist/Ulyanov_MT/io_status.h new file mode 100644 index 0000000..5e603af --- /dev/null +++ b/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; diff --git a/2025.04.04/dist/Ulyanov_MT/matrix.c b/2025.04.04/dist/Ulyanov_MT/matrix.c new file mode 100644 index 0000000..af34d91 --- /dev/null +++ b/2025.04.04/dist/Ulyanov_MT/matrix.c @@ -0,0 +1,63 @@ +#include +#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; + } diff --git a/2025.04.04/dist/Ulyanov_MT/matrix.h b/2025.04.04/dist/Ulyanov_MT/matrix.h new file mode 100644 index 0000000..93ce19f --- /dev/null +++ b/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); diff --git a/2025.04.04/dist/Ulyanov_MT/task20.c b/2025.04.04/dist/Ulyanov_MT/task20.c new file mode 100644 index 0000000..9743ba2 --- /dev/null +++ b/2025.04.04/dist/Ulyanov_MT/task20.c @@ -0,0 +1,176 @@ +#include +#include +#include +#include +#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; + } diff --git a/2025.04.18/Tasks09.pdf b/2025.04.18/Tasks09.pdf new file mode 100644 index 0000000..5e21bf5 Binary files /dev/null and b/2025.04.18/Tasks09.pdf differ