Last active
December 3, 2020 10:36
-
-
Save hoenirvili/9ac9b8178284f7b739d31a5cee9cb81f to your computer and use it in GitHub Desktop.
euclidian algorithm, linear diophantine and linear_congruential equation solver
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
#include <stdio.h> | |
#include <assert.h> | |
#include <stdlib.h> | |
#include <string.h> | |
#include <stdbool.h> | |
#include <errno.h> | |
typedef struct pair{ | |
int alpha; | |
int beta; | |
}pair; | |
typedef struct er { | |
int gcd; | |
pair pair; | |
}er; | |
static inline int mod(int a, int b) { | |
int m = a % b; | |
if (m < 0) | |
// m += (b < 0) ? -b : b; // avoid this form: it is UB when b == INT_MIN | |
m = (b < 0) ? m - b : m + b; | |
return m; | |
} | |
static inline int max(int a, int b) { return (a>b)? a : b;} | |
// Extended Euclidian algorithm | |
static er euclid(int a, int b) | |
{ | |
er e = {0}; | |
if ((a == 0) || (b == 0)) { | |
a = abs(a); | |
b = abs(b); | |
e.gcd = max(a, b); | |
e.pair.alpha = 1; | |
return e; | |
} | |
bool sw = false; | |
if (a < b) { | |
a ^= b; | |
b ^= a; | |
a ^= b; | |
sw = true; | |
} | |
assert(a > b && a > 0 && b > 0); | |
pair v = {0}; | |
pair v0 = {.alpha = 1, .beta = 0}; | |
pair v1 = {.alpha = 0, .beta = 1}; | |
int q, r = 0; | |
while (b != 0) { | |
q = a / b; | |
r = mod(a, b); | |
a = b; | |
b = r; | |
v = v0; | |
v0 = v1; | |
v1.alpha = v.alpha - q * v1.alpha; | |
v1.beta = v.beta - q * v1.beta; | |
} | |
if (sw) { | |
v0.alpha ^= v0.beta; | |
v0.beta ^= v0.alpha; | |
v0.alpha ^= v0.beta; | |
} | |
e.gcd = abs(a); e.pair = v0; | |
return e; | |
} | |
typedef struct solution{ | |
int x; | |
int y; | |
#define MAX_BUFF 256 | |
char err[MAX_BUFF]; | |
}solution; | |
// Linear Diophantine equations solution finder | |
static solution diophantic(int a, int b, int c) | |
{ | |
solution s = {0}; | |
er e = euclid(a,b); | |
// if gcd | c, exits a k from Z that c = gcd * k | |
if (c % e.gcd) { | |
snprintf(s.err, MAX_BUFF, | |
"No integer solution for %d and %d\n", a, b); | |
return s; | |
} | |
c = c / e.gcd; | |
s.x = e.pair.alpha * c; | |
s.y = e.pair.beta * c; | |
return s; | |
} | |
#define MAX_SOL 50 | |
static void linear_congruential(int a, int b, int m) | |
{ | |
int s = 0; | |
er e = euclid(a, m); | |
if (m % e.gcd) | |
printf("No Linear congurential solutions for a = %d b = %d m = %d", | |
a, b, m); | |
int b0 = b / e.gcd; | |
int x0 = e.pair.alpha * b0; | |
printf("For our linear congurential a = %d, b = %d, m = %d Solutions: ", | |
a,b,m); | |
for (int i = 0; i < e.gcd; i++) { | |
s = mod((x0 + (i * m) / e.gcd), m); | |
printf("%d ", s); | |
} | |
printf("\n"); | |
} | |
typedef struct it { | |
int a; | |
int b; | |
}it; | |
typedef struct dit{ | |
int a; | |
int b; | |
int c; | |
}dit; | |
#define ARRAY_SIZE(arr) \ | |
(sizeof(arr) / sizeof((arr)[0])) | |
int main(void) | |
{ | |
it its[] = { | |
{0, 5}, | |
{13, 0}, | |
{21, 13}, | |
{311, 11}, | |
{6239871, 71164}, | |
{3252, 32}, | |
{31, 3}, | |
{9999, 33}, | |
{125612,2351}, | |
{821,512}, | |
{71,11}, | |
{17, 11}, | |
{195, 42}, | |
}; | |
printf("\nExtended Euclidian\n"); | |
it one = {0}; | |
for (size_t i = 0; i < ARRAY_SIZE(its); i++) { | |
one = its[i]; | |
er e = euclid(one.a, one.b); | |
printf("a = %d, b = %d, gcd = %d, alpha = %d, beta = %d\n", | |
one.a, one.b, e.gcd, e.pair.alpha, e.pair.beta); | |
} | |
dit d[] = { | |
{21, 13, 5}, | |
{97, 35, 13}, | |
{195, 42, 12}, | |
{47, 30, 1}, | |
}; | |
printf("\nLinear Diophantine equations\n"); | |
dit done = {0}; | |
for (size_t i = 0; i < ARRAY_SIZE(d); i++) { | |
done = d[i]; | |
solution s = diophantic(done.a, done.b, done.c); | |
if (strlen(s.err) != 0) { | |
printf("%s\n",s.err); | |
continue; | |
} | |
printf("Solution for a = %d, b = %d, c = %d, are x = %d, y = %d\n", | |
done.a, done.b, done.c, s.x, s.y); | |
} | |
printf("\nLinear congurential\n"); | |
linear_congruential(5, 25, 10); | |
linear_congruential(17, 3, 29); | |
linear_congruential(6, 4, 10); | |
linear_congruential(17, 1, 43); | |
printf("\n"); | |
} |
If you can even do that on gist. Cuz I can't open a pull request.
@VladuZ
Wait a second, let me review the changes first.
@VladuZ
Thanks for the patch, it's now merged.
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Merge my fork <3