Compare commits

..

No commits in common. "main" and "recursive_polynomial_simplification" have entirely different histories.

11 changed files with 42 additions and 149 deletions

View File

@ -1,89 +1,3 @@
# Polynomial Interpolation # Polynomial Interpolation
ANSI C program which composes polynomial of n - 1 degree that passes through n dots. It presents it in Newton interpolation polynomial and monic form. ANSI C program which composes polynomial of n - 1 degree that passes through n dots.
## Interface
Application accepts as standart input decimal below 2147483647 `n` as number of dots, followed by n dots in format: `<x> (space) <y>` on each line, where `x` is an abscisse and `y` is an ordinate of single dot. Dot coordinates must fit [2.22507e-308;1.79769e+308] range by modulo.
Result will be printed to standart output in the following format:
Newton polynomial form:
$$f_0 - f_1*(x-x_0) + ... + f_n(x-x_0)*(x-x_1)*...*(x-x_{n-1})$$
Simplified coefficients array (starting from 0 upto n-1 power):
$$a_0 a_1 ... a_{n-1} a_n$$
Polynomial in monic form:
$$a_0 - a_1*x + ... + a_{n-1}*x^(n-2) + a_n*x^(n-1)$$
Where $f_i$ is a divided difference of $y_1,...,y_i$, $a_i$ are coefficients of resulting monic polynomial
## Data structure
- `n` is an `unsigned int` variable, that is used to input and store number of dots
- `x` is a pointer to array of `n` `double`s, that is used to store abscisses of dots
- `y` is a pointer to array of `n` `double`s, that is used to store ordinates of dots
- `coefficients` is a pointer to array of `n` `double`s, that is used to store coefficients of monic interpolation polynomial
- `i`, `j` are `int` variables, those are used in loops as iterators
- `tmp_polynomial` is a pointer to array of `n` `double`s, that is used to store coefficients of polynomial, resulting during simplification of interpolation polynomial summands.
## Example
Build and run application:
```bash
gcc main.c
./a.out
```
### Input/output
For input n = 3 and the following dots
```plain
1 5
2 3
4 8
```
Output is
```plain
Newton polynomial form:
5 - 2*(x-1) + 1.5*(x-1)*(x-2)
Simplified coefficients array (starting from 0 upto n-1 power):
10 -6.5 1.5
Polynomial in standart form:
10 - 6.5*x + 1.5*x^2
```
### Illustrations
#### Example 1
<img src="./img/console.png" />
<img src="./img/wolfram.png" />
<img src="./img/plot.png" />
#### Example 2
<img src="./img/console2.png" />
or
<img src="./img/console3.png" />
<img src="./img/wolfram2.png" />
<img src="./img/plot2.png" />

Binary file not shown.

Before

Width:  |  Height:  |  Size: 35 KiB

Binary file not shown.

Before

Width:  |  Height:  |  Size: 53 KiB

Binary file not shown.

Before

Width:  |  Height:  |  Size: 51 KiB

Binary file not shown.

Before

Width:  |  Height:  |  Size: 78 KiB

Binary file not shown.

Before

Width:  |  Height:  |  Size: 206 KiB

Binary file not shown.

Before

Width:  |  Height:  |  Size: 83 KiB

Binary file not shown.

Before

Width:  |  Height:  |  Size: 145 KiB

View File

@ -1,5 +1,4 @@
import sys import sys
import math
try: try:
n = int(sys.argv[1]) n = int(sys.argv[1])
@ -23,4 +22,4 @@ def f(x: int) -> int:
return res return res
for i in range(n): for i in range(n):
print(i, math.sin(i)) print(i, f(i))

76
main.c
View File

@ -16,16 +16,16 @@ double fabs(double x)
array x stands for x array x stands for x
number i stands for index of evaluated difference (from 0) number i stands for index of evaluated difference (from 0)
number d stands for order of difference (from 0) number d stands for order of difference (from 0)
example: https://en.wikipedia.org/wiki/Newton_polynomial#Examples */ example: https://shorturl.at/tBCPS */
double div_diff(double *y, double *x, unsigned i, unsigned 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]); return (y[i] - y[i - 1]) / (x[i] - x[i - d]);
} }
/* Evaluates divided differences of n values - array of some kind of derivatives with big enough dx /* Evaluates divided differences of n values - array of some kind of derivatives with big enough dx
Example: https://en.wikipedia.org/wiki/Newton_polynomial#Examples Example: https://shorturl.at/tBCPS
Warning: result is evaluated in `double *y` array */ Warning: result is evaluated in `double *y` array */
double *div_diff_es(double *x, double *y, unsigned 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 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 for (int j = n - 1; j >= i; j--) // evaluate from the end of array, decreacing number of step every repeation
@ -38,35 +38,26 @@ double *div_diff_es(double *x, double *y, unsigned n)
Coeficients of simplified polynomial computation Coeficients of simplified polynomial computation
*/ */
/* Simplifies Newton polynomial with `el_coef` array of divided differences, void simplify_polynomial(double *res, double *el_coef, double *x, unsigned int n)
and `x` as array of x coordinates of dots,
and `n` is number of elements of this sum */
void simplify_polynomial(double *res, double *el_coef, double *x, unsigned n)
{ {
double *tmp_polynomial // Temporary array for storage of coefficients of multiplication of (x-x_i) polynomial double *tmp_polynomial = (double *)malloc(sizeof(double) * n);
= (double *)malloc(sizeof(double) * n); tmp_polynomial[0] = 1;
for (int i = 1; i < n; i++)
tmp_polynomial[i] = 0;
tmp_polynomial[0] = 1; // Set polynomial to 1 to start multiplication with it
for (int i = 0; i < n; i++) // For each elemnt of sum for (int i = 0; i < n; i++)
if (el_coef[i])
{ {
if (i > 0) // Start multiplication from second element of sum if (i > 0)
mult_by_root(tmp_polynomial, x[i - 1], i - 1); mult_by_root(tmp_polynomial, x[i - 1], i - 1);
for (int j = 0; j <= i; j++) // For each cumputed coefficient of i'th polynomial of sum for (int j = 0; j <= i; j++)
res[j] += el_coef[i] * tmp_polynomial[j]; // Add it, multiplied with divided difference, to sum res[j] += el_coef[i] * tmp_polynomial[j];
} }
free(tmp_polynomial);
} }
/* `res` is an array of coefficients of polynomial, which is multiplied with (x - `root`) polynomial. void mult_by_root(double *res, double root, unsigned int step)
`power` is the power of `res` polynomial */
void mult_by_root(double *res, double root, unsigned power)
{ {
for (int j = power + 1; j >= 0; j--) for (int j = step + 1; j >= 0; j--)
res[j] = (j ? res[j - 1] : 0) - (root * res[j]); // coefficient is k_i-1 - root * k_i res[j] = (j ? res[j - 1] : 0) - (root * res[j]);
} }
/* /*
@ -74,7 +65,7 @@ void mult_by_root(double *res, double root, unsigned power)
*/ */
/* Prints interpolation polynomial in Newton notation */ /* Prints interpolation polynomial in Newton notation */
void print_newton_poly(double *f, double *x, unsigned n) void print_newton_poly(double *f, double *x, unsigned int n)
{ {
printf("Newton polynomial form:\n"); printf("Newton polynomial form:\n");
for (int i = 0; i < n; i++) for (int i = 0; i < n; i++)
@ -92,16 +83,16 @@ void print_newton_poly(double *f, double *x, unsigned n)
else if (f[i] < 0) // If it is the first summond and coefficient is below zero else if (f[i] < 0) // If it is the first summond and coefficient is below zero
printf("-"); printf("-");
printf("%g", fabs(f[i])); // Print coefficient without sign printf("%lf", fabs(f[i])); // Print coefficient without sign
for (int j = 0; j < i; j++) // For each (x-xi) bracket 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]) // If summond is not zero, print it
{ {
if (x[j] > 0) if (x[j] > 0)
printf("*(x-%g)", x[j]); printf("*(x-%lf)", x[j]);
else else
printf("*(x+%g)", -x[j]); printf("*(x+%lf)", -x[j]);
} }
else else
printf("*x"); printf("*x");
@ -114,18 +105,16 @@ void print_newton_poly(double *f, double *x, unsigned n)
printf("\n"); printf("\n");
} }
/* Returns inputed by user number of dots */ unsigned int insert_n()
unsigned insert_n()
{ {
printf("Insert number of dots: "); printf("Insert number of dots: ");
unsigned n = 0; unsigned int n = 0;
scanf("%u", &n); scanf("%u", &n);
return n; return n;
} }
/* Reads pairs of x'es and y'es of n dots to corresponding array */ void insert_coords(double *xes, double *yes, unsigned int n)
void insert_coords(double *xes, double *yes, unsigned n)
{ {
printf("Insert dots coordinates in the following format:\n<x> (space) <y>\nEach dot on new line\n"); printf("Insert dots coordinates in the following format:\n<x> (space) <y>\nEach dot on new line\n");
@ -139,22 +128,19 @@ void insert_coords(double *xes, double *yes, unsigned n)
} }
} }
/* Prints array of n doubles */ void print_array(double *arr, unsigned int n)
void print_array(double *arr, unsigned n)
{ {
printf("Simplified coefficients array (starting from 0 upto n-1 power):\n"); printf("Simplified coefficients array (starting from 0 upto n-1 power):\n");
for (int i = 0; i < n; i++) for (int i = 0; i < n; i++)
printf("%g ", arr[i]); printf("%lf ", arr[i]);
printf("\n"); printf("\n");
} }
/* Prints interpolation polynomial in standart form void print_poly(double *coef, unsigned int n)
e.g. a*x^2 + b*x + c */
void print_poly(double *coef, unsigned n)
{ {
printf("Polynomial in standart form:\n"); printf("Simplified polynom:\n");
for (int i = 0; i < n; i++) for (int i = 0; i < n; i++)
{ {
@ -165,10 +151,10 @@ void print_poly(double *coef, unsigned n)
printf("+ "); printf("+ ");
else else
printf("- "); printf("- ");
else if (coef[i] < 0) else
printf("-"); printf("-");
printf("%g", fabs(coef[i])); printf("%lf", fabs(coef[i]));
if (i > 0) if (i > 0)
printf("*x"); printf("*x");
if (i > 1) if (i > 1)
@ -199,8 +185,6 @@ int main()
print_newton_poly(f, x, n); print_newton_poly(f, x, n);
double *coefficients = (double *)malloc(sizeof(double) * n); double *coefficients = (double *)malloc(sizeof(double) * n);
for (unsigned i = 0; i < n; i++)
coefficients[i] = 0;
simplify_polynomial(coefficients, f, x, n); simplify_polynomial(coefficients, f, x, n);
@ -208,9 +192,5 @@ int main()
print_poly(coefficients, n); print_poly(coefficients, n);
free(x);
free(y);
free(coefficients);
return 0; return 0;
} }

View File

@ -13,24 +13,24 @@ double fabs(double x);
Business logic Business logic
*/ */
double div_diff(double *y, double *x, unsigned i, unsigned d); double div_diff(double *y, double *x, unsigned int i, unsigned int d);
double *div_diff_es(double *x, double *y, unsigned n); double *div_diff_es(double *x, double *y, unsigned int n);
/* /*
User interface User interface
*/ */
unsigned insert_n(); unsigned int insert_n();
void print_newton_poly(double *f, double *x, unsigned n); void print_newton_poly(double *f, double *x, unsigned int n);
void insert_coords(double *x, double *y, unsigned n); void insert_coords(double *x, double *y, unsigned int n);
void print_array(double *arr, unsigned n); void print_array(double *arr, unsigned int n);
void print_poly(double *coef, unsigned n); void print_poly(double *coef, unsigned int n);
/* /*
Coeficients of simplified polynomial computation Coeficients of simplified polynomial computation
*/ */
void simplify_polynomial(double *res, double *el_coef, double *x, unsigned n); void simplify_polynomial(double *res, double *el_coef, double *x, unsigned int n);
void mult_by_root(double *res, double root, unsigned step); void mult_by_root(double *res, double root, unsigned int step);
#endif #endif