149 lines
3.6 KiB
C++
149 lines
3.6 KiB
C++
#include <math.h>
|
|
#include <stdio.h>
|
|
int p = 0;
|
|
int m = 0, n = 0;
|
|
double A[5][10000], Z[10000];
|
|
|
|
double get_content(int row, int column) {
|
|
if (abs(row - column) > p / 2) {
|
|
return 0;
|
|
}
|
|
if (row < 0 || row >= n) {
|
|
return 0;
|
|
}
|
|
if (column < 0 || column >= n) {
|
|
return 0;
|
|
}
|
|
return A[(row - column + p / 2)][column];
|
|
}
|
|
|
|
int set_content(int row, int column, double val) {
|
|
if (abs(row - column) > p / 2) {
|
|
return 0;
|
|
}
|
|
if (row < 0 || row >= n) {
|
|
return 0;
|
|
}
|
|
if (column < 0 || column >= n) {
|
|
return 0;
|
|
}
|
|
A[(row - column + p / 2)][column] = val;
|
|
return 0;
|
|
}
|
|
|
|
int output() {
|
|
for (int i = 0; i < p; i++) {
|
|
for (int j = 0; j < n; j++) {
|
|
printf("%lf ", get_content(i, j));
|
|
}
|
|
printf("\n");
|
|
}
|
|
printf("\n");
|
|
return 0;
|
|
}
|
|
|
|
int outputZ() {
|
|
for (int i = 0; i < n; i++) {
|
|
printf("%lf ", Z[i]);
|
|
}
|
|
printf("\n");
|
|
return 0;
|
|
}
|
|
|
|
int main() {
|
|
scanf("%d", &p);
|
|
scanf("%d %d", &n, &m);
|
|
for (int i = 0; i < p / 2; i++) {
|
|
for (int j = p / 2 - i; j < n; j++) {
|
|
scanf("%lf", &A[i][j]);
|
|
}
|
|
}
|
|
for (int i = p / 2; i < p; i++) {
|
|
for (int j = 0; j < n - (i - p / 2); j++) {
|
|
scanf("%lf", &A[i][j]);
|
|
}
|
|
}
|
|
|
|
// printf("Readin\n");
|
|
// output();
|
|
|
|
for (int col = 0; col < n; col++) {
|
|
double diag_elem = get_content(col, col);
|
|
// make everything below the diag to 0 row by row
|
|
for (int row = col + 1; row < n; row++) {
|
|
double base = get_content(row, col);
|
|
if (base != 0) {
|
|
double coefficient = base / diag_elem;
|
|
set_content(row, col, -coefficient);
|
|
for (int l = col + 1; l < n; l++) {
|
|
set_content(row, l, get_content(row, l) - coefficient * get_content(col, l));
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// printf("L ready\n");
|
|
// output();
|
|
|
|
for (int row = 0; row < n; row++) {
|
|
double diag_elem = get_content(row, row);
|
|
set_content(row, row, 1 / diag_elem);
|
|
for (int col = row + 1; col < n; col++) {
|
|
set_content(row, col, -get_content(row, col) / diag_elem);
|
|
}
|
|
}
|
|
|
|
// printf("All ready\n");
|
|
// output();
|
|
|
|
for (int i = 0; i < m; i++) {
|
|
// read in one col of z
|
|
for (int j = 0; j < n; j++) {
|
|
scanf("%lf", &Z[j]);
|
|
}
|
|
|
|
// printf("Readin\n");
|
|
// outputZ();
|
|
|
|
// calculate Z = L-1 * Z
|
|
// do this column by column
|
|
for (int j = 0; j < n; j++) {
|
|
// diag of U is always 1, z[j] * U[j][j] = z[j], so k doesn't have to equal to j
|
|
for (int k = j + 1; k < j + 3 && k < n; k++) {
|
|
Z[k] += Z[j] * get_content(k, j);
|
|
}
|
|
}
|
|
|
|
// printf("L-1 * Z\n");
|
|
// outputZ();
|
|
|
|
// calculate Z = D-1 * Z
|
|
for (int j = 0; j < n; j++) {
|
|
Z[j] = Z[j] * get_content(j, j);
|
|
}
|
|
|
|
// printf("D-1 L-1 * Z\n");
|
|
// outputZ();
|
|
|
|
// calculate Z = U-1 * Z
|
|
// do this column by column
|
|
// printf("Starting U-1\n");
|
|
for (int j = n - 1; j >= 0; j--) {
|
|
// diag of U is always 1, z[j] * U[j][j] = z[j], so k doesn't have to equal to j
|
|
for (int k = j - 1; k > j - 3 && k >= 0; k--) {
|
|
Z[k] += Z[j] * get_content(k, j);
|
|
}
|
|
// outputZ();
|
|
}
|
|
|
|
// printf("U-1 D-1 L-1 * Z\n");
|
|
// outputZ();
|
|
|
|
// output
|
|
for (int j = 0; j < n; j++) {
|
|
printf("%.4lf ", Z[j]);
|
|
}
|
|
printf("\n");
|
|
}
|
|
return 0;
|
|
} |