From 8933e4b2208e7e51f9d74320e4eabd9d25458649 Mon Sep 17 00:00:00 2001 From: dm1sh Date: Sat, 30 Oct 2021 13:12:39 +0300 Subject: [PATCH 1/6] Fixed resize function, separated array index convertion to convert_pos function, added separated get function --- main.c | 40 +++++++++++++++++++++++-------------- polynominal_interpolation.h | 2 ++ 2 files changed, 27 insertions(+), 15 deletions(-) diff --git a/main.c b/main.c index e007adf..cd3423f 100644 --- a/main.c +++ b/main.c @@ -41,12 +41,12 @@ arr *resize(arr *a, int new_size) double *new_p = (double *)malloc(sizeof(double) * new_size); for (int i = 0; i < min(new_size, a->size); i++) - new_p[i] = a->p[i]; + new_p[i] = get(a, i); free(a->p); for (int i = a->size; i < new_size; i++) - new_p = 0; + new_p[i] = 0; a->p = new_p; a->size = new_size; @@ -54,19 +54,29 @@ arr *resize(arr *a, int new_size) return a; } +int convert_pos(int size, int pos) +{ + pos = pos % size; + if (pos < 0) + pos = size + pos; + + return pos; +} + void insert(arr *a, int pos, double val) { - pos = pos % a->size; - if (pos < 0) - pos = a->size + pos; + a->p[convert_pos(a->size, pos)] = val; +} - a->p[pos] = val; +double get(arr *a, int pos) +{ + return a->p[convert_pos(a->size, pos)]; } arr *add(arr *a, arr *b) { for (int i = 0; i < a->size; i++) - insert(a, i, a->p[i] + b->p[i]); + insert(a, i, get(a, i) + get(b, i)); return a; } @@ -76,7 +86,7 @@ arr *mult(arr *a, double mul) arr *res = init(a->size); for (int i = 0; i < a->size; i++) - insert(res, i, a->p[i] * mul); + insert(res, i, get(a, i) * mul); return res; } @@ -90,7 +100,7 @@ void printa(arr *a) printf("\n"); for (int i = 0; i < a->size; i++) - printf("%5.2f ", a->p[i]); + printf("%5.2f ", get(a, i)); printf("\n"); } @@ -101,7 +111,7 @@ arr *arr_without_el(arr *a, int ex_pos) { if (i == ex_pos) continue; - insert(res, pos, a->p[i]); + insert(res, pos, get(a, i)); pos++; } @@ -112,7 +122,7 @@ arr *reverse(arr *a) { arr *res = init(a->size); for (int i = 0; i < a->size; i++) - insert(res, i, a->p[a->size - 1 - i]); + insert(res, i, get(a, a->size - 1 - i)); return res; } @@ -161,7 +171,7 @@ int mult_by_index(arr *a, int *coords, int n) { double res = 1; for (int i = 0; i < n; i++) - res = res * a->p[coords[i]]; + res = res * get(a, coords[i]); return res; } @@ -173,7 +183,7 @@ int sum_of_mult_of_n_combinations(arr *a, int n) if (a->size == 1) { - return a->p[0]; + return get(a, 0); } double acc = 0; @@ -197,7 +207,7 @@ double compose_denominator(arr *a, int pos) if (i == pos) continue; - res = res * (a->p[pos] - a->p[i]); + res = res * (get(a, pos) - get(a, i)); } return res; } @@ -211,7 +221,7 @@ arr *compose_interpolation_polynomial(arr *xes, arr *ys) { int minus = !(xes->size % 2); double denominator = compose_denominator(xes, j); - double multiplicator = ys->p[j]; + double multiplicator = get(ys, j); arr *xis = arr_without_el(xes, j); diff --git a/polynominal_interpolation.h b/polynominal_interpolation.h index 983cf2f..f702dcf 100644 --- a/polynominal_interpolation.h +++ b/polynominal_interpolation.h @@ -20,7 +20,9 @@ typedef struct arr *init(int n); arr *resize(arr *a, int new_size); +int convert_pos(int size, int pos); void insert(arr *a, int pos, double val); +double get(arr *a, int pos); arr *add(arr *a, arr *b); arr *mult(arr *a, double mul); void printa(arr *a); From 44281e14b31640e5ce8efe2503fdb0ec8b835af3 Mon Sep 17 00:00:00 2001 From: dm1sh Date: Sat, 30 Oct 2021 13:15:40 +0300 Subject: [PATCH 2/6] Updated readme --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index 973dcb6..d78d042 100644 --- a/README.md +++ b/README.md @@ -1,3 +1,3 @@ # Polynomial Interpolation -ANSI C program which composes polynomial of n - 1 degree for n dots. +ANSI C program which composes polynomial of n - 1 degree that passes through n dots. From 49619c27e0b5d43f91d877f2e4dc5e109d6ae6f9 Mon Sep 17 00:00:00 2001 From: dm1sh Date: Sat, 30 Oct 2021 15:35:40 +0300 Subject: [PATCH 3/6] First Newton interpolation polynomial implementation. No dots input, no polynomial simplification --- main.c | 295 ++++-------------------------------- polynominal_interpolation.h | 36 +---- 2 files changed, 33 insertions(+), 298 deletions(-) diff --git a/main.c b/main.c index cd3423f..5634325 100644 --- a/main.c +++ b/main.c @@ -1,281 +1,46 @@ -#include -#include - #include "./polynominal_interpolation.h" -/* - Utils -*/ - -int min(int a, int b) +/* Divided difference is evaluated for: + array y stands for f(x) + array x stands for x + number i stands for index of evaluated difference (from 0) + number d stands for order of difference (from 0) + example: https://shorturl.at/tBCPS */ +double div_diff(double *y, double *x, int i, int d) { - return (a + b - abs(a - b)) / 2; + return (y[i] - y[i - 1]) / (x[i] - x[i - d]); } -int max(int a, int b) +/* Evaluates divided differences of n values - array of some kind of derivatives with big enough dx + Example: https://shorturl.at/tBCPS + Warning: result is evaluated in `double *y` array */ +double *div_diff_es(double *x, double *y, int n) { - return (a + b + abs(a - b)) / 2; + for (int i = 1; i < n; i++) // first element remains unchanged + for (int j = n - 1; j >= i; j--) // evaluate from the end of array, decreacing number of step every repeation + y[j] = div_diff(y, x, j, i); + + return y; } -/* - Array utils -*/ - -arr *init(int n) +/* Prints interpolation polynomial in Newton notation */ +void print_newton_poly(double *f, double *x, int n) { - arr *a = (arr *)malloc(sizeof(arr)); + for (int i = 0; i < n; i++) + { + printf("(%lf)", f[i]); + for (int j = 0; j < i; j++) + printf("*(x-(%lf))", x[j]); - a->size = n; - - a->p = (double *)malloc(sizeof(double) * n); - for (int i = 0; i < n; i++) - insert(a, i, 0); - - return a; -} - -arr *resize(arr *a, int new_size) -{ - if (a->size == new_size) - return a; - - double *new_p = (double *)malloc(sizeof(double) * new_size); - for (int i = 0; i < min(new_size, a->size); i++) - new_p[i] = get(a, i); - - free(a->p); - - for (int i = a->size; i < new_size; i++) - new_p[i] = 0; - - a->p = new_p; - a->size = new_size; - - return a; -} - -int convert_pos(int size, int pos) -{ - pos = pos % size; - if (pos < 0) - pos = size + pos; - - return pos; -} - -void insert(arr *a, int pos, double val) -{ - a->p[convert_pos(a->size, pos)] = val; -} - -double get(arr *a, int pos) -{ - return a->p[convert_pos(a->size, pos)]; -} - -arr *add(arr *a, arr *b) -{ - for (int i = 0; i < a->size; i++) - insert(a, i, get(a, i) + get(b, i)); - - return a; -} - -arr *mult(arr *a, double mul) -{ - arr *res = init(a->size); - - for (int i = 0; i < a->size; i++) - insert(res, i, get(a, i) * mul); - - return res; -} - -void printa(arr *a) -{ - printf("Array of size %d:\n", a->size); - - for (int i = 0; i < a->size; i++) - printf("%5d ", i + 1); - printf("\n"); - - for (int i = 0; i < a->size; i++) - printf("%5.2f ", get(a, i)); - printf("\n"); -} - -arr *arr_without_el(arr *a, int ex_pos) -{ - arr *res = init(a->size - 1); - for (int i = 0, pos = 0; i < a->size; i++) - { - if (i == ex_pos) - continue; - insert(res, pos, get(a, i)); - pos++; - } - - return res; -} - -arr *reverse(arr *a) -{ - arr *res = init(a->size); - for (int i = 0; i < a->size; i++) - insert(res, i, get(a, a->size - 1 - i)); - - return res; -} - -void free_arr(arr *a) -{ - free(a->p); - free(a); -} - -/* - Business logic -*/ - -int has_comb(int *arr, int n, int k) -{ - if (n == k) - return 0; - - int pos = k - 1; - - if (arr[pos] == n - 1) - { - if (k == 1) - return 0; - - while ((pos > 0) && arr[pos] == n - 1) - { - pos--; - arr[pos]++; - } - - for (int i = pos + 1; i < k; i++) - arr[i] = arr[i - 1] + 1; - - if (arr[0] > n - k) - return 0; - } - else - arr[pos]++; - - return 1; -} - -int mult_by_index(arr *a, int *coords, int n) -{ - double res = 1; - for (int i = 0; i < n; i++) - res = res * get(a, coords[i]); - - return res; -} - -int sum_of_mult_of_n_combinations(arr *a, int n) -{ - if (n == 0) - return 1; - - if (a->size == 1) - { - return get(a, 0); - } - - double acc = 0; - - int coords[n]; - for (int i = 0; i < n; i++) - coords[i] = i; - - acc += mult_by_index(a, coords, n); - while (has_comb(coords, a->size, n)) - acc += mult_by_index(a, coords, n); - - return acc; -} - -double compose_denominator(arr *a, int pos) -{ - double res = 1; - for (int i = 0; i < a->size; i++) - { - if (i == pos) - continue; - - res = res * (get(a, pos) - get(a, i)); - } - return res; -} - -arr *compose_interpolation_polynomial(arr *xes, arr *ys) -{ - arr *res = init(xes->size); - - arr *jcoef = init(xes->size); - for (int j = 0; j < xes->size; j++) - { - int minus = !(xes->size % 2); - double denominator = compose_denominator(xes, j); - double multiplicator = get(ys, j); - - arr *xis = arr_without_el(xes, j); - - for (int i = 0; i < xes->size; i++) - { - double k_sum = sum_of_mult_of_n_combinations(xis, xes->size - 1 - i); - insert(jcoef, i, (minus ? -1 : 1) * (multiplicator * k_sum) / denominator); - minus = !minus; - } - - res = add(res, jcoef); - - free_arr(xis); - } - - free_arr(jcoef); - - return res; + if (i != n - 1) + printf("+"); + } } int main() { - printf("Insert number of dots: "); - int n = 0; - scanf("%d", &n); + double x[] = {0, 1, 2, 3}, + y[] = {-2, -5, 0, -4}; - printf("Insert dots coordinates in the following format:\n (space) \nEach dot on new line\n"); - - arr *xes = init(n); - arr *ys = init(n); - - for (int i = 0; i < n; i++) - { - double x, y; - scanf("%lf %lf", &x, &y); - - insert(xes, i, x); - insert(ys, i, y); - } - - printf("Inserted the following doths:\n"); - printa(xes); - printa(ys); - - arr *res = compose_interpolation_polynomial(xes, ys); - - printf("Resulting polynomial will have such coeficients:\n"); - arr *reversed = reverse(res); - printa(reversed); - - free_arr(reversed); - free_arr(res); - free_arr(xes); - free_arr(ys); - - return 0; + print_newton_poly(div_diff_es(x, y, 4), x, 4); } \ No newline at end of file diff --git a/polynominal_interpolation.h b/polynominal_interpolation.h index f702dcf..729d043 100644 --- a/polynominal_interpolation.h +++ b/polynominal_interpolation.h @@ -1,42 +1,12 @@ #ifndef POLYNOMIAL_INTERPOLATION_H #define POLYNOMIAL_INTERPOLATION_H -/* - Utils -*/ - -int min(int a, int b); -int max(int a, int b); - -/* - Array utils -*/ - -typedef struct -{ - int size; - double *p; -} arr; - -arr *init(int n); -arr *resize(arr *a, int new_size); -int convert_pos(int size, int pos); -void insert(arr *a, int pos, double val); -double get(arr *a, int pos); -arr *add(arr *a, arr *b); -arr *mult(arr *a, double mul); -void printa(arr *a); -arr *arr_without_el(arr *a, int ex_pos); -arr *reverse(arr *a); - /* Business logic */ -int has_comb(int *arr, int n, int k); -int mult_by_index(arr *a, int *coords, int n); -int sum_of_mult_of_n_combinations(arr *a, int n); -double compose_denominator(arr *a, int pos); -arr *compose_interpolation_polynomial(arr *xes, arr *ys); +double div_diff(double *y, double *x, int i, int d); +double *div_diff_es(double *x, double *y, int n); +void print_newton_poly(double *f, double *x, int n); #endif \ No newline at end of file From c83cc2894c09a8bda139032de16d37b0820aad2a Mon Sep 17 00:00:00 2001 From: dm1sh Date: Sat, 30 Oct 2021 16:16:56 +0300 Subject: [PATCH 4/6] Added polynomial input interface --- main.c | 33 ++++++++++++++++++++++++++++++--- polynominal_interpolation.h | 10 ++++++++++ 2 files changed, 40 insertions(+), 3 deletions(-) diff --git a/main.c b/main.c index 5634325..d92ea80 100644 --- a/main.c +++ b/main.c @@ -37,10 +37,37 @@ void print_newton_poly(double *f, double *x, int n) } } +unsigned int insert_n() +{ + printf("Insert number of dots: "); + unsigned int n = 0; + scanf("%u", &n); + + return n; +} + +void insert_coords(double *xes, double *yes, unsigned int n) +{ + printf("Insert dots coordinates in the following format:\n (space) \nEach dot on new line\n"); + + for (int i = 0; i < n; i++) + { + double x, y; + scanf("%lf %lf", &x, &y); + + xes[i] = x; + yes[i] = y; + } +} + int main() { - double x[] = {0, 1, 2, 3}, - y[] = {-2, -5, 0, -4}; + unsigned n = insert_n(); - print_newton_poly(div_diff_es(x, y, 4), x, 4); + double *x = (double *)malloc(sizeof(double) * n), + *y = (double *)malloc(sizeof(double) * n); + + insert_coords(x, y, n); + + print_newton_poly(div_diff_es(x, y, n), x, n); } \ No newline at end of file diff --git a/polynominal_interpolation.h b/polynominal_interpolation.h index 729d043..4ec01bd 100644 --- a/polynominal_interpolation.h +++ b/polynominal_interpolation.h @@ -1,12 +1,22 @@ #ifndef POLYNOMIAL_INTERPOLATION_H #define POLYNOMIAL_INTERPOLATION_H +#include +#include + /* Business logic */ double div_diff(double *y, double *x, int i, int d); double *div_diff_es(double *x, double *y, int n); + +/* + User interface +*/ + +unsigned int insert_n(); void print_newton_poly(double *f, double *x, int n); +void insert_coords(double *x, double *y, unsigned int n); #endif \ No newline at end of file From 95a1c4cb8194e92e0c5105ce7e641ce354a031b3 Mon Sep 17 00:00:00 2001 From: dm1sh Date: Sun, 31 Oct 2021 02:06:56 +0300 Subject: [PATCH 5/6] Added simpified polynomial coefficients computation --- main.c | 104 ++++++++++++++++++++++++++++++++++-- polynominal_interpolation.h | 15 ++++-- 2 files changed, 112 insertions(+), 7 deletions(-) diff --git a/main.c b/main.c index d92ea80..0ea1d65 100644 --- a/main.c +++ b/main.c @@ -1,12 +1,16 @@ #include "./polynominal_interpolation.h" +/* + Newton interpolation polynomial +*/ + /* Divided difference is evaluated for: array y stands for f(x) array x stands for x number i stands for index of evaluated difference (from 0) number d stands for order of difference (from 0) example: https://shorturl.at/tBCPS */ -double div_diff(double *y, double *x, int i, int d) +double div_diff(double *y, double *x, unsigned int i, unsigned int d) { return (y[i] - y[i - 1]) / (x[i] - x[i - d]); } @@ -14,7 +18,7 @@ double div_diff(double *y, double *x, int i, int d) /* Evaluates divided differences of n values - array of some kind of derivatives with big enough dx Example: https://shorturl.at/tBCPS Warning: result is evaluated in `double *y` array */ -double *div_diff_es(double *x, double *y, int n) +double *div_diff_es(double *x, double *y, unsigned int n) { for (int i = 1; i < n; i++) // first element remains unchanged for (int j = n - 1; j >= i; j--) // evaluate from the end of array, decreacing number of step every repeation @@ -23,8 +27,78 @@ double *div_diff_es(double *x, double *y, int n) return y; } +/* + Coeficients of simplified polynomial computation +*/ + +void simplify_polynomial(double *res, double *rev_el_coef, double *x, unsigned int n) +{ + for (int i = 0; i < n; i++) + if (rev_el_coef[i]) + for (int j = 0; j <= i; j++) + res[i - j] += (j % 2 ? -1 : 1) * rev_el_coef[i] * compute_sum_of_multiplications_of_k(x, j, i); +} + +double compute_sum_of_multiplications_of_k(double *arr, unsigned int k, unsigned int n) +{ + if (k == 0) + return 1; + + if (k == 1 && n == 1) + return arr[0]; + + unsigned int *selected = (unsigned int *)malloc(sizeof(unsigned int) * k); // Indexes of selected for multiplication elements + + int i = 0, // index of `arr` array + j = 0; // index of `selected` array + + double sum = 0; + + while (j >= 0) + { + if (i <= (n + (j - k))) + { + selected[j] = i; + + if (j == k - 1) + { + sum += mult_by_indexes(arr, selected, k); + i++; + } + else + { + i = selected[j] + 1; + j++; + } + } + else + { + j--; + if (j >= 0) + i = selected[j] + 1; + } + } + + free(selected); + + return sum; +} + +double mult_by_indexes(double *arr, unsigned int *indexes, unsigned int size) +{ + double res = 1; + for (int i = 0; i < size; i++) + res *= arr[indexes[i]]; + + return res; +} + +/* + User Interface +*/ + /* Prints interpolation polynomial in Newton notation */ -void print_newton_poly(double *f, double *x, int n) +void print_newton_poly(double *f, double *x, unsigned int n) { for (int i = 0; i < n; i++) { @@ -60,6 +134,18 @@ void insert_coords(double *xes, double *yes, unsigned int n) } } +void print_array(double *arr, unsigned int n) +{ + for (int i = 0; i < n; i++) + printf("%lf ", arr[i]); + + printf("\n"); +} + +/* + Main +*/ + int main() { unsigned n = insert_n(); @@ -69,5 +155,15 @@ int main() insert_coords(x, y, n); - print_newton_poly(div_diff_es(x, y, n), x, n); + double *f = div_diff_es(x, y, n); + + print_newton_poly(f, x, n); + + double *coeficients = (double *)malloc(sizeof(double) * n); + + simplify_polynomial(coeficients, f, x, n); + + print_array(coeficients, n); + + return 0; } \ No newline at end of file diff --git a/polynominal_interpolation.h b/polynominal_interpolation.h index 4ec01bd..2815c53 100644 --- a/polynominal_interpolation.h +++ b/polynominal_interpolation.h @@ -8,15 +8,24 @@ Business logic */ -double div_diff(double *y, double *x, int i, int d); -double *div_diff_es(double *x, double *y, int n); +double div_diff(double *y, double *x, unsigned int i, unsigned int d); +double *div_diff_es(double *x, double *y, unsigned int n); /* User interface */ unsigned int insert_n(); -void print_newton_poly(double *f, double *x, int n); +void print_newton_poly(double *f, double *x, unsigned int n); void insert_coords(double *x, double *y, unsigned int n); +void print_array(double *arr, unsigned int n); + +/* + Coeficients of simplified polynomial computation +*/ + +void simplify_polynomial(double *res, double *rev_el_coef, double *x, unsigned int n); +double compute_sum_of_multiplications_of_k(double *x, unsigned int k, unsigned int n); +double mult_by_indexes(double *arr, unsigned int *indexes, unsigned int size); #endif \ No newline at end of file From 9155c87b65b58712b45c75a82ed832aec318f4e1 Mon Sep 17 00:00:00 2001 From: dm1sh Date: Sun, 31 Oct 2021 03:28:40 +0300 Subject: [PATCH 6/6] Updated (shortened) Newton polynomial printing, added simplified polynomial printing --- main.c | 82 +++++++++++++++++++++++++++++++++---- polynominal_interpolation.h | 6 +++ 2 files changed, 80 insertions(+), 8 deletions(-) diff --git a/main.c b/main.c index 0ea1d65..0f13dc1 100644 --- a/main.c +++ b/main.c @@ -1,5 +1,12 @@ #include "./polynominal_interpolation.h" +/* Utils */ + +double fabs(double x) +{ + return x > 0 ? x : -x; +} + /* Newton interpolation polynomial */ @@ -100,15 +107,42 @@ double mult_by_indexes(double *arr, unsigned int *indexes, unsigned int size) /* Prints interpolation polynomial in Newton notation */ void print_newton_poly(double *f, double *x, unsigned int n) { + printf("Newton polynomial form:\n"); for (int i = 0; i < n; i++) { - printf("(%lf)", f[i]); - for (int j = 0; j < i; j++) - printf("*(x-(%lf))", x[j]); + if (f[i]) // If coefficient != 0 + { + /* Coefficient sign and sum symbol */ + if (i > 0 && f[i - 1]) // If it's not the first summond + { + if (f[i] > 0) + printf("+ "); + else + printf("- "); + } + else if (f[i] < 0) // If it is the first summond and coefficient is below zero + printf("-"); - if (i != n - 1) - printf("+"); + printf("%lf", fabs(f[i])); // Print coefficient without sign + + for (int j = 0; j < i; j++) // For each (x-xi) bracket + { + if (x[j]) // If summond is not zero, print it + { + if (x[j] > 0) + printf("*(x-%lf)", x[j]); + else + printf("*(x+%lf)", -x[j]); + } + else + printf("*x"); + } + + printf(" "); + } } + + printf("\n"); } unsigned int insert_n() @@ -136,12 +170,42 @@ void insert_coords(double *xes, double *yes, unsigned int n) void print_array(double *arr, unsigned int n) { + printf("Simplified coefficients array (starting from 0 upto n-1 power):\n"); + for (int i = 0; i < n; i++) printf("%lf ", arr[i]); printf("\n"); } +void print_poly(double *coef, unsigned int n) +{ + printf("Simplified polynom:\n"); + + for (int i = 0; i < n; i++) + { + if (coef[i]) + { + if (i > 0 && coef[i - 1]) + if (coef[i] > 0) + printf("+ "); + else + printf("- "); + else + printf("-"); + + printf("%lf", fabs(coef[i])); + if (i > 0) + printf("*x"); + if (i > 1) + printf("^%d ", i); + else printf(" "); + } + } + + printf("\n"); +} + /* Main */ @@ -159,11 +223,13 @@ int main() print_newton_poly(f, x, n); - double *coeficients = (double *)malloc(sizeof(double) * n); + double *coefficients = (double *)malloc(sizeof(double) * n); - simplify_polynomial(coeficients, f, x, n); + simplify_polynomial(coefficients, f, x, n); - print_array(coeficients, n); + print_array(coefficients, n); + + print_poly(coefficients, n); return 0; } \ No newline at end of file diff --git a/polynominal_interpolation.h b/polynominal_interpolation.h index 2815c53..1cca3ec 100644 --- a/polynominal_interpolation.h +++ b/polynominal_interpolation.h @@ -4,6 +4,11 @@ #include #include +/* + Utils +*/ +double fabs(double x); + /* Business logic */ @@ -19,6 +24,7 @@ unsigned int insert_n(); void print_newton_poly(double *f, double *x, unsigned int n); void insert_coords(double *x, double *y, unsigned int n); void print_array(double *arr, unsigned int n); +void print_poly(double *coef, unsigned int n); /* Coeficients of simplified polynomial computation