#include #include 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("%.4lf ", 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 < col + 3; 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"); outputZ(); } return 0; }