AZEN-SGG/2nd_Sem_Bogachev
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Добавил ускорение в 14ю

Browse source
This commit is contained in:
AZEN-SGG 2025-04-15 21:09:52 +03:00
parent 9ac038df34
commit 92a3a951c2
8 changed files with 100092 additions and 71 deletions
Download patch file Download diff file Expand all files Collapse all files

1
2025.04.04/dist/Krivoruchenko_SK/Makefile vendored
View file

@ -35,6 +35,7 @@ gdb: clean $(TARGET)
# Профилировочная сборка (gprof)
prof: CFLAGS = -mfpmath=sse -std=gnu99 -pg -O3
prof: LFLAGS = -lm =pg
prof: clean $(TARGET)
clean:

122
2025.04.04/dist/Krivoruchenko_SK/solve.c vendored
View file

@ -7,69 +7,32 @@
#define EPS 1.2e-16
// Прямой ход Гауса
void gauss_inverse(const int n, const int k, double * restrict A, double * restrict X)
{
const double inv_akk = 1./A[k*n + k];
for (int j = 0; j < k+1; ++j)
X[k*n + j] *= inv_akk;
for (int j = k+1; j < n; ++j)
{
A[k*n + j] *= inv_akk;
X[k*n + j] *= inv_akk;
}
for (int i = k+1; i < n; ++i)
{
const double aik = A[i*n + k];
for (int j = 0; j < k+1; ++j)
X[i*n + j] -= X[k*n + j] * aik;
for (int j = k+1; j < n; ++j)
{
A[i*n + j] -= A[k*n + j] * aik;
X[i*n + j] -= X[k*n + j] * aik;
}
}
}
// Обратный ход метода Гаусса
void gauss_back_substitution(const int n, double * restrict A, double * restrict X)
{
// Идём с последней строки и вычитаем её из последующих
for (int k = n-1; k > 0; --k)
for (int i = 0; i < k; ++i)
{
const double aik = A[i*n + k];
for (int j = 0; j < n; ++j)
X[i*n + j] -= X[k*n + j] * aik;
}
}
// 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;
double maximum = -1.;
int max_i = 0, max_j = 0;
// Ищем максимальный элемент минора
for (int i = 0; i < n; ++i)
for (int j = 0; j < n; ++j)
{
const double aij = fabs(A[i*n + j]);
if (aij > maximum) {
maximum = aij;
max_i = i;
max_j = j;
}
}
// Проходимся по главным минорам
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;
}
for (int k = 0; k < n; ++k)
{
double inv_akk;
// Если максимальный по модулю элемент равен нулю, значит матрица вырождена
if (fabs(maximum) <= eps)
return SINGULAR;
@ -111,7 +74,40 @@ int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c)
}
}
gauss_inverse(n, k, A, X);
// Прямой ход метода Гаусса
inv_akk = 1./A[k*n + k];
for (int j = 0; j < k+1; ++j)
X[k*n + j] *= inv_akk;
for (int j = k+1; j < n; ++j)
{
A[k*n + j] *= inv_akk;
X[k*n + j] *= inv_akk;
}
maximum = -1.;
for (int i = k+1; i < n; ++i)
{
const double aik = A[i*n + k];
for (int j = 0; j < k+1; ++j)
X[i*n + j] -= X[k*n + j] * aik;
for (int j = k+1; j < n; ++j)
{
const double aij = A[i*n + j] - A[k*n + j] * aik;
if (fabs(aij) > maximum)
{
maximum = fabs(aij);
max_i = i;
max_j = j;
}
A[i*n + j] = aij;
X[i*n + j] -= X[k*n + j] * aik;
}
}
}
gauss_back_substitution(n, A, X);
@ -140,5 +136,17 @@ int t14_solve(int n, double * restrict A, double * restrict X, int * restrict c)
return 0;
}
// Обратный ход метода Гаусса
void gauss_back_substitution(const int n, double * restrict A, double * restrict X)
{
// Идём с последней строки и вычитаем её из последующих
for (int k = n-1; k > 0; --k)
for (int i = 0; i < k; ++i)
{
const double aik = A[i*n + k];
for (int j = 0; j < n; ++j)
X[i*n + j] -= X[k*n + j] * aik;
}
}

3
2025.04.04/dist/Krivoruchenko_SK/solve.h vendored
View file

@ -1,8 +1,7 @@
#ifndef SOLVE_H
#define SOLVE_H
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);
void gauss_back_substitution(const int n, double * restrict A, double * restrict X);
#endif

2
2025.04.18/01Ex/main.c
View file

@ -17,7 +17,7 @@ int main(int argc, char *argv[])
if (
!((argc == 4) &&
sscanf(argv[1], "%lf", &x_0) == 1 &&
sscanf(argv[2], "%d", &n) == 1 &&
((sscanf(argv[2], "%d", &n) == 1) && n > 0) &&
((name = argv[3]) && name))
) {
fprintf(stderr, "Usage: %s x_0 n filename\n", argv[0]);

16
2025.04.18/02Ex/solve.c
View file

@ -7,8 +7,7 @@
// the Newton interpolation polynomial
double t2_solve (const double x_0, const int n, const double * restrict X, double * restrict Y)
{
double last_x = X[n-1];
double last_y = Y[n-1];
double last_x, last_y;
double value = 0;
double start_value = 0;
@ -16,25 +15,22 @@ double t2_solve (const double x_0, const int n, const double * restrict X, doubl
for (int k = 0; k < n-1; ++k)
{
printf ("------- K = %d -------\n", k);
last_y = X[n-1];
for (int i = n-2; i >= k; --i)
{
const double x_i = X[i-k];
const double y_i = Y[i];
last_x = X[i+1];
if (fabs(last_x - x_i) < DBL_EPSILON)
return DBL_MAX;
// printf ("Y[%d] = (%lf - %lf) / (%lf - %lf)\n", i+1, last_y, y_i, last_x, x_i);
Y[i+1] = (last_y - y_i) / (last_x - x_i);
printf ("I = %d, f(x%d, ... , x%d) = %lf\n", i, i-k+1, i+2, Y[i]);
printf ("I = %d, f(x%d, ... , x%d) = %lf\n", i, i-k+1, i+2, Y[i+1]);
last_x = x_i;
last_y = y_i;
}
last_x = X[n-1];
last_y = Y[n-1];
}
last_x = X[0];
@ -46,8 +42,6 @@ double t2_solve (const double x_0, const int n, const double * restrict X, doubl
start_value *= (x_0 - last_x);
if (fabs(Y[i]) > DBL_EPSILON)
value += Y[i] * start_value;
printf ("i = %d, start_value = %lf, value = %lf, last_x = %lf\n", i, start_value, value, last_x);
printf ("Y[%d] = %lf\n", i, Y[i]);
last_x = X[i];
}

100000
2025.04.18/tests/big_linear.txt Normal file
View file

File diff suppressed because it is too large Load diff

13
2025.04.18/tests/sinx.txt Normal file
View file

@ -0,0 +1,13 @@
0 0
0.523598 0.5
1.047197 0.866025
1.570796 1
2.094395 0.866025
2.617993 0.5
3.141592 0
3.665191 -0.5
4.188790 -0.866025
4.712388 -1
5.235987 -0.866025
5.759586 -0.5
6.283185 0

6
2025.04.18/tests/x_2.txt Normal file
View file

@ -0,0 +1,6 @@
1 1
2 4
3 9
4 16
5 25
6 36