Last active
July 13, 2022 17:45
-
-
Save kala13x/6341e241d9b11f66ee5bc65484b7e557 to your computer and use it in GitHub Desktop.
Lagrange polynomial and coefficient calculations
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
/* | |
* lagint.c | |
* 2019-2020 Sun Dro (f4tb0y@protonmail.com) | |
* | |
* Calculate Lagrange interpolated polynomials and | |
* displays the coefficient of each polynomials | |
*/ | |
#include <stdio.h> | |
#include <errno.h> | |
#include <string.h> | |
#include <stdlib.h> | |
#include <unistd.h> | |
#include <signal.h> | |
#include <fcntl.h> | |
#include <math.h> | |
#include <sys/types.h> | |
#include <sys/wait.h> | |
#include <sys/stat.h> | |
#include <sys/file.h> | |
/* | |
Example input: | |
6.0,8.0,1.0,12.0,5.0,5.0,9.0,7.0,8.0,2.0,7.0,4.0,4.0,6.0,5.5,2.0 | |
4.0,9.0,6.0,3.0,2.0,9.0,7.0,3.0,9.0,2.0,1.0,6.0,8.0,7.0,5.8,4.0 | |
1.0,8.0,7.0,7.0,5.0,9.0,2.0,3.0,9.0,1.0,17.0,0.0,3.0,3.0,8.0,5.0 | |
9.0,10.0,7.0,7.0,8.0,9.0,5.0,3.0,2.0,1.0,3.0,1.0,4.0,8.0,6.0,7.0 | |
1.0,1.0,3.0,4.0,5.0,8.0,7.0,3.0,15.0,1.0,4.0,4.0,8.0,10.0,9.0,7.0 | |
8.0,7.0,1.0,0.0,18.0,8.0,2.0,9.0,4.0,4.0,3.0,6.0,6.0,8.0,10.0,9.0 | |
5.0,8.0,8.0,3.0,7.0,7.0,3.0,2.0,15.0,1.0,10.0,6.0,6.0,6.0,4.1,6.0 | |
1.0,9.0,3.0,3.0,4.0,4.0,9.0,4.0,6.0,2.0,7.0,1.0,5.0,3.0,3.5,3.0 | |
*/ | |
#define BUFF_MAX 8192 | |
#define SIZE_ARR 32 | |
struct Node { | |
double coeff; | |
int power; | |
struct Node* next; | |
}; | |
struct polynomial { | |
struct Node *p[10]; | |
}; | |
void exit_failure() | |
{ | |
/* Terminate every process before exit */ | |
killpg(getpgrp(), SIGUSR2); | |
exit(EXIT_FAILURE); | |
} | |
int read_file(int fd, char *buffer, size_t size) | |
{ | |
lseek(fd, SEEK_SET, 0); | |
int bytes = read(fd, buffer, size - 1); | |
int length = (bytes > 0) ? bytes : 0; | |
if (length == 0 && errno == 0) | |
return read_file(fd, buffer, size); | |
buffer[length] = '\0'; | |
return length; | |
} | |
int get_line_count(char *buffer, size_t length) | |
{ | |
size_t i, lastp = 0; | |
int count = 0; | |
/* Count new line characters */ | |
for (i = 0; i < length; i++) | |
{ | |
if (buffer[i] == '\n') | |
{ | |
lastp = i; | |
count++; | |
} | |
} | |
/* Check if last line doesnot contain new line character */ | |
if ((i - lastp > 1) && buffer[i] != '\n') count++; | |
return count; | |
} | |
int get_line(const char *path, int num, char *output, size_t size) | |
{ | |
int fd = open(path, O_RDONLY, 644); | |
if (fd < 0) | |
{ | |
fprintf(stderr, "Can not open file: %s\n", strerror(errno)); | |
return 0; | |
} | |
char buffer[BUFF_MAX]; | |
struct flock lock; | |
lock.l_type = F_RDLCK; | |
lock.l_whence = SEEK_SET; | |
lock.l_start = 0; | |
lock.l_len = 0; | |
if (fcntl(fd, F_SETLKW, &lock) == -1) | |
{ | |
fprintf(stderr, "Can not lock the file: %s\n", strerror(errno)); | |
exit_failure(); | |
} | |
int bytes = read_file(fd, buffer, sizeof(buffer)); | |
if (bytes <= 0) | |
{ | |
fprintf(stderr, "Can not read file: %s\n", strerror(errno)); | |
return 0; | |
} | |
/* Also releases all locks */ | |
close(fd); | |
int i, len = 0, line_num = 1; | |
for (i = 0; i < bytes; i++) | |
{ | |
if (num == line_num) | |
{ | |
output[len++] = buffer[i]; | |
if (len == size) break; | |
} | |
if (buffer[i] == '\n') line_num++; | |
} | |
output[len] = '\0'; | |
return len; | |
} | |
int put_line(const char *path, int num, const char *line, size_t length) | |
{ | |
int fd = open(path, O_RDWR, 644); | |
if (fd < 0) | |
{ | |
fprintf(stderr, "Can not open file: %s\n", strerror(errno)); | |
return 1; | |
} | |
char buffer[BUFF_MAX]; | |
struct flock lock; | |
lock.l_type = F_WRLCK; | |
lock.l_whence = SEEK_SET; | |
lock.l_start = 0; | |
lock.l_len = 0; | |
if (fcntl(fd, F_SETLKW, &lock) == -1) | |
{ | |
fprintf(stderr, "Can not lock the file: %s\n", strerror(errno)); | |
exit_failure(); | |
} | |
int bytes = read_file(fd, buffer, sizeof(buffer)); | |
if (bytes <= 0) | |
{ | |
fprintf(stderr, "Can not read file: %s\n", strerror(errno)); | |
return 0; | |
} | |
char front_buffer[BUFF_MAX]; | |
char back_buffer[BUFF_MAX]; | |
int front_len = 0; | |
int front_pos = 0; | |
int back_len = 0; | |
int i, line_num = 1; | |
for (i = 0; i < bytes; i++) | |
{ | |
if (line_num < num) | |
{ | |
front_buffer[front_len++] = buffer[i]; | |
front_pos = i; | |
} | |
else if (line_num > num) | |
{ | |
back_buffer[back_len++] = buffer[i]; | |
} | |
if (buffer[i] == '\n') line_num++; | |
} | |
int full_len = front_len + back_len + length; | |
int front_offset = 0, back_offset = 0, line_offset = 0; | |
int insert_done = 0; | |
for (i = 0; i < full_len; i++) | |
{ | |
if (front_pos && i <= front_pos) | |
{ | |
buffer[i] = front_buffer[front_offset++]; | |
} | |
else if (!insert_done) | |
{ | |
buffer[i] = line[line_offset++]; | |
if (line_offset == length) insert_done = 1; | |
} | |
else | |
{ | |
buffer[i] = back_buffer[back_offset++]; | |
} | |
} | |
buffer[full_len] = '\0'; | |
/* | |
It is necessary to truncate the file manually because the O_TRUNC flag | |
will damage the result of read() function and due to the synchronization | |
the reading and writing must be done in one call of the open() function. | |
*/ | |
ftruncate(fd, 0); | |
/* Move offset in the begining of the file */ | |
lseek(fd, SEEK_SET, 0); | |
bytes = write(fd, buffer, full_len); | |
close(fd); /* Also releases all locks */ | |
return bytes; | |
} | |
int get_values(const char *line, int len, float (*x)[SIZE_ARR], float (*y)[SIZE_ARR]) | |
{ | |
char value[16]; | |
int i, offset = 0, isx = 0; | |
int ix = 0, iy = 0; | |
for (i = 0; i < len; i++) | |
{ | |
value[offset] = line[i]; | |
if (value[offset] == ',' || value[offset] == '\n') | |
{ | |
value[offset] = '\0'; | |
offset = 0; | |
isx = !isx; | |
if (isx) (*x)[ix++] = atof(value); | |
else (*y)[iy++] = atof(value); | |
if (ix >= SIZE_ARR || iy >= SIZE_ARR) | |
{ | |
fprintf(stderr, "Too big input\n"); | |
return 0; | |
} | |
continue; | |
} | |
offset++; | |
} | |
return ix; | |
} | |
/* Function add a new node at the end of list */ | |
struct Node* add_node(struct Node* start, double coeff, int power) | |
{ | |
/* Create a new node */ | |
struct Node* newnode = (struct Node*)malloc(sizeof(struct Node)); /* new Node; */ | |
newnode->coeff = coeff; | |
newnode->power = power; | |
newnode->next = NULL; | |
/* If linked list is empty */ | |
if (start == NULL) | |
return newnode; | |
/* If linked list has nodes */ | |
struct Node* ptr = start; | |
while (ptr->next != NULL) | |
ptr = ptr->next; | |
ptr->next = newnode; | |
return start; | |
} | |
/* Functionn to display coefficients from linked list */ | |
void print_coeff(struct Node* ptr, int id) | |
{ | |
printf("Polynomial %d: ", id); | |
while (ptr->next != NULL) | |
{ | |
printf("%s%0.1f,", | |
(ptr->coeff >= 0) ? | |
"+" : "", ptr->coeff); | |
ptr = ptr->next; | |
} | |
printf("%0.1f\n", ptr->coeff); | |
} | |
void clear_list(struct Node* head) | |
{ | |
while (head != NULL) | |
{ | |
struct Node *tmp = head; | |
head = head->next; | |
free(tmp); | |
} | |
} | |
/* Function to add coefficients of two elements having same powerer */ | |
void remove_duplicates(struct Node* start) | |
{ | |
struct Node *ptr1; | |
struct Node *ptr2; | |
struct Node *dup; | |
ptr1 = start; | |
/* Pick elements one by one */ | |
while (ptr1 != NULL && ptr1->next != NULL) | |
{ | |
ptr2 = ptr1; | |
/* | |
Compare the picked element | |
with rest of the elements | |
*/ | |
while (ptr2->next != NULL) | |
{ | |
/* If powerer of two elements are same */ | |
if (ptr1->power == ptr2->next->power) | |
{ | |
/* Add their coefficients and put it in 1st element */ | |
ptr1->coeff = ptr1->coeff + ptr2->next->coeff; | |
dup = ptr2->next; | |
ptr2->next = ptr2->next->next; | |
/* remove the 2nd element */ | |
free(dup); | |
} | |
else | |
ptr2 = ptr2->next; | |
} | |
ptr1 = ptr1->next; | |
} | |
} | |
/* Function to multiply two polynomial numbers */ | |
struct Node* multiply(struct Node* poly1,struct Node* poly2, | |
struct Node* poly3) | |
{ | |
/* | |
Create two pointer and store the | |
address of 1st and 2nd polynomials | |
*/ | |
struct Node *ptr1 = poly1; | |
struct Node *ptr2 = poly2; | |
while (ptr1 != NULL) | |
{ | |
while (ptr2 != NULL) | |
{ | |
double coeff; | |
int power; | |
/* | |
Multiply the coefficient of both | |
polynomials and store it in coeff | |
*/ | |
coeff = ptr1->coeff * ptr2->coeff; | |
power = ptr1->power + ptr2->power; | |
/* | |
Invoke add_node function to create a | |
newnode by passing three parameters | |
*/ | |
poly3 = add_node(poly3, coeff, power); | |
/* | |
move the pointer of 2nd polynomial | |
two get its next term | |
*/ | |
ptr2 = ptr2->next; | |
} | |
/* | |
Move the 2nd pointer to the | |
starting point of 2nd polynomial | |
*/ | |
ptr2 = poly2; | |
/* move the pointer of 1st polynomial */ | |
ptr1 = ptr1->next; | |
} | |
/* | |
this function will be invoke to add | |
the coefficient of the elements | |
having same powerer from the resultant linked list | |
*/ | |
remove_duplicates(poly3); | |
return poly3; | |
} | |
/* | |
This is the function to Calculate Lagrange and displays the coefficient of lagrange | |
*/ | |
struct Node* func(struct Node* poly1, struct Node* poly2, double val) | |
{ | |
/* | |
Create two pointer and store the | |
address of 1st and 2nd polynomials | |
*/ | |
struct Node *ptr1 = poly1; | |
while (ptr1 != NULL) | |
{ | |
double coeff; | |
int power; | |
coeff = ptr1->coeff * val; | |
power = ptr1->power; | |
poly2 = add_node(poly2, coeff, power); | |
ptr1 = ptr1->next; | |
} | |
return poly2; | |
} | |
struct Node* func1(struct Node* poly1, struct Node* poly2) | |
{ | |
/* | |
Create two pointer and store the | |
address of 1st and 2nd polynomials | |
*/ | |
struct Node *ptr1 = poly1; | |
while (ptr1 != NULL) | |
{ | |
double coeff; | |
int power; | |
coeff = ptr1->coeff; | |
power = ptr1->power; | |
poly2 = add_node(poly2, coeff, power); | |
ptr1 = ptr1->next; | |
} | |
return poly2; | |
} | |
void multiplay_polynomial(struct polynomial *poly[], int count, double* arr, int id) | |
{ | |
struct Node *h[count]; | |
struct Node *h1[count]; | |
int i, val = 0, counter = 0; | |
for (i = 0; i < count; i++) | |
{ | |
h[i] = NULL; | |
h1[i] = NULL; | |
} | |
for (i = 0; i < count; i++) | |
{ | |
int j; | |
for (j = 0; j < count-2; j++) | |
{ | |
struct Node *hold = NULL; | |
hold = multiply(poly[i]->p[j], poly[i]->p[j+1], hold); | |
clear_list(poly[i]->p[j+1]); | |
poly[i]->p[j+1] = hold; | |
} | |
h[i] = poly[i]->p[j]; | |
val++; | |
} | |
for (i = 0; i < val; i++) h1[i] = func(h[i], h1[i], arr[i]); | |
for (i = 1; i < val; i++) h1[0] = func1(h1[i],h1[0]); | |
remove_duplicates(h1[0]); | |
print_coeff(h1[0], id); | |
for (i = 0; i < count; i++) | |
clear_list(h1[i]); | |
} | |
double calculate_lagrange(float x[SIZE_ARR], float y[SIZE_ARR], int count, double xp, int display_id) | |
{ | |
/* Store the Result */ | |
double result_sum = 0; | |
/* Array to store value of Lagrange coefficients */ | |
double *coeff_array = (double*)malloc(sizeof(double)*count); | |
double *denominator = (double*)malloc(sizeof(double)*count); | |
struct polynomial *poly[count]; | |
int i, j, counter[count]; | |
for (i = 0; i < count; i++) | |
poly[i] = (struct polynomial*)malloc(sizeof(struct polynomial)); | |
for (i = 0; i < count; i++) | |
{ | |
counter[i] = 0; | |
double result = 1; | |
double c1 = 1.0; | |
for (j = 0; j < count; j++) | |
{ | |
/* Using Formula */ | |
if (j != i) | |
{ | |
struct Node *poly1 = NULL; | |
poly1 = add_node(poly1, 1, 1); | |
poly1 = add_node(poly1, -x[j], 0); | |
poly[i]->p[counter[i]++] = poly1; | |
result = result * (double)((xp - x[j]) / (x[i] - x[j])); | |
c1 *= (x[i] - x[j]); | |
} | |
} | |
denominator[i] = y[i] / c1; | |
coeff_array[i] = result; | |
result_sum += y[i] * coeff_array[i]; | |
} | |
if (display_id >= 0) | |
multiplay_polynomial(poly, count, denominator, display_id); | |
free(coeff_array); | |
free(denominator); | |
for (i = 0; i < count; i++) | |
{ | |
for (j = 0; j < counter[i]; j++) | |
clear_list(poly[i]->p[j]); | |
free(poly[i]); | |
} | |
return result_sum; | |
} | |
float estimate_error(float y_real, float y_estimated) { | |
return (y_real-y_estimated) / y_real * 100; | |
} | |
void signal_handler(int sig) | |
{ | |
if (sig == SIGINT) | |
{ | |
printf("Interrupt signal: %d (CTRL+C)\n", sig); | |
killpg(getpgrp(), SIGINT); /* Interrupt every child before exit */ | |
exit(EXIT_SUCCESS); /* Interrupt is not an error */ | |
} | |
else if (sig == SIGUSR2) | |
{ | |
printf("Terminated by process: %d\n", sig); | |
exit(EXIT_FAILURE); | |
} | |
else if (sig == SIGUSR1) | |
{ | |
/*printf("Received SIGUSR1: %d\n", getpid()); /* debug */ | |
/* | |
Do nothing, Signal handler for SIGUSR1 is registered only to avoid | |
the default action for this signal which will terminate the program. | |
*/ | |
} | |
} | |
int register_signals() | |
{ | |
struct sigaction sact; | |
sigemptyset(&sact.sa_mask); | |
sact.sa_flags = 0; | |
sact.sa_handler = signal_handler; | |
if (sigaction(SIGUSR1, &sact, NULL) != 0) | |
{ | |
fprintf(stderr, "Failed to setup SIGUSR1 handler: %s\n", strerror(errno)); | |
return 0; | |
} | |
if (sigaction(SIGUSR2, &sact, NULL) != 0) | |
{ | |
fprintf(stderr, "Failed to setup SIGUSR2 handler: %s\n", strerror(errno)); | |
return 0; | |
} | |
if (sigaction(SIGINT, &sact, NULL) != 0) | |
{ | |
fprintf(stderr, "Failed to setup SIGINT handler: %s\n", strerror(errno)); | |
return 0; | |
} | |
/* | |
Explicitly setting the disposition of SIGCHLD to SIG_IGN causes | |
any child process that subsequently terminates to be immediately | |
removed from the system instead of being converted into a zombie. | |
*/ | |
sact.sa_handler = SIG_IGN; | |
if (sigaction(SIGCHLD, &sact, NULL) != 0) | |
{ | |
fprintf(stderr, "Failed to setup SIGCHLD handler: %s\n", strerror(errno)); | |
return 0; | |
} | |
return 1; | |
} | |
void mask_and_suspend() | |
{ | |
sigset_t mask, oldmask; | |
sigemptyset(&mask); | |
sigaddset(&mask, SIGUSR1); | |
sigprocmask(SIG_BLOCK, &mask, &oldmask); | |
sigsuspend(&oldmask); | |
sigprocmask(SIG_UNBLOCK, &mask, NULL); | |
} | |
void notify_process(pid_t pid) | |
{ | |
if (kill(pid, SIGUSR1) == -1) | |
{ | |
fprintf(stderr, "Failed to send signal: %s\n", strerror(errno)); | |
exit_failure(); | |
} | |
} | |
void start_calculations(const char *path, int num) | |
{ | |
/* Wait access to start calculation */ | |
mask_and_suspend(); | |
char line[BUFF_MAX]; | |
float x[SIZE_ARR], y[SIZE_ARR]; | |
int length = get_line(path, num, line, sizeof(line)); | |
if (length <= 0) return; | |
int count = get_values(line, length, &x, &y); | |
if (count == 0) | |
{ | |
fprintf(stderr, "Invalid line: %d\n", num); | |
return; | |
} | |
/* First round */ | |
float p = calculate_lagrange(x, y, 6, x[6], -1); | |
snprintf(line, sizeof(line), | |
"%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f," | |
"%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f\n", | |
x[0],y[0],x[1],y[1],x[2],y[2],x[3],y[3], | |
x[4],y[4],x[5],y[5],x[6],y[6],x[7],y[7],p); | |
if (!put_line(path, num, line, strlen(line))) | |
{ | |
fprintf(stderr, "Can not put line: %d (%s)", num, strerror(errno)); | |
return; | |
} | |
/* Notify parent process about finish of first round */ | |
notify_process(getppid()); | |
/* Wait parent to accept second round */ | |
mask_and_suspend(); | |
/* Second round (also print polynomial’s 7 coefficients) */ | |
float r = calculate_lagrange(x, y, 7, x[6], num-1); | |
snprintf(line, sizeof(line), | |
"%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f," | |
"%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f,%.1f\n", | |
x[0],y[0],x[1],y[1],x[2],y[2],x[3],y[3], | |
x[4],y[4],x[5],y[5],x[6],y[6],x[7],y[7],p,r); | |
if (!put_line(path, num, line, strlen(line))) | |
fprintf(stderr, "Can not put line: %d (%s)", num, strerror(errno)); | |
} | |
pid_t spawn_process(const char *path, int num) | |
{ | |
/* | |
On success, the PID of the child process is returned in the parent, | |
and 0 is returned in the child. On failure, -1 is returned in the | |
parent, no child process is created and errno is set appropriately. | |
*/ | |
pid_t pid = fork(); | |
if (pid == -1) | |
{ | |
fprintf(stderr, "Failed to spawn child process: %d (%s)\n", | |
num, strerror(errno)); | |
exit_failure(); | |
} | |
else if (pid == 0) | |
{ | |
/* We are in a child process */ | |
start_calculations(path, num); | |
exit(EXIT_SUCCESS); | |
} | |
return pid; | |
} | |
int main(int argc, char* argv[]) | |
{ | |
if (argc < 2) | |
{ | |
fprintf(stderr, "Please specify file name\n"); | |
return 1; | |
} | |
const char *path = argv[1]; | |
char buffer[BUFF_MAX]; | |
int fd = open(path, O_RDWR, 644); | |
if (fd < 0) | |
{ | |
fprintf(stderr, "Can not open file: %s\n", strerror(errno)); | |
return 1; | |
} | |
int length = read_file(fd, buffer, sizeof(buffer)); | |
if (length <= 0) | |
{ | |
fprintf(stderr, "Can not read file: %s\n", strerror(errno)); | |
return 1; | |
} | |
close(fd); | |
int count = get_line_count(buffer, length); | |
if (count == 0) | |
{ | |
fprintf(stderr, "Invalid or empty file: %s\n", path); | |
return 1; | |
} | |
if (!register_signals()) | |
{ | |
fprintf(stderr, "Failed to register signal handlers: %s\n", path); | |
return 1; | |
} | |
pid_t pids[count]; | |
int i; | |
/* Use only first 8 row */ | |
count = count > 8 ? 8 : count; | |
/* Spawn all child processes */ | |
for (i = 1; i <= count; i++) | |
pids[i-1] = spawn_process(path, i); | |
/* Wait for first round to be completed */ | |
for (i = 0; i < count; i++) | |
{ | |
notify_process(pids[i]); | |
mask_and_suspend(); | |
} | |
int used_count = 0; | |
float err_sum = 0.; | |
for (i = 1; i <= count; i++) | |
{ | |
char line[BUFF_MAX]; | |
float x[SIZE_ARR], y[SIZE_ARR]; | |
int length = get_line(path, i, line, sizeof(line)); | |
if (length <= 0) return 1; | |
int x_count = get_values(line, length, &x, &y); | |
if (x_count == 0) | |
{ | |
fprintf(stderr, "Invalid line: %d\n", i); | |
return 1; | |
} | |
if (x_count > 8) | |
{ | |
/* x[8] is calculated p(x_7) with degree of 5 */ | |
err_sum += estimate_error(y[6], x[8]); | |
used_count++; | |
} | |
} | |
float average_err = (err_sum / used_count); | |
average_err = (average_err < 0) ? -average_err : average_err; | |
printf("Error of polynomial of degree 5: %.1f\n", average_err); | |
/* Run second round in childs */ | |
for (i = 0; i < count; i++) notify_process(pids[i]); | |
/* Wait for second round to be completed */ | |
for (i = 0; i < count; i++) wait(NULL); | |
used_count = 0; | |
err_sum = 0.; | |
for (i = 1; i <= count; i++) | |
{ | |
char line[BUFF_MAX]; | |
float x[SIZE_ARR], y[SIZE_ARR]; | |
int length = get_line(path, i, line, sizeof(line)); | |
if (length <= 0) return 1; | |
int x_count = get_values(line, length, &x, &y); | |
if (x_count == 0) | |
{ | |
fprintf(stderr, "Invalid line: %d\n", i); | |
return 1; | |
} | |
if (x_count > 8) | |
{ | |
/* Actully y[8] is calculated p(x_7) with degree of 6 */ | |
err_sum += estimate_error(y[6], y[8]); | |
used_count++; | |
} | |
} | |
average_err = (err_sum / (float)used_count); | |
average_err = (average_err < 0) ? -average_err : average_err; | |
printf("Error of polynomial of degree 6: %.1f\n", average_err); | |
/* That's all */ | |
return 0; | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment